acer

32313 Reputation

29 Badges

19 years, 315 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are replies submitted by acer

I have deleted a duplicate of this Question.

You can add followup details or followup queries here, for this problem. Please don't submit new and separate Question threads for it.

I advise you to provide code that reproduces the problem. You could upload and attach a worksheet using the green up-arrow in the Mapleprimes editor.

What does it mean to you, for a plot to have a transparent background?

Do you mean it in the context of exporting, eg. as .GIF file?

What would the contour represent?

(Once you answer that first query) How many such contours do you think there are?

@AmirHosein Sadeghimanesh 

This user repeatedly asks Questions here without mentioning that he is using the older version Maple 2019. (I have now marked this Question with that version number.)

In Maple 2019 the Student:-Basics package did not contain the FactorSteps, SolveSteps, or SimplifiySteps commands.

Using a single left-quoted name as first argument to nprintf is technically weaker because it can be broken by the obscure possibility of prior assignment to that name.  As unlikely as that is possibility is, it is more defensive programming to prevent that breakage.

A string is the simplest, safer alternative. The alternative of a left single-quoted reference to the global right single-quoted is more complicated. In your example with curry that defensive alternative may require nested single left-quotes for such protective safety, at which point the process is far more complicated than simply using double-quotes. Eg,

restart;
`#mo("%s",mathcolor="%a");`:=oof:

cprint:= curry(nprintf, '':-`#mo("%s",mathcolor="%a");`''):

cprint("This is green", green);
cprint("This is green", red);

note: I know that Carl is well aware of all this. I mention it for the general readership.

@nm You might consider using the "Branch" icon at the bottom of your Question body (next to Reply/Answer.).

That should create a separate Question thread, but with cross-reference link added automatically.

The OP's attached worksheet executes without error for me in both Maple 2022.0 and Maple 2022.1 for Linux.

@ijuptilk You can augment that list, conts, with any value you choose. Below I include the min and max values from the generated Array of data.

In the case that the surface lies very close to its minimum (or maximum) over an area (not just points, or along curves) then there may be visual artefacts when the contour plot is rendered. I accomodate that by making a small offset to the contour values, which seems to satisfactorily make the artefact vanish.

As for an Array of the data, you already asked me about that. It was op([1,3],...) of the structure returned by plot3d. You can assign those to names of your choice. I do this below.

Recall that I mentioned before that you could construct an Array (or Matrix, or list of lists) of data either using a nested seq call, or using the plot3d command. I used the latter, as I think that makes the associated ranges easier for you to understand than for the corresponding increments for the seq calls..

By the way, given the Array of data, a you could also call listcontplot and obtain a similar result as contourplot (with the corresponding grid and ranges). You can transform the result to the corresponding ranges. See attached.

restart:

doCalc:= proc(xi, u)  #u is the \bar(H): normalize magnetic field magnitude,
                          # where H = bar(H)*H__a
                 option remember;
                 # Import Packages
                 uses ArrayTools, Student:-Calculus1, LinearAlgebra,
                      ListTools, RootFinding, plots, ListTools:
                 local gamma__1:= .1093,
                       alpha__3:= -0.1104e-2,
                       k__1:= 6*10^(-12),
                       d:= 0.2e-3,
                       theta0:= 0.1e-3,
                       eta__1:= 0.240e-1, chi:= 1.219*10^(-6),
                       alpha:= 1-alpha__3^2/(gamma__1*eta__1),
                       theta_init:= theta0*sin(Pi*z/d),
                       H__a:= Pi*sqrt(k__1/chi)/d,
                       H := u*H__a,     
                       c:=alpha__3*xi*alpha/(eta__1*(4*k__1*q^2/d^2-alpha__3*xi/eta__1 - chi*H^2)),
                       w := chi*H^2*eta__1*alpha/(4*k__1*q^2/d^2-alpha__3*xi/eta__1 - chi*H^2),
                       n:= 20,
                       g, f, b1, b2, qstar, OddAsymptotes, ModifiedOddAsym,
                       qstarTemporary, indexOfqstar2, qstar2, AreThereComplexRoots,
                       soln1, soln2, qcomplex1, qcomplex2, gg, qq, m, pp, j, i,
                       AllAsymptotes, p, Efun, b, aa, F, A, B, Ainv, r, theta_sol, v, Vfun, v_sol,minp,nstar;

                 if not [xi,u]::[numeric,numeric] then
                    return 'procname'(args);
                 end if;

# Assign g for q and plot g, Set q as a complex and Compute the Special Asymptotes

                 g:= q-(1-alpha)*tan(q)- (w*q + c)*tan(q):
                 f:= subs(q = x+I*y, g):
                 b1:= evalc(Re(f)) = 0:
                 b2:= evalc(Im(f)) = 0:
                 qstar:= (fsolve(1/c = 0, q = 0 .. infinity)):
                 OddAsymptotes:= Vector[row]([seq(evalhf(1/2*(2*j + 1)*Pi), j = 0 .. n)]);

# Compute Odd asymptote

                 ModifiedOddAsym:= abs(`-`~(OddAsymptotes, qstar));
                 qstarTemporary:= min(ModifiedOddAsym);
                 indexOfqstar2:= SearchAll(qstarTemporary, ModifiedOddAsym);
                 qstar2:= OddAsymptotes(indexOfqstar2);

# Compute Odd asymptote

                 AreThereComplexRoots:= type(true, 'truefalse');
                 try
                      soln1:= fsolve({b1, b2}, {x = min(qstar2, qstar) .. max(qstar2, qstar), y = 0 .. infinity});
                      soln2:= fsolve({b1, b2}, {x = min(qstar2, qstar) .. max(qstar2, qstar), y = -infinity .. 0});
                      qcomplex1:= subs(soln1, x+I*y);
                      qcomplex2:= subs(soln2, x+I*y);
                 catch:
                       AreThereComplexRoots:= type(FAIL, 'truefalse');
                 end try;

# Compute the rest of the Roots

                 #OddAsymptotes:= Vector[row]([seq(evalhf((1/2)*(2*j+1)*Pi), j = 0 .. n)]);
                 AllAsymptotes:= sort(Vector[row]([OddAsymptotes, qstar]));
                 if   AreThereComplexRoots
                 then
                      gg := [ qcomplex1,
                              qcomplex2,
                              op(Roots(g, q = 0.1e-3 .. AllAsymptotes[n], numeric))
                            ];
                 elif not AreThereComplexRoots
                 then gg:= [ op(Roots(g, q = 0.1e-3 .. AllAsymptotes[n], numeric))
                           ];
                 end if:

# Remove the repeated roots if any & Redefine n

                 qq:= MakeUnique(gg):
                 m:= numelems(qq):

## Compute the `τ_n`time constants

                 for i to m do
                 p[i] := evalf(gamma__1*alpha/(4*k__1*qq[i]^2/d^2 - alpha__3*xi/eta__1- chi*H^2)):
if not p[i]::complex(numeric) then print(p[i], [xi,u], qq[i]); end if;
                 end do:

                Digits := 15;
## Compute θ_n from initial conditions

                for i to m do
                Efun[i] := cos(qq[i]-2*qq[i]*z/d)-cos(qq[i]):
                end do:

## Compute integral coefficients for off-diagonal elements θ_n matrix

                printlevel := 2:
                for i to m do
                    for j from i+1 to m do
                        b[i, j] := evalf(Int(Efun[i]*Efun[j], z = 0 .. d)):
                        if not b[i, j]::complex(numeric) then
                            b[i, j] := evalf[15](Int(Efun[i]*Efun[j], z = 0 .. d,
                                                 digits=15, method=_Dexp, epsilon=1e-12)):
                            if not b[i, j]::complex(numeric) then
                               print("A",[Efun[i]*Efun[j], z = 0 .. d]);
                            end if;
                        end if;
                        b[j, i] := b[i, j]:
                        aa[i, j] := b[i, j]:
                        aa[j, i] := b[j, i]:
                    end do :
                end do:

## Calculate integral coefficients for diagonal elements in theta_n matrix

                for i to m do
                   aa[i, i] := evalf(Int(Efun[i]^2, z = 0 .. d)):
                   if not aa[i, i]::complex(numeric) then
                      aa[i, i] := evalf[15](Int(Efun[i]^2, z = 0 .. d,
                                                digits=15, epsilon=1e-12));
                      if not aa[i, i]::complex(numeric) then
                         print("B",[Efun[i]^2, z = 0 .. d]);
                      end if;
                   end if;
                end do:

## Calculate integrals for RHS vector

               for i to m do
               F[i] := evalf(Int(theta_init*Efun[i], z = 0 .. d));
               if not F[i]::complex(numeric) then
                  F[i] := evalf[15](Int(theta_init*Efun[i], z = 0 .. d,
                                     digits=15, method=_Dexp, epsilon=1e-12));
                  if not F[i]::complex(numeric) then
                     print("C",[theta_init*Efun[i], z = 0 .. d]);
                  end if;
               end if;
               end do:

## Create matrix A and RHS vector B

               A := Matrix([seq([seq(aa[i, j], i = 1 .. m)], j = 1 .. m)]):
               B := Vector([seq(F[i], i = 1 .. m)]):

## Calculate solve A*x=B

              r := LinearSolve(A,B);

## Define Theta(z,t)
             theta_sol := add(r[i]*Efun[i]*exp(-t/p[i]), i = 1 .. m):

## Compute v_n for n times constant

             for i to m do
             v[i] := (-2*k__1*alpha__3*qq[i])*(1/(d^2*eta__1*alpha*gamma__1))+ alpha__3^2*xi/(2*(eta__1)^2*qq[i]*alpha*gamma__1)+xi/(2*eta__1*qq[i]) + alpha__3*chi*H^2/(2*eta__1*qq[i]*gamma__1*alpha):
             end do:

## Compute v(z,t) from initial conditions
             for i to m do
             Vfun[i] := d*sin(qq[i]-2*qq[i]*z/d)+(2*z-d)*sin(qq[i]):
             end do:

## Define v(z,t)
             v_sol := add(v[i]*r[i]*Vfun[i]*exp(-t/p[i]), i = 1 .. m):

##
             minp:=min( seq( Re(p[i]), i=1..m) ):
             member(min(seq( Re(p[i]), i=1..m)),[seq( Re(p[i]), i=1..m)],'nstar'):

## Return all the plots
                 return theta_init, theta_sol, v_sol, minp, eval(p), nstar, p[nstar], g, H, H__a;
                 end proc:

# Run the calculation for supplied value of 'xi'

# Put the returned quantities in a simple list

ans:=CodeTools:-Usage([doCalc(-0.06, 0.001)]):
evalf(ans[9]);
evalf(ans[10]);

memory used=454.19MiB, alloc change=366.15MiB, cpu time=4.71s, real time=4.53s, gc time=486.75ms

0.3484926715e-1

34.84926715

with(plots):

d:= 0.2e-3:

testproc := proc(j, u, k) option remember;
   local calcvals;
   if not [j,u,k]::[numeric,numeric,posint] then
      return 'procname'(args);
   end if;
   calcvals:=doCalc(j,u);
   evalf(calcvals[k]);
end proc:

NN:=8; # number of values (excluding minv and maxv)

8

(gridji,rngj,rngi):=[11,21],1..3,-2.0..-0.0;
PP[gridji,4]:=CodeTools:- Usage(plot3d(max(-1e2,testproc(i,j,4)), j=rngj, i=rngi, grid=gridji)):

[11, 21], 1 .. 3, -2.0 .. -0.

memory used=105.53GiB, alloc change=160.00MiB, cpu time=18.98m, real time=16.00m, gc time=2.84m

PP[gridji,4];

data[gridji,4]:=op([1,3],PP[gridji,4]);

data[[11, 21], 4] := Vector(4, {(1) = ` 1..11 x 1..21 `*Array, (2) = `Data Type: `*float[8], (3) = `Storage: `*rectangular, (4) = `Order: `*C_order})

(minv,maxv):=[min,max](data[gridji,4])[];

-100., -2.78860304200000009

(color1,color2):="Red","Yellow":
conts:=[minv, seq(minv+(i-1)*(maxv-minv)/(NN+2-1),i=2..NN+1), maxv];

[-100., -89.19873367, -78.39746734, -67.59620101, -56.79493469, -45.99366836, -35.19240202, -24.39113570, -13.58986937, -2.78860304200000009]

colorlist:=ColorTools:-Color~(((c1,c2,N)->[seq(c1+(i-1)*(c2-c1)/(N-1),
                                               i=1..N)])([ColorTools:-Color(color2)[]],
                                                         [ColorTools:-Color(color1)[]],
                                                         nops(conts))):
#re-using values computed during plot3d on same grid.
CodeTools:- Usage(contourplot(max(-1e2,testproc(i,j,4)), j=rngj, i=rngi, grid=gridji,
                              contours=conts, thickness=0, coloring=[color1,color2],
                              filledregions, axes=boxed)):
PC:=subsindets(%,specfunc(anything,THICKNESS),u->THICKNESS(0.2)):
display(PC,
        seq(plot(-2.0,1..1,legend=evalf[4](conts[i]),thickness=15,
                 color=colorlist[nops(conts)-i+1],legendstyle=[location=right]),
            i=1..nops(conts)),
        size=[500,400]);

memory used=1.76MiB, alloc change=0 bytes, cpu time=59.00ms, real time=60.00ms, gc time=0ns

subsindets(plots:-listcontplot(Matrix(data[gridji,4]),
                               contours=conts-~1e-9,filledregions,thickness=0),
           specfunc(anything,THICKNESS),u->THICKNESS(0.2)):
LCP[gridji,4]:=plottools:-transform((x,y)->[lhs(rngj)+(rhs(rngj)-lhs(rngj))*(x-1)/(gridji[1]-1),
                                            lhs(rngi)+(rhs(rngi)-lhs(rngi))*(y-1)/(gridji[2]-1)])(%):
display(LCP[gridji,4],
        seq(plot(-2.0,1..1,legend=evalf[4](conts[i]),thickness=15,
                 color=colorlist[nops(conts)-i+1]),
            i=1..nops(conts)),
        size=[500,400],legendstyle=[location=right], axes=boxed);

(gridji,rngj,rngi):=[11,21],1..3,-2.0..-0.0;
PP[gridji,7]:=CodeTools:- Usage(plot3d(Im(testproc(i,j,7)), j=rngj, i=rngi, grid=gridji)):

[11, 21], 1 .. 3, -2.0 .. -0.

memory used=301.40KiB, alloc change=0 bytes, cpu time=3.00ms, real time=4.00ms, gc time=0ns

PP[gridji,7];

data[gridji,7]:=op([1,3],PP[gridji,7]);

data[[11, 21], 7] := Vector(4, {(1) = ` 1..11 x 1..21 `*Array, (2) = `Data Type: `*float[8], (3) = `Storage: `*rectangular, (4) = `Order: `*C_order})

(minv,maxv):=[min,max](data[gridji,7])[];

-36.4108769999999993, 0.

(color1,color2):="Red","Yellow":
conts:=[minv, seq(minv+(i-1)*(maxv-minv)/(NN+2-1),i=2..NN+1), maxv];

[-36.4108769999999993, -32.36522400, -28.31957100, -24.27391800, -20.22826500, -16.18261200, -12.13695900, -8.09130600, -4.04565300, 0.]

colorlist:=ColorTools:-Color~(((c1,c2,N)->[seq(c1+(i-1)*(c2-c1)/(N-1),
                                               i=1..N)])([ColorTools:-Color(color2)[]],
                                                         [ColorTools:-Color(color1)[]],
                                                         nops(conts))):
#re-using values computed during plot3d on same grid.
CodeTools:- Usage(contourplot(Im(testproc(i,j,7)), j=rngj, i=rngi, grid=gridji,
                              contours=conts-~1e-10, thickness=0, coloring=[color1,color2],
                              filledregions, axes=boxed)):
PC:=subsindets(%,specfunc(anything,THICKNESS),u->THICKNESS(0.2)):
display(PC,
        seq(plot(-2.0,1..1,legend=evalf[4](conts[i]),thickness=15,
                 color=colorlist[nops(conts)-i+1],legendstyle=[location=right]),
            i=1..nops(conts)),
        size=[500,400]);

memory used=1.49MiB, alloc change=0 bytes, cpu time=49.00ms, real time=51.00ms, gc time=0ns

subsindets(plots:-listcontplot(Matrix(data[gridji,7]),
                               contours=conts-~1e-9,filledregions,thickness=0),
           specfunc(anything,THICKNESS),u->THICKNESS(0.2)):
LCP[gridji,7]:=plottools:-transform((x,y)->[lhs(rngj)+(rhs(rngj)-lhs(rngj))*(x-1)/(gridji[1]-1),
                                            lhs(rngi)+(rhs(rngi)-lhs(rngi))*(y-1)/(gridji[2]-1)])(%):
display(LCP[gridji,7],
        seq(plot(-2.0,1..1,legend=evalf[4](conts[i]),thickness=15,
                 color=colorlist[nops(conts)-i+1]),
            i=1..nops(conts)),
        size=[500,400],legendstyle=[location=right], axes=boxed);

 

NULL

Download Case_3_IjuptilK_300522_ac_data.mw

@Rouben Rostamian  A friendly comment: one doesn't need the expense of unapply, merely to perform this substitution.

T := expr->eval(expr,[x=(sqrt(3)*y-x)/2,
                      y=(-sqrt(3)*y-x)/2]):

I have deleted seven duplicates of this Question, all posted within the last 45 minutes. Please stop spamming the forum with these.

If you have followup queries or details for this problem then please add them as Comments/Replies in this Question thread, instead of spawing new and separate Question threads.

@Tokoro Intersection points at the boundaries should also be included in the candidate points (ie. as well as what extrema returns). And the results from extrema can be filtered, according to the desired ranges.

I left out the end-points for the previous computation of the maximum because I could see by the graph that it did not occur at a boundary. But ideally you could use the a robust approach for both min and max.

restart;

m := 1: n := 2:

Zc := m*cos(1/180*x*Pi) + n*cos(1/180*y*Pi);

cos((1/180)*x*Pi)+2*cos((1/180)*y*Pi)

Zs := m*sin(1/180*x*Pi) + n*sin(1/180*y*Pi);

sin((1/180)*x*Pi)+2*sin((1/180)*y*Pi)



From the graph we can observe that the curve passes through the boundary
only along x=0 or x=180.
Programmatically we could also solve along y=0 and y=180 (and get no further
results there).

endpts := [map(yy->[x=0,y=yy],
               select(y->is(y>=0 and y<=180),
                      [solve(eval(Zc=Zs,x=0))]))[],
           map(yy->[x=180,y=yy],
               select(y->is(y>=0 and y<=180),
                      [solve(eval(Zc=Zs,x=180))]))[]]:
endpts := simplify(endpts);

[[x = 0, y = (45*arctan(3*7^(1/2))+45*Pi)/Pi], [x = 180, y = (-45*arctan(3*7^(1/2))+45*Pi)/Pi]]

extrema(Zs, {Zc=Zs}, {x,y}, 's'):


Filter from `s` the points which are in the desired ranges.
Augment that with the end-points.

news := select(sol->eval(x,sol)>=0 and eval(x,sol)<=180
               and eval(y,sol)>=0 and eval(y,sol)<=180,
               s) union {endpts[]};

{[x = 0, y = (45*arctan(3*7^(1/2))+45*Pi)/Pi], [x = 180, y = (-45*arctan(3*7^(1/2))+45*Pi)/Pi], {x = 45, y = 45}}

objs := eval~(Zs, news);

{(3/2)*2^(1/2), 2*cos((1/4)*arctan(3*7^(1/2))+(1/4)*Pi), 2*sin((1/4)*arctan(3*7^(1/2))+(1/4)*Pi)}

minobj := min(objs);

2*cos((1/4)*arctan(3*7^(1/2))+(1/4)*Pi)

minpts := select(u->eval(Zs,u)=minobj, news);

{[x = 180, y = (-45*arctan(3*7^(1/2))+45*Pi)/Pi]}

evalf(minobj), evalf(minpts);

.8228756542, {[x = 180., y = 24.29518893]}

maxobj := max(objs);

(3/2)*2^(1/2)

maxpts := select(u->eval(Zs,u)=maxobj, news);

{{x = 45, y = 45}}

evalf(maxobj), evalf(maxpts);

2.121320343, {{x = 45., y = 45.}}

 

Download plot-15-iPlot_ac.mw

@pik1432 Only with the interpretation ( eg. of p*(Delta*delta) ) as multiplication is your input equivalent to your stated desired result.

With the function call interpretation your input expression is not equivalent to your desired target expression.

@ijuptilk The contourplot command had not been updated in Maple 18 to make use of a thickness value less than 1.

And (funny as it sounds) thickness=0 produces a line actually thicker than it could be.

But the Maple GUI knows how to render even thinner lines. The trick is to forcibly shoehorn the smaller thickness value into the constructed PLOT structure.

The contourplot command doesn't make use of, say, thickess=0.2 as an option. But we can forcibly replace with what we want in the PLOT structure, using subsindets.

@ijuptilk Here is some explanation of those code snippets that you mentioned:

P:=plot3d(x*y,x=1..2,y=1..2,grid=[3,4]):


This is the z-data, ie. height, or values of x*y.

op([1,3],P);

Matrix(3, 4, {(1, 1) = 1., (1, 2) = 1.33333333333333, (1, 3) = 1.66666666666667, (1, 4) = 2.00000000000000, (2, 1) = 1.50000000000000, (2, 2) = 2., (2, 3) = 2.50000000000000, (2, 4) = 3.00000000000000, (3, 1) = 2., (3, 2) = 2.66666666666667, (3, 3) = 3.33333333333333, (3, 4) = 4.00000000000000})


The min and max values from those z-values, then assigned to minv and maxv respectively..

(minv,maxv) := [min,max](op([1,3],P))[];

1., 3.99999999999999956

nconts:=4:


This is a variant of a standard technique, and generates a sequence of nconts evenly spaced values. But I chose to omit the lower value, so that in the contour-plot the colored/filled regions denote areas less than the next higher contour. (In that manner, there's no point in showing the lowest contour corresponding to z=x*y=minv .)


conts:=[seq(minv+(u-1)*(maxv-minv)/(nconts+2-1),u=1..nconts+2)][2..nconts+1];

[1.600000000, 2.200000000, 2.800000000, 3.400000000]

(color1,color2):="Red","Yellow":


This turns each of those two colors into the corresponding two lists of RGB values.

temp := [ColorTools:-Color(color2)[]], [ColorTools:-Color(color1)[]];

[1.00000000, 1.00000000, 0.], [1.00000000, 0., 0.]


This generates a sequence of lists that transition linearly between the first and the second of the above RGB representations of color1 and color2. This is that standard technique, once again: generate N results on the (here, 3-dimensional) line between c1 and c2 inclusive.

((c1,c2,N)->[seq(c1+(u-1)*(c2-c1)/(N-1),u=1..N)])(temp, nops(conts));

[[1.00000000, 1.00000000, 0.], [1.00000000, .6666666667, 0.], [1.00000000, .3333333333, 0.], [1.00000000, 0., 0.]]

Download some_expl.mw

First 103 104 105 106 107 108 109 Last Page 105 of 591