acer

32333 Reputation

29 Badges

19 years, 326 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are answers submitted by acer

I looked at the other responses only after getting the following last night. I am not surprised that the various approaches have elements in common.

ee := abs(z)*(z-4-I)+2*I = (5-I)*z:

ans1:=eval~(z=a+b*I,
            [solve(evalc([Re,Im](subs(z=a+I*b,ee))),explicit)]):

lprint(%);
[z = -1, z = 83/410*(-1/2*(25+2*I*26^(1/2))^(1/3)-9/2/(25+2*I*26^(1/2))^(1/3)+3-1/2*I*3^(1/2)*
((25+2*I*26^(1/2))^(1/3)-9/(25+2*I*26^(1/2))^(1/3)))^2+323/820*(25+2*I*26^(1/2))^(1/3)+2907/
820/(25+2*I*26^(1/2))^(1/3)-1063/410+323/820*I*3^(1/2)*((25+2*I*26^(1/2))^(1/3)-9/(25+2*I*26^(
1/2))^(1/3))+I*(9/410*(-1/2*(25+2*I*26^(1/2))^(1/3)-9/2/(25+2*I*26^(1/2))^(1/3)+3-1/2*I*3^(1/2
)*((25+2*I*26^(1/2))^(1/3)-9/(25+2*I*26^(1/2))^(1/3)))^2+119/820*(25+2*I*26^(1/2))^(1/3)+1071/
820/(25+2*I*26^(1/2))^(1/3)-219/410+119/820*I*3^(1/2)*((25+2*I*26^(1/2))^(1/3)-9/(25+2*I*26^(1
/2))^(1/3))), z = 83/410*((25+2*I*26^(1/2))^(1/3)+9/(25+2*I*26^(1/2))^(1/3)+3)^2-323/410*(25+2
*I*26^(1/2))^(1/3)-2907/410/(25+2*I*26^(1/2))^(1/3)-1063/410+I*(9/410*((25+2*I*26^(1/2))^(1/3)
+9/(25+2*I*26^(1/2))^(1/3)+3)^2-119/410*(25+2*I*26^(1/2))^(1/3)-1071/410/(25+2*I*26^(1/2))^(1/
3)-219/410)]

length(%);
                              1233

ans2:=simplify(evalc(ans1)):

lprint(%);
[z = -1, z = (-747/205-81/205*I)*cos(1/3*arctan(2/25*26^(1/2)))^2+1/410*(-525+195*I-(1494+162*
I)*3^(1/2)*sin(1/3*arctan(2/25*26^(1/2))))*cos(1/3*arctan(2/25*26^(1/2)))+385/82+21/82*I+(105/
82-39/82*I)*3^(1/2)*sin(1/3*arctan(2/25*26^(1/2))), z = -158/205-69/205*I+(1494/205+162/205*I)
*cos(1/3*arctan(2/25*26^(1/2)))^2+(105/41-39/41*I)*cos(1/3*arctan(2/25*26^(1/2)))]

length(%);
                               524

simplify(convert(eval~((rhs-lhs)(ee),ans1),expln));
                             [0, 0, 0]

simplify(convert(eval~((rhs-lhs)(ee),ans2),expln));
                             [0, 0, 0]

Here are a few ways to deal with this issue.

The first is to use formula for one of the domain ranges. This requires joining two surfaces, and a slight seam is visible.

The second and third involve using two variants of cylindrical coordinates. One of these has an quite visible seam. Fortunately the latter (made using an easy stock coordinate known to plot3d) is the nicest looking.

3dplot_domain.mw

Here is an image of the third plot.

Another less satisfying approach is to increase the grid size as an option to the plot3d command. This barely works since manual rotation of the plot slows down as the plot structure gets large, and the jaggedness is only reduced and smaller.

Another (advanced approach) is to write one's own MESH generator that looks for the boundary where the expression becomes non-real, and adjusts adjacent data points to match. I've been working on this (as one of many evening projects) for a while now.

[edit] I forgot that you were trying to make an Explore call. Using Carl's clever idea of sqrt(max(0,...)) then my third variant (usual cylindrical coordinates) could be explored as,

Explore( plot3d([sqrt(max(0,z^2+c))], theta=0..2*Pi, z=-2..2,
                coords=cylindrical,
                style=surface, color=green),
         c=-1.0 .. 1.0 );

The difference being that Carl has used the parametric calling-sequence of plot3d while I've used the expression calling-sequence (while we both used coords=cylindrical). Also, Carl noticed that only one of the two signed pieces is needed in this variant.

And Carl chose a nice fixed view for the plot3d call inside the Explore call. Using that nuance, the above looks even better with the extra option view=[-2..2,-2..2,-sqrt(2)..sqrt(2)] in the plot3d call.

[Wait, where did Carl's fine Answer go?! It seems to have disappeared.]

restart;

ee := RootOf(x^2+1,x):


allvalues(ee);
                             I, -I

If the expressions containing RootOf were the result from solve then you could try passing that command its explicit option.

The answer is that it is an ill-judged and misguided idea to treat an immutable data structure as if it were mutable, since that is relatively bad for performance and needlessly wastes time and memory resources.

It's also a bad idea to have an easy syntax for doing just that, especially if it would obscure the fact that a full copy of the whole structure is being done upon every "update" to an entry.

A set is not a mutable data structure. It would be misguided to treat it as if it were.

A list also is not a mutable data structure. About 20 years ago, or so, the people in charge of Maple development had a serious lapse in judgement and allowed the easy syntax of indexed assignment of a list element, eg, A[i]:=somthing where A is a list. And ever since then people have had an easier task of reducing code performance by inducing full copies of their lists upon each element "update".

 

Do you want a 2D or a 3D contour plot, and if the latter then do you want it separate from the surface plot or together?

Do you want the data interpolated?

surfcont.mw

If you want a 2D contour plot with the data interpolated (to get smoother contour lines) then you could also look at the Interpolation package new to Maple 2018.

[edited] Since you have Maple 2018 here's how you can use the new Interpolation package to obtain a single procedure that can be used to obtain the surface, or the 3D contour plot, or the 2D contour plot -- each smoothed by interpolation if you wish.

Apart from the data Matrix, all you need more are two Vectors to denote the tacit ranges for x and y that the Matrix data represents.

surfcont_Interp.mw

 

[edit] I've removed part of my answer which was relevant or correct because I failed to notice why this change below occurs. [end-of-edit]

It seems to have changed between Maple 16.02 and Maple 17.02. That is, using the CLI,

/usr/local/maple/maple17.02/bin/maple -q -s -c "print(kernelopts(version));" -c "print(evalf(frac(Pi^20)));" -c "quit"

                               Maple 17.02, X86 64 LINUX, Sep 5 2013, Build ID 872941

                                                         23.

/usr/local/maple/maple16.02/bin/maple -q -s -c "print(kernelopts(version));" -c "print(evalf(frac(Pi^20)));" -c "quit"

                               Maple 16.01, X86 64 LINUX, May 6 2012, Build ID 744592

                                                         0.

 

The following computes quickly in my Maple 2018.0 version.

restart;

ee := sqrt(b)*sqrt(1-4*p/b)-2*arctan(sqrt((9*p/b-22201/10000)/(9/4-9*p/b))):

plots:-implicitplot([Re(ee)=0, Im(ee)=0], b=0..5, p=0..5,
                    gridrefine=1, crossingrefine=7);

Your expression is 2*Pi periodic, so we can call minimize or maximize with the restricted range x=0..2*Pi.

Or we can do the general computation, using diff, solve, etc.

restart;

f := (cos(x)+sqrt(3)*sin(x))/(cos(x)+sin(x)+2);

(cos(x)+3^(1/2)*sin(x))/(cos(x)+sin(x)+2)

(1)

simplify([maximize(f, x=0..2*Pi,'location')]);
evalf(%);

[4*(12+2*3^(1/2))^(1/2)/(3^(1/2)*(12+2*3^(1/2))^(1/2)+(12+2*3^(1/2))^(1/2)+2*3^(1/2)+12), {[{x = arctan((3^(1/2)*(12+2*3^(1/2))^(1/2)-1+3^(1/2))/(3^(1/2)-3+(12+2*3^(1/2))^(1/2)))}, 4*(12+2*3^(1/2))^(1/2)/(3^(1/2)*(12+2*3^(1/2))^(1/2)+(12+2*3^(1/2))^(1/2)+2*3^(1/2)+12)]}]

 

[.6001956932, {[{x = 1.231247593}, .6001956932]}]

(2)

simplify([minimize(f, x=0..2*Pi,'location')]);
evalf(%);

[4*(12+2*3^(1/2))^(1/2)/(3^(1/2)*(12+2*3^(1/2))^(1/2)+(12+2*3^(1/2))^(1/2)-2*3^(1/2)-12), {[{x = arctan((3^(1/2)*(12+2*3^(1/2))^(1/2)+1-3^(1/2))/(-3^(1/2)+3+(12+2*3^(1/2))^(1/2)))+Pi}, 4*(12+2*3^(1/2))^(1/2)/(3^(1/2)*(12+2*3^(1/2))^(1/2)+(12+2*3^(1/2))^(1/2)-2*3^(1/2)-12)]}]

 

[-3.332246505, {[{x = 4.004740164}, -3.332246505]}]

(3)

sols := map(s->simplify(radnormal(s)),
          [solve(diff(f,x), x, allsolutions)]);

[arctan(2+(1/6)*(3-3^(1/2))*(12+2*3^(1/2))^(1/2))+2*Pi*_Z3, arctan(2+(1/6)*(3^(1/2)-3)*(12+2*3^(1/2))^(1/2))+(2*_Z3-1)*Pi]

(4)

maxf := radnormal(expand( eval(f, x=sols[1]) ));
evalf(%);

(3^(1/2)*(12+2*3^(1/2))^(1/2)-(12+2*3^(1/2))^(1/2)+4*3^(1/2)+2)/(3^(1/2)*(12+2*3^(1/2))^(1/2)+4*3^(1/2)+(12+2*3^(1/2))^(1/2)+2)

 

.6001956940

(5)

minf := radnormal(expand( eval(f, x=sols[2]) ));
evalf(%);

(3^(1/2)*(12+2*3^(1/2))^(1/2)-(12+2*3^(1/2))^(1/2)-4*3^(1/2)-2)/(3^(1/2)*(12+2*3^(1/2))^(1/2)-4*3^(1/2)+(12+2*3^(1/2))^(1/2)-2)

 

-3.332246501

(6)

plot(f, x=-10..10, view=-4..1);

 

 


Download max_min_example.mw

 

Here is an example, using Q__outflow and e as parameters (whose earlier assignments I commented out).

I changed Q[outflow] to Q__outflow , although in 2D Input they both appear as the subscripted name.

You could comment out additonal parameters and add them into the Explore call.

382483three_reservoir_prob_ac.mw

 

Try using calls to time[real]() before and after, instead of time() .

(Perhaps also consider CodeTools:-Usage )

Here is what there was, up until the point of corruption (in a Plot Output).

This attached is essentially the same as what the Standard GUI of either Maple 2016.2 or Maple 2018.0 are able to recover, plus one additonal trailing paragraph of 2D Input that I could salvage.

But the next Output onwards, to the end of the file, is corrupt and lacking information.

I've left it as 2D Input and a Document.

PLOTS-SOLID_ac.mw

ps. Why do you not install the 2016.2 point-release , which is a free download if you already have a licensed copy of Maple 2016.0? It may not have prevented this corruption example, but perhaps it's a good idea regardless.

Why are you not telling us specifically what you intend on doing with the intermediate expression (with the additional abs calls in it)?

It is not helpful to hold back such additonal information, which can make someone's suggested approach work -- or not -- in your broader endeavor.

Apart from not wasting responders' time, there's also a good chance that they might suggest a better way to approach the larger problem.

In the absence of those extra details about your intermediary computations, I'll throw in the following suggestion below.

Note that you can make use of a terse calling sequence to the op command, in lieu of the nested call(s).

restart;

Y:=sqrt(1-y^2)/x + x/sqrt(1-y^2);

(-y^2+1)^(1/2)/x+x/(-y^2+1)^(1/2)

(1)

Y:=subsindets(Y,anything^({1/2,-1/2}),ee->%abs(op(1,ee))^(op(2,ee)));     

%abs(-y^2+1)^(1/2)/x+x/%abs(-y^2+1)^(1/2)

(2)

subsindets(Y,%abs(anything)^({1/2,-1/2}),ee->op(op(1,ee))^op(2,ee));

(-y^2+1)^(1/2)/x+x/(-y^2+1)^(1/2)

(3)

subsindets(Y,%abs(anything)^({1/2,-1/2}),ee->op([1,..],ee)^op(2,ee));

(-y^2+1)^(1/2)/x+x/(-y^2+1)^(1/2)

(4)

 

Download subsindets_abs_sqrt.mw

More to the point: why can't you use some other name entirely, instead of abs?

Done in Maple 2018.0 on 64bit Linux on a fast machine.

restart;

plots:-setoptions(gridlines=false);

Omega := 2*Pi*N;
R0 := a*tanh((a^2-mu)/(2*T_c))*ln((2*a^2+2*a*q+q^2-2*mu-I*Omega)/(2*a^2-2*a*q+q^2-2*mu-I*Omega))/q-2:                     
T_c := 0.169064e-1; mu := .869262; #N := 0;
R1re := Q->Int(subs(q=Q,unapply(Re(R0),a)), 0.1e-3 .. 100, method=_d01ajc):
R1im := Q->Int(subs(q=Q,unapply(Im(R0),a)), 0.1e-3 .. 100, method=_d01ajc):
R1 := unapply('R1re'(q)+'I*R1im(q)', q, numeric):

2*Pi*N

 

0.169064e-1

 

.869262

(1)

N := 0;

st := time():

  plot(abs(R1(q)), q = 0.1e-2 .. 10);

(time()-st)*'seconds';

0

 

 

.822*seconds

(2)

N := 10;

st := time():

  plot(abs(R1(q)), q = 0.1e-2 .. 10);

(time()-st)*'seconds';

10

 

 

.740*seconds

(3)

st := time():

count := 0;
for NN from 0 to 10 by 0.25 do
  count := count+1;
  N := NN;
  PP[count]:=plot(abs(R1(q)), q = 0.1e-2 .. 10);
end do:

count;

(time()-st)*'seconds';

0

 

41

 

31.282*seconds

(4)

plots:-display(convert(PP,list), insequence=true);

 

 

Download plot_abs_ReIm.mw

 

restart;

expr := 7*ln(arcsin(x))-(1/2)*ln(x-1)*sin(x)-(1/2)*ln(x+1)+f;

7*ln(arcsin(x))-(1/2)*ln(x-1)*sin(x)-(1/2)*ln(x+1)+f

(1)

subsindets(expr, 'specfunc(ln)', ee->ln(abs(op(ee))));

7*ln(abs(arcsin(x)))-(1/2)*ln(abs(x-1))*sin(x)-(1/2)*ln(abs(x+1))+f

(2)

expr:=sqrt(x)+1/sqrt(x);

x^(1/2)+1/x^(1/2)

(3)

subsindets(expr, anything^({1/2,-1/2}), ee->abs(op(1,ee))^(op(2,ee)));

abs(x)^(1/2)+1/abs(x)^(1/2)

(4)

 

Download subsindets_type.mw

 

orthopoly is a table-based package, not a module.

You could do it this way:

f := proc(n)
   uses H=orthopoly[':-H'];
   H(n,z);
end proc;

If your source for f is in a plaintext file (being read of loaded to a .mla Library) then you could also utilize the $define directive instead of uses. Eg.

$define H=orthopoly[':-H']

f := proc(n)
   H(n,z);
end proc;

$undef H

In this way the source code for f is still more legible.

[edited: a side-note for some explanation of the original issue] When the table-based package orthopoly is loaded using with the global name :-H is assigned the procedure orthopoly[H]. That's in contrast to module-based packages, where invoking with will rebind the name of an export but not assign to the global name. Eg,

with(orthopoly):
eval(:-H);
                     proc(n, x)  ...  end;

with(LinearAlgebra):
eval(:-MatrixNorm);
                           MatrixNorm

A table-based package is not a module, and doesn't have exports per se.

But the uses package statement sets up a temporary scope for the binding of the names of package exports. As such it works with module-based packages (in examples with the terse syntax like uses LinearAlgebra, where the names of all exports get rebound within scope), whose functionality relates to name binding. But it doesn't work with table-based packages (using a syntax like your uses orthopoly attempt), because Maple doesn't see exports to rebind.

exports(orthopoly);
Error, invalid input: exports expects its 1st argument, m,
to be of type `module`, but received orthopoly

It's an understatement to say that the Help page documentation of uses is thin.

First 176 177 178 179 180 181 182 Last Page 178 of 336