acer

32405 Reputation

29 Badges

19 years, 351 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are answers submitted by acer

Using odeplot and a single dsolve (numeric) call -- with parameters -- this can be done somewhat efficiently.

(I compute it also with Carl's Answer's phaseportrait code, to illustrate that odeplot is reasonably quick. There's about a 10-fold timing speedup on my Maple 2020 on Linux, for this number of different Initial Conditions. The performance speedup increases with the number of ICs.)

 

First, using plots:-odeplot and a dsolve solution with parameters.

 

restart;

str:=time[real]():
dsol:=dsolve({diff(r(t),t)=r(t)^2*(r(t)+sin(theta(t))),
              diff(theta(t),t)=r(t)^3*(r(t)+theta(t)+cos(theta(t))),
              theta(0)=partheta, r(0)=parr}, numeric,
              parameters=[partheta, parr]):

odeplotpar:=proc(thetapar,rpar,i)
  dsol(parameters=[thetapar,rpar]);
  plots:-odeplot(dsol, [:-theta(t),:-r(t)], t=0..0.5,
                 numpoints=ceil(0.5/0.0005),
                 color=COLOR(HSV,0.85*i/16,0.85,0.7));
end proc:

P:=plottools:-transform((t,r)->[r*cos(t),r*sin(t)])(
 plots:-display(
   plots:-fieldplot([r^3*(r+theta+cos(theta)),r^2*(r+sin(theta))],
                    theta=0..2*Pi, r=0.7..3, fieldstrength=fixed(0.3),
                    color=COLOR(HSV,0,0,0.85)),
   seq(odeplotpar(i*Pi/8,0.75,i), i=0..15),
   axiscoordinates=polar, size=[600,600],
   axis[2]=[tickmarks=[subticks=0]],
   axis[1]=[gridlines=false, tickmarks=[[seq(0..3,0.5)]]])):
time[real]()-str;

.200

P;

 

For a timing comparison, a similar plot done with DEtools:-phaseportrait.
(One can also compare with DEplot, etc).

 

restart;

str:=time[real]():
P:=plots:-display( #needed to get polar axes
    plottools:-transform((t,r)-> r*~[cos,sin](t))( #transform to polar
        DEtools:-phaseportrait(
            [
                diff(r(t),t) = r(t)^2*(r(t)+sin(theta(t))),
                diff(theta(t),t) = r(t)^3*(r(t)+theta(t)+cos(theta(t)))
            ],
            [theta,r](t),
            t= 0..0.5, #independent variable range for all trajectories
            (i0:= map2(`~`[`=`], [theta,r](0), `[]`~(Pi/8*~[$0..15], .75))), #inits
            linecolor= [seq(COLOR(HSV, .85*i/ni, .85, .7), i= 0..(ni:= nops(i0)-1))],
            thickness= 1, #trajectory thickness (min. is 0, not 1)
            method= rkf45, stepsize= 0.0005,
            color= COLOR(HSV, 0, 0, .85), #light-grey arrows
            arrowsize= .3
        )
    ),
    axiscoordinates= polar, size= [600$2],
    axis[2]= [tickmarks= [subticks= 0]], axis[1]= [gridlines= false]
 ):
time[real]()-str;

1.949

#P;

 

Download polar_odeplot.mw

Here are two possible avenues ideas: 1) percent as a Unit, and 2) percent as an evalf'able constant. (I'll leave out the potential rabbit hole of some kind of object.)

2D Input mode seems to be a problem. And if that's the principle goal then this won't be worth much.

Note that in the Maple 2020 GUI the relevant output below renders OK for me, with an actual percent symbol even for the alternate. The Mapleprimes inline display doesn't seem to handle it, so a blank box appears below.

I could also have used the entity `%` for the Units variant. I used the so-called "Arabic Percent Symbol" only because it looked a little different. Both seem to be limited, and not work as 2D Input.

restart;

 

Entering % in 2D Input mode is problematic. If I use `%` or another

encoding of the usual % symbol then as 2D Input it gets parsed as the

so-called ditto (last result).

 

It's difficult even if I use an encoding for Arabic Percent Sign. (There are

other alternates, but my Maple 2020.0 on Linux doesn't render them

in output. (It shows only an empty box for the symbol.) In decimal

that can be referred to as the entity `٪` and in hex `٪`.

This too does not serve for 2D Input because the Punctuation palette

has only the problematic regular %. And cut&paste of this alternative

turns into the problematic usual % symbol upon parsing.

 

So, I don't know how to use a pretty-printed % for 2D Input. If that is

solved then one of these two mechanisms might be improved.

 

Units:-AddUnit('percent', 'symbol' = '`٪`',
               'prefix' = 'SI', 'conversion' = 1/100);

P := 45*Unit(percent);

45*Units:-Unit(Ùª)

convert( P * 123/7, units, 1 );

1107/140

convert( P * 123.0/7, units, 1 );

7.907142856

 

with(Units:-Simple):

simplify( P * 123/7, units );

1107/140

simplify( P * 123.0/7, units );

7.907142857

123/7*Unit(meter) * 45*Unit(percent);

evalf(%);

(1107/140)*Units:-Unit(m)

7.907142857*Units:-Unit(m)

restart;

`evalf/constant/%` := proc(x) 0.01; end proc:

`evalf/constant/percent` := proc(x) 0.01; end proc:

 

P := 45 * `%`;

45*%

evalf( P );

.45

evalf( P * 123/7 );

7.907142857

Q := 45 * percent;

45*percent

evalf( Q );

.45

evalf( Q * 123/7 );

7.907142857

 

Download percent.mw

[edited] As far as I know, Units and evalf-extensions either do some inner manipulations with locals or rely upon global names. So I don't see offhand and attractive way to make either of these two mechanisms work with a top-level local % declaration.

For example, the following allows for easy 2D Input, but in order to get the final float value it also needs the convert to global (which makes it awkward to use).

restart;

local %:

`evalf/constant/%` := proc(x) 0.01; end proc:

 

 

P := 45*%

45*%

evalf( convert(P,`global`) * 123/7 );

7.907142857

 

Download percent_local_evalf.mw

The awkwardness might all be hidden behind an object, so that entry could be something easy like, say %() as an empty function call. And then the prettyprint could be made nice and clean, and a static evalf export could do the conversion via global, and regular prettyprint output could hide the function call brackets. Oof.

I notice that you cut down your y-range to 0.0..0.7 in your original 3D plot. If you had used the range y=0..1 then in your Maple 2018 there would be a ragged edge along one side of the surface.

For y>1/sqrt(2)=0.707... you can still get a decent rendering of the surface or the contours.  See below, for the higher y-contour values. (Some of this is automatically better in Maple 2020, btw.)

You can also get a little tighter near to x=1, where you had used 0.99 for the upper edge. (This aspect is not fully better, automatically, in Maple 2020.0.)

The effect of those changes can also bring down the surface edge towards z=0 along boundaries, ie. a more faithful representation IMO.

For the 2D curves it's a little trickier to get tight up to the lower x-value (within x=0..1) for y-values less than 1/sqrt(2). One can fiddle with the lower plot range (eg. a conditional and a formula), but roundoff error makes even that tricky. Easier is to allow adaptive 2D plotting and higher numpoints handle it.

Also, you can raise the 2D curves up to their contour values, to include in a 3D plot, if you'd like.

You can also get a variant of that using the contourplot3d command (with axis relabelling). Note that with the parametric calling sequence you can rearrange which formulas get used for each coordinate, and that can be matched by corresponding axis relabelling and reorientation.

restart

kernelopts(version)

`Maple 2018.2, X86 64 LINUX, Nov 16 2018, Build ID 1362973`

with(plots)

z := proc (x, y) options operator, arrow; (x^2+x*y+y^2)*sqrt((1-y^2)^2/(x+y)^2-(x-y)^2) end proc

S := [solve({((-y^2+1)/(x+y))^2-(x-y)^2 = 0}, {y})]; S1 := eval(y, S[1])

eps := 0.1e-8; A := plot3d([x, y, z(x, y)], x = 0 .. 1-eps, y = 0 .. S1-eps, style = surface)

Colors := [yellow, orange, red, green, blue, khaki, violet, cyan, pink, magenta]; P3d := plots:-display(seq((plottools:-transform(proc (a, b) options operator, arrow; [a, y, b] end proc))(plot([x, z(x, y), x = 0 .. 1-eps], color = Colors[1+trunc(y/(.1))], thickness = 3, numpoints = 1000)), y = 0 .. .9, .1))

plots:-display(A, P3d, labels = [x, y, z])

plot([seq([x, z(x, y), x = 0 .. 1-eps], y = 0 .. .9, .1)], color = Colors, thickness = 2, numpoints = 1000, legend = [seq(y = C, C = 0 .. .9, .1)], labels = [x, z])

display(plot3d([x, z(x, y), y], x = 0 .. 1-eps, y = 0 .. S1-eps, style = surface), contourplot3d([x, z(x, y), y], x = 0 .. 1-eps, y = 0 .. S1-eps, contours = [seq(C, C = 0 .. .9, .1)], grid = [151, 151], thickness = 3), orientation = [-100, -15, 35], labels = [x, z, y])

 

Download Contour_Plot_MaplePrimes_ac.mw

Try passing the option  style=surface  to the plot3d call in question.

(Or, if your Maple version is truly ancient,  style=patchnogrid ).

The visual effect that you are calling "grid" is sometimes referred to as a wireframe. Its resolution details (when present in the rendering) are affected by the choice of the grid option, yes. But whether it appears at all -- or not -- is controlled by the style option.

See the Help page with Topic  plot3d,options  for more detail on the possible values that can be used for the style option.

If you search this site then you can find previous questions with the same topic.

For example, here.

I found that recent item by looking at the questions with the phaseportrait tag. It was right at the top.

There are other ways to search the site. Entering the terms in this sites Search bar produces a collection.

You can use multiple initial-conditions with the DEplot command. There are examples on its Help page.

You could also use the odeplot command.

I'm sure that someone will be quite willing to do all your problems in full (even coursework questions). But perhaps you could also show how far you've gotten so far.

Is this something like the look you're going for? (You can adjust the style, transparency, etc). Q1_ac.mw

If you don't need the solid boundary of the leading side (that moves) then it's a great deal simpler to do. I just prefer this look. I find that too much transparency makes these kinds of rotations overly difficult to visualize (and poorly so).

restart;
Lift3D:=plottools[transform]((x,y,z)->[x,z+sqrt(6-(sqrt(x^2+y^2))),y]):
F:=plottools[transform]((x,y)->[x,y+sqrt(6-x),0])(plot(sqrt(4-r)-sqrt(6-r),r=0..4,
                                                   filled,transparency=0.0)):
P := A -> plots:-display(orientation=[-90,15,-160],
            Lift3D(plot3d([[r,theta,sqrt(4-r)-sqrt(6-r)]],
                   r=0..4, theta=0..A, coords=cylindrical, filled,
                   orientation=[-90,15,-160], lightmodel=light1,
                   grid=[40,max(2,ceil(evalf(49*A/(2*Pi))))],
                   scaling=constrained, axes=normal, style=patch)),
            `if`(A<6.28,[F,plottools[rotate](F,-A,[[0,0,0],[0,1,0]])][],NULL),
            labels=[x,y,z]):
plots[animate](P, [A], A=0..2*Pi, frames=17, paraminfo=false);

This may be similarly resource intensive as an animated use of 3D plots:-shadebetween (for a similar number of frames). I could probably improve it by ditching the filled option in the plot3d call (inside P), and using a grid=[...,2] cylinder instead. That's done pretty well before -- the cylinder shape can have varying height, to joint the two surfaces where they don't intersect within the desired domain. Let me know if memory/GUI performance is a huge issue.

Here is a revision, which works.

I put a few comments into the Action Code, to explain how the substitution (for n) in the formula can be done.

I did not otherwise change the algorithm. It could be done more efficiently. I left it in the spirit of your prior work.

Iterative_components_ac.mw

[Edit] This is better than my original example (which passed n=3.0 and while relevant was not the most precise as explanation). I re-instated this Answer as it was accidentally deleted by someone else.

restart;
Digits := 20:
eval( ((n+1.)/n)^(n+1.), n=3 );

         3.1604938271604938268

restart;
Digits := 20:
evalf( eval( ((n+1)/n)^(n+1), n=3 ) );

         3.1604938271604938272

The difference between seq and $ is one of evaluation. (See also `||` versus cat.)

By the way, in the construction of your a list there is a single call to procedure u, into the result of which are placed the ten different integer arguments. In the construction of your b list there are ten calls to procedure u (each of which applies evalf).

This is difficult to test properly and carefully because of the following:
1) So-called "plot persistance", where an Output region of the GUI retails details of some "global plotting options" unless the output is fully removed.
2) Some additonal ordering issue that is hard to pin down -- possibly GUI internal state.

Having said that... I believe that the problem might be overcome by forcing the order in which the ORIENTATION and PROJECTION substructures appear within the PLOT3D structure.

Execute the following within a completely new and blank Worksheet. Note the use of separate Execution Groups. proj_issue_2019.mw

restart;

# new Execution Group
redo := proc(P::specfunc(PLOT3D)) local s,r;
  (s,r):=:-selectremove(:-type,[op(P)],'specfunc({:-PROJECTION,:-ORIENTATION})');
  :-PLOT3D(r[],:-sort(s)[]);
end proc:

# new Execution Group
P1 := plot3d(x*y, x=0..10, y=0..10, grid=[2,2], projection = 0.2, orientation = [55, 75, 0]):
P1;

# new Execution Group
P2 := redo(P1):
P2;

# new Execution Group
op(P1);
op(P2);

For me, the forced reordering fixed the issue.

However, I found that with sufficient jiggery-pokery-editing-rerunning I could somehow fool the GUI into rendering the previously incorrect plot as OK, even though I was taking care about point 1) above. And I know what I'm doing. Hence I conclude that there is some additional weirdness in the GUI state that can also make it magically appear fixed for all Exec Groups (ie, point 2 above).

I am seeing the bug in Maple 2017.2 to 2020.0 on 64bit Linux. But it works properly in Maple 17.02, 18.02, 2015.2 and 2016.2, in which versions it does not make any difference in which order the substructures appear.

I have submitted a bug report (Software Change Request).

Finally, it appears to not be sufficient to simply switch the order in which the projection and orientation options are supplied to the plot3d command. Sometimes that works, and sometimes not -- but do not be deceived into thinking it is sufficient, because you may be fooled by points 1) and 2) above.

You could look at combinat:-nextcomb.

That allows you to process/walk the combinations one by one without generating them all up front. But it also allows you to choose the size of the set.

If you have an original set of N things which are not simply integers then you can apply the nextcomb command to produce the combinations of the ordinals. You can use that to index into the original set of things, to get the appropriate subsets.

restart;
with(combinat):

S := {a,b,c,d,e,f,g}:

N := nops(S);
                             N := 7

Q := firstcomb(N, 3);
                         Q := {1, 2, 3}
S[[Q[]]];
                           {a, b, c}

Q := nextcomb( Q, N );
                         Q := {1, 2, 4}
S[[Q[]]];
                           {a, b, d}

Q := nextcomb( Q, N );
                         Q := {1, 2, 5}
S[[Q[]]];
                           {a, b, e}

You could also look at the Iterator package (available in your Maple 2017), with which you may get improved performance (of the generator, rather than of what you do with the iterates, naturally). If the cost of what you compute with the iterates dwarfs the generator cost then you might as well use whichever you find more convenient.

I have to state that your "problem" is not something to do with simplify specifically. It's to do with evaluation.

And the behavior is completely normal. Your expression that is assigned to T contains the global l1 and l2 names, and if you subsequently assign values to those same names then it should be expected that any full evaluation of T will result in the values of l1,l2 being utilized. Perhaps you know all that, ok.

Yes, you could pass an unevaluated T into a modified version (using the uneval modifier) of your procedure simplify2 and then evaluate the passed argument only 2 levels. Note that any inadvertant full evaluation of the procedure's 2-level evaluated version of the expression would break the effect.

Note that you could also pass a 1-level evaluation eval(T,1) into your original simplify2 and get a simlar result. But again, any inadvertant full evaluation within simplify2 would break the effect. simplify2_orig.mw

You could also robustify the route of using the uneval modifier, by forcibly replacing all the names with generated locals. In the attached I done that, and made a few other edits. It also shows a few other "length" metrics. cost_simplify.mw

But now I come to my main point. Your Question stated: `The variable "l1" stands for "length 1" in a mechanic problem as well as "length of the expression"'.  The core problem is use of the same global name l1 for different purposes! Don't use the name l1 (of the same scope) within an expression to stand for "length 1" (of a beam or whatever) while also using or assigning l1 (at the same scope) to mean "length of the expression". That is wrong programming. Trying strange programming hacks to work around this fundamental programming mistake is not the right way to go.

I realize that you've supplied a smaller example of your full work, to represent this issue. But much better would be to back up and figure out how to use different names for different purposes.

There are examples of plotting multiple expressions on the Help page for Topic plot.

Here is is with f and df as expressions. Notice that in this case the plot command's second argument is x=-5..5 (ie. name=range).

f := ( (x^3) / (2) ) - 5*x + 2;

               1  3          
          f := - x  - 5 x + 2
               2             

eval(f, x=4);

                 14

df := diff(f, x);

                   3  2    
             df := - x  - 5
                   2       

plot([f, df], x=-5..5);

You could also add various options, for example,

plot([f, df], x=-5..5,
     thickness=3,
     color=["Purple", "Green"],
     legend=['f','diff'('f',x)]);

See the Help page for Topic plot,options for more details.

And here it is below with f and df as procedures (arrow operators) instead. The plot will look the same. Notice that in this case the plot command's second argument is just the range -5..5 .

f := x -> ( (x^3) / (2) ) - 5*x + 2;

      f := x -> (1/2)*x^3-5*x+2

f(4);

                 14

df := D(f);

       df := x -> (3/2)*x^2-5

plot([f, df], -5 .. 5);

The very best way to make progress with Maple is to learn the programming language, bit by bit. It will reward you, going forward, and it will allow you to deal with more complicated examples than the purely "clickable"/tutor approaches could ever handle.

Make sure that you are running the latest point-release for your major version, ie. Maple 2019.2 .

The strange floats in your results hint that your license verification might be having trouble (either a cracked version or a rare false positive detection of such). The later point-release may do better.

This works in Maple 11 (as well as later versions):

expr := -(729*beta*(1/2*(-1/9*(-5/27*beta-1/9)*lambda^6-1/9*(1/9*beta^2*TT+
(10/9*TE+2/3*TT)*beta-1/9*TE)*lambda^4+5/27*(2/5*TT*(3/2*TT+TE)*beta+TE*
(TE-2/5*TT))*beta*lambda^2-1/9*TE*beta^2*TT*(-TT+TE))*p1(m,t)^2+(beta*
TT-1/3*lambda^2)^2*(-1/3*lambda^2+TE)^3*((p2(m,t))/(lambda^2-3*TE)+3/2*(
(lambda^2+TE)*p1(m,t)^2)/((lambda^2-3*TE)^3))))/((3*beta*TT-lambda^2)^3*
(3*TE-lambda^2)^2):

collect(expr, [p2,p1], simplify);

                                       2                  2         2
     3 beta p2(m, t)         (-5 lambda  + 3 beta TT) beta  p1(m, t)
   ------------------- + 3/2 ----------------------------------------
                     2                                   2 3
   3 beta TT - lambda                 (3 beta TT - lambda )

The result is very close to your given target expression, up to the leading signs or the two terms. You gave the target as,

-3*beta*p2(m,t)/(lambda^2-3*beta*TT)                                            
+3*beta^2*(5*lambda^2-3*beta*TT)*p1(m,t)^2/(2*(lambda^2-3*beta*TT)^3);

                                  2          2                      2
      3 beta p2(m, t)         beta  (5 lambda  - 3 beta TT) p1(m, t)
  - ------------------- + 3/2 ---------------------------------------
          2                                  2             3
    lambda  - 3 beta TT               (lambda  - 3 beta TT)

Calling simplify(expr) in Maple 2020 produces a result similar to this next result in Maple 11 (in which version the size option of simplify is not enabled by default).

simplify(simplify(expr),size);

     //                2\2                /          2          \               \                                                                  
     ||          lambda |                 |  5 lambda           |         2     |                                                                  
  27 ||beta TT - -------|  p2(m, t) + 1/6 |- --------- + beta TT| p1(m, t)  beta| beta                                                             
     \\             3   /                 \      3              /               /                                                                  
  ------------------------------------------------------------------------------------                                                             
                                               2 3                                                                                            
                            (3 beta TT - lambda )

Here is something from the first file: Beregningsdokument_ac.mw

There seemed to be something mixed up with the XML Font tags and the Code Edit Region. (You had two copies of the TgPakken  source within it -- I deleted the duplicate which had Font irregularities, as well as a misplaced closing </Font> tag. I'm not sure if it related to UTF-8 and strings in the Code Edit Region, or something related to an attempt at marking up a duplicate copy of the source in the C.E.R.) Please check it.

Here is something from the second: Beregningsdokument_2_ac.mw   It was the same issue, so please check this file.

This problem is not the same kind as in the older Question to which you linked, and not fixable by Joe's procedure in that thread. Also, neither Document is recovered (even partially) by opening in Maple 2020.0.

You might want to make regular backup snapshots from now on.

By the way, you seem to be running Maple 17.00, without the (free) point-release update 17.02.

First 127 128 129 130 131 132 133 Last Page 129 of 336