acer

9 years, 149 days


These are answers submitted by acer

As Carl points out, the performance regression seems due to plot3d no longer using evalhf for procedure J. That is a shame, and I don't see a rationale for it so I'd call this regression a bug.

I'm not sure what Yiannis means by saying the values like lim, lambda, etc, have to be globals. It's easy enough to burn them into a generated procedure, though I note that doing so does not by itself make Maple 18.02 utilize evalhf.

By the way, in Maple 18 the colors from Maple 13 can be duplicated by using the lightmodel=none option, though it looks ugly when inlined into this sight, for some strange reason. I notice also that in the session that it's computed it has the matching colors, but after saving then when reopened in Maple 18 it initially appears different and applies a default lightmodel (until recomputed, sigh).

On my Windows 7 machine with 64bit Maple 18.02 it takes about 14 seconds to produce this plot using evalhf around a wrapped procedure call.

Download Euler_modif.mw

Using the Compiler package it is possible to make both the 3D plot and the 2D densityplot compute even faster. I also separated out the real and imaginary parts, while doing so. (That separation, and even some amount of shoehorning custom formula into a template of procedure J can be programmatic, but making it automatic can be tricky. Here, I did some manual copy and pasting after a few evalc computations.)

The end 3D plot appears the same as above, but computes in a about 3 seconds on my machine.

Euler_compiled.mw

It would be even faster to write an evalhf'ble or compilable procedure which would act inplace on a float[8] Array argument and populate it inplace using the iteration scheme. Such an Array could then be stuffed into a GRID substructure of a PLOT3D structure. It's often much faster to populate an Array all within a single procedure call than to have to invoke a scalar-valued procedure for each entry. (That's how Maple 18's Fractals:-EscapeTime package's commands do it.)

acer

My guess is that you are in the habit of sometimes reexecuting out of order some commands in your Worksheet/Document. That is particularly problematic when you utilize assume. The code below mimics the kind of trouble you can get into when you repeatedly execute calls to assume and assignments with the assumed names.

restart:

interface(showassumed=1):

U := x*P^2/3;
                                      1    2
                                 U := - x P 
                                      3     

assume(P>0);

diff(U,P);

                                    2    
                                    - x P
                                    3    

U := x*P^2/3;

                                      1    2
                                 U := - x P 
                                      3     

assume(P>0);

diff(U,P);

                                      0

Such endeavors are often easier and less problematic if you utilize assuming instead of assume. You can even make an assignment like,

   conds := P>0, Q>0:

and then later on do computations to take that into effect.  Eg,

   some_commmand(something) assuming conds;

 

acer

value( sum( Sum(...) ) )

February 25 2015 acer 10867

You mentioned at one point that you prefered not so use uneval-quoted 'sum'(...) for the inner summation because then you wouldn't get the summation symbol. I suppose that you might feel similarly about using a nested call to add.

It seems that for this example another workaround variant is to use inert Sum for the inner summation. That typesets with a  summation symbol, albeit in gray.

In 1D Maple Notation that could look like this.

X := [0, 1, 2, 3]:
n := 2:
value( sum(Sum(a[k]*X[i]^k, i = 1 .. nops(X)), k = 0 .. n) );

                           4 a[0] + 6 a[1] + 14 a[2]

In 2D Input it could look as follows. Note that I used command-completion (Ctrl-Shift-Spacebar on my Linux) after typing the word Sum and chose the command "template" for inert definite summation from the popup menu offerings.

I also keystrokes Ctl-Shift-underscore in my Maple to get the subscripted display of X[i] and a[k]. You could also just enter X[i] and a[k] to be displayed as indexed names, if you prefer. In earlier version of Maple getting 2D Input for an indexed name involved a different set of keystrokes.

X := [0, 1, 2, 3]:

n := 2:

value(sum(Sum(a[k]*X[i]^k, i = 1 .. nops(X)), k = 0 .. n))

4*a[0]+6*a[1]+14*a[2]

Download typesetsum.mw

You also wondered why you couldn't place uneval-quotes around the entire inner summation, when using sum. For this example it seems that you can, if you handle nops(X) slightly differently. (The nops(X) looks ugly as upper index value on the summation symbol, anyway, in my opinion.) I used Maple 18.01 here too. In 1D Maple Notation,

X := [0, 1, 2, 3]:
n := 2:
r := nops(X):
eval( sum('sum(a[k]*X[i]^k, i = 1 .. r)', k = 0 .. n) );

                           4 a[0] + 6 a[1] + 14 a[2]

And in 2D Input using command-completions (or palette for typeset definite sum) the following has both summation symbols in black. Of course the uneval-quotes detract somewhat.

X := [0, 1, 2, 3]:

n := 2:

r := nops(X):

eval(sum('sum(a[k]*X[i]^k, i = 1 .. r)', k = 0 .. n))

4*a[0]+6*a[1]+14*a[2]

Download typesetsum2.mw

 

acer

something

February 25 2015 acer 10867


restart:

f := x^2+y^2 = 1:

# One way is to select only those solutions which `is` can determine
# produces a positive `y`, given assumptions on `x` (that we observe
# through manual inspection).

op(select(e->is(eval(y,e),positive), [solve( {f}, y )])) assuming x>-1, x<1;

{y = (-x^2+1)^(1/2)}

# Another way, but somewhat similar.

solve( {f, y>0}, y, useassumptions ) assuming x>-1, x<1;

{y = (-x^2+1)^(1/2)}

# The following gets us information about key x values,
# without having to deduce and supply them manually.
#
# The resulting piecewise is awkward. The empty list is presumably
# there in order to distinguish between NULL and `undefined` (the
# latter having some special connotations in Maple).

S := solve( {f, y>0}, y, parametric );

S := piecewise(x <= -1, [], x < 1, [[y = sqrt(-x^2+1)]], 1 <= x, [])

# This next has a rather artifical look to it,
# since it involves typing in what one observes about
# the conditions of the piecewise S.

op(simplify(S)) assuming x>-1, x<1;

[y = (-x^2+1)^(1/2)]

# The following picks out values from piecewise S and associated
# conditions under which the values would hold.
# This particular code to do that is graceless and doubtless
# not robust. I expect it could be done much better. The idea is
# to associate with the ith value in the pieceise the logical
# conjunction of the ith condition and the negation of conditions
# that preceded it.

K:=convert(subsindets(S,listlist,op),list):
W:=seq( `if`(K[2*i]<>[],
             [op(K[2*i]),
             (%And(seq(subsindets(
                        subsindets(
                          subsindets(
                            subsindets(K[2*j-1],
                                       `<=`,u->(op(1,u) &> op(2,u))),
                            `<`,u->(op(1,u)>=op(2,u))),
                          specfunc(`&>`),u->(op(1,u)>op(2,u))),
                        `=`,u->(op(1,u)<>op(2,u))),
                      j=1..i-1),
                  K[2*i-1]))],
             NULL), i=1..nops(K)/2 );

[y = (-x^2+1)^(1/2), %And(-1 < x, x < 1)]

lprint(%);

[y = (-x^2+1)^(1/2), %And(-1 < x, x < 1)]

op(simplify(S)) assuming value(W[2]);

[y = (-x^2+1)^(1/2)]

 


Download solvegrin2.mw

acer

undocumented

February 24 2015 acer 10867

In Maple 18 a new 2D plot option size was introduced. This can be used with the 2D plot command and with 2D commands in the plots package, and when using the plots:-display command on 2D plots.

In Maple 17 the Standard Java GUI already understood an earlier form of the relavent substructure (or a PLOT structure) related to sizing. But the usual documented 2D plotting commands mentioned above did not have any means of inserting this substructure. But it could still be accomplished in Maple 17 (for just the choice of posint pixel numbers width-by-height) using the undocumented Plot:-Structure:-SetSize command.

And there are limitations. The result below of calling SetSize cannot subsequently be combined with other 2D plots or otherwise generally adjusted using the plots:-display command. You'd have to first do all the combining, and then make the call to SetSize be the very last action before finally doing a simple print call on the resized PLOT structure.

Also, I don't think that the File>Print mechanism of the Maple 17's GUI's menubar understands the structure. But pdf Export might...

 

restart:

kernelopts(version);

`Maple 17.02, X86 64 WINDOWS, Sep 5 2013, Build ID 872941`

P := plot(sin(x), x=-2*Pi..2*Pi, gridlines=false):

P;

Plot:-Structure:-SetSize(P,600,200);

 


Download m17setsize.mw

In Maple 18 one should definitely not go this route. Use the documented size option instead.

acer

Grid

February 20 2015 acer 10867

The Grid package would provide more insulation.

acer

LibraryTools

February 20 2015 acer 10867

Put the procedures and computed results into a Library archive (.mla extension).

Build that, and put things into it, and extract from it, using the LibraryTools package.

acer

angles

February 20 2015 acer 10867

In a followup Comment to Carl's Answer you asked about minimal and maximal angles which would attain the right-center field home run, given the two initial velocities as mentioned in your original Question.


restart:

interface( warnlevel = 0 ):

m := 0.145:      #mass of baseball

g := 9.81:       #gravity

C_d := 0.35:     #drag coefficient

r := 0.037:      #radius

y0 := 1.5:       #y(0): height from which ball is hit

yf := 2.5:       #height of home-run fence

xf := 114:       #right-center field home-run distance

v_x := diff( x(t), t ):

v_y := diff( y(t), t ):

eqx := diff( v_x,t) = -((C_d)*rho*Pi*(r^2)*(v_x)*sqrt((v_x)^2 +(v_y)^2))/(2*m):

eqy := diff( v_y,t) = -((C_d)*rho*Pi*(r^2)*(v_y)*sqrt((v_x)^2 +(v_y)^2))/(2*m)-g:

ics := x(0) = 0, y(0) = y0, D(x)(0) = v0*cos(alpha), D(y)(0) = v0*sin(alpha):

Sol := dsolve( { eqx, eqy, ics  }, numeric, parameters = [ v0, rho, xhalt, alpha ],

               events = [ [x(t) - xhalt, halt] ], output = listprocedure ):

Y := eval( y(t), Sol ):

HeightatHalt := proc( v0, rho, xtarget, angle )

   Y( parameters = [ :-v0 = v0, :-rho = rho, :-xhalt = xtarget, :-alpha = angle ] );

   Y( 200 ); # 200 second airborne limit

end proc:

# smallest v0 which will attain right-center home run at
# sea level air resistance and angle alpha = Pi/4.

v_sealevel := fsolve( 'HeightatHalt'( v0, 1.2, xf, Pi/4 ) = yf, v0 = 30..100 );

46.65980695

# smallest v0 which will attain right-center home run at
# Denver air resistance and angle alpha = Pi/4.

v_denver := fsolve( 'HeightatHalt'( v0, 1.2*0.9, xf, Pi/4 ) = yf, v0 = 30..100 );

44.92447287

 

with(Optimization):

# smallest angle to get home run with v0 = v_sealevel
# and sea level air resitance.

alpha_low[v_sealevel,1.2]
   := fsolve('HeightatHalt'( v_sealevel, 1.2, xf, alpha ) - yf = 0,
             alpha=0..Pi/4 );
#  := rhs( Minimize( -alpha
#                    + ('HeightatHalt'( v_sealevel, 1.2, xf, alpha ) - yf)^2,
#                    method = nonlinearsimplex )[2,1] );

.6216930153

# smallest angle to get home run with v0 = v_sealevel
# and Denver air resitance.

alpha_low[v_sealevel,1.2*0.9]
   := fsolve('HeightatHalt'( v_sealevel, 1.2*0.9, xf, alpha ) - yf = 0,
             alpha=0..Pi/4 );
#  := rhs( Minimize( -alpha
#                    + ('HeightatHalt'( v_sealevel, 1.2*0.9, xf, alpha ) - yf)^2,
#                    method = nonlinearsimplex )[2,1] );

.5285125608

# largest angle to get home run with v0 = v_sealevel
# and Denver air resitance.

alpha_high[v_sealevel,1.2*0.9]
   := fsolve('HeightatHalt'( v_sealevel, 1.2*0.9, xf, alpha ) - yf = 0,
             alpha=Pi/4..Pi/2 );

.8947974228

# smallest angle to get home run with v0 = v_denver
# and Denver air resitance.

alpha_low[v_denver,1.2*0.9]
   := fsolve('HeightatHalt'( v_denver, 1.2*0.9, xf, alpha ) - yf = 0,
             alpha=0..Pi/4 );
#  := rhs( Minimize( -alpha
#                    + ('HeightatHalt'( v_denver, 1.2*0.9, xf, alpha ) - yf)^2,
#                    method = nonlinearsimplex )[2,1] );

.6406089003

with(plots):

display(

  # Pi/4 was originally used to find v_sealevel, v_denver
  plot( Pi/4, v0=43..50, color=black ),

  plot( [v_denver, alpha, alpha=0..Pi/4],
        color=blue, linestyle=dot ),

  plot( alpha_low[v_denver,1.2*0.9], v0=43..v_denver,
        color=blue, linestyle=dash ),

  implicitplot( 'HeightatHalt'( v0, 1.2*0.9, xf, alpha ) - yf,
                v0=43..50, alpha=0..Pi/2, gridrefine=2,
                color=blue ),

  plot( [v_sealevel, alpha, alpha=0..alpha_high[v_sealevel,1.2*0.9]],
        color=red, linestyle=dot ),

  plot( alpha_low[v_sealevel,1.2], v0=43..v_sealevel,
        color=red, linestyle=dash ),

  plot( alpha_low[v_sealevel,1.2*0.9], v0=43..v_sealevel,
        color=red, linestyle=dash ),

  plot( alpha_high[v_sealevel,1.2*0.9], v0=43..v_sealevel,
        color=red, linestyle=dash ),

  implicitplot( 'HeightatHalt'( v0, 1.2, xf, alpha ) - yf,
                v0=43..50, alpha=0..Pi/2, gridrefine=2,
                color=red ),

  axis[2] = [tickmarks=piticks], # or use decimalticks
  view = [43..50, 0..3*Pi/8],
  labels = [v0, alpha]
              );

 


Download bballfull.mw

Sorry about the gridlines in that plot, they are an unfortunate artefact of the MapleNet server currently used by MaplePrimes. They don't appear in the attachment.

acer

mechanism

February 20 2015 acer 10867

As far as I know Maple doesn't have a prebuilt mechanism by which general computation can be done via the GPU.

Do you have some reason to think that it does?

acer

identity

February 19 2015 acer 10867


restart:

baseident := isolate(convert(FunctionAdvisor(LaguerreL,"identities",quiet)[6],`global`),
                     LaguerreL(a-2,b,z));

LaguerreL(a-2, b, z) = (LaguerreL(a, b, z)*a-(b-z+2*a-1)*LaguerreL(a-1, b, z))/(1-a-b)

ident := subs(a=n,b=l+1/2,z=y,baseident);

LaguerreL(n-2, l+1/2, y) = (LaguerreL(n, l+1/2, y)*n-(l-1/2-y+2*n)*LaguerreL(n-1, l+1/2, y))/(1/2-n-l)

ee := (2*l+2*n-1)*LaguerreL(n-2,l+1/2,y)-(2*l+4*n-2*y-1)*LaguerreL(n-1,l+1/2,y);

(2*l+2*n-1)*LaguerreL(n-2, l+1/2, y)-(2*l+4*n-2*y-1)*LaguerreL(n-1, l+1/2, y)

simplify(subs(ident,ee));

-2*LaguerreL(n, l+1/2, y)*n

simplify(ee,{ident});

-2*LaguerreL(n, l+1/2, y)*n

kernelopts(version);

`Maple 16.02, X86 64 LINUX, Nov 18 2012, Build ID 788210`

 


Download LaguerreL_ident.mw

acer

reinstall

February 19 2015 acer 10867

There have been several reports of the same error message, specific to Maple 18 (or an update of same). It seems that a reinstall usually fixes the problem.

I don't know whether the problem is with corruption of the installer when it's downloaded (in which case you'll have to download again), or whether the problem is at the installation step from a viable download (in which case you just need to retry with the same installer).

acer

something

February 19 2015 acer 10867

It's a bit of a chore, but one of the identities that FunctionAdvisor reports for LegendreP seems to be useful. I extracted that and assigned it as identityA below.

Also, to convert back from LegendreP to SphericalY I had to invert the reverse conversion, which I assigned as identityB below.

It's possible that there is a relation similar to identityA which could be applied to the expression while still in SphericalY form. And that might remove the need for identityB. But the SphericalY form seems to need even more combinations of simplification steps.

Not all the simplify(...,size) calls wrapped around simplify calls may be necessary, but it keeps some of the intermediate expressions a little tighter and easier to understand I think.

I ran this on Maple 18.01.

restart:

with(VectorCalculus):
SetCoordinates(spherical[r,theta,phi]):

# Originally I used these assumptions for all simplification steps.
# But it seems as if they are not needed...
#conds := lambda::posint, mu::integer, abs(mu)<=lambda:

ee := Laplacian(SphericalY(lambda,mu,theta,phi)):

identityA := subs(a=lambda+2,b=mu,z=cos(theta),
                  convert(FunctionAdvisor(LegendreP,"identities",quiet)[6],`global`)):

identityB := isolate(SphericalY(lambda,mu,theta,phi)
                       = convert(SphericalY(lambda,mu,theta,phi),LegendreP),
                     LegendreP(lambda,mu,cos(theta))):

combine(simplify(simplify(convert(ee,LegendreP)),size),radical):
simplify(simplify(%),size):

simplify(simplify(subs(identityA,%),size)):

normal(simplify(subs(identityB,%)));

            lambda (lambda + 1) SphericalY(lambda, mu, theta, phi)
          - ------------------------------------------------------
                                       2                          
                                      r                                            

acer

pairs := convert(M,listlist);

acer

multiplication

February 18 2015 acer 10867

Two instances of an adjacent pair of bracketed terms are being interpreted as function application, so if you intended them as products then you left out two instances of multiplcation syntax.

The simplest thing, in that case, would be to insert a pair of multiplication symbols.

solve({(1-x)/(1+b) = b*y, (1-y)*(1-a)/(-a*b+1) = b*z, (1-z)*(1-b)/(-a*b+1) = a*x}, {x, y, z});

acer

Maple 13

February 17 2015 acer 10867

If the problem is that the grid shown on the rendered surface is too dark in your Maple 13 then how about using one of these option choices for your call to plots:-display?

style = patchnogrid

or

style = wireframeopaque

acer

1 2 3 4 5 6 7 Last Page 1 of 107