acer

with assuming...

Answered: acer 32373

September 10 2022

1 0

alpha^m*alpha^n*n/alpha

Download Q_9-9-22_limit_in_integral_ac.mw

type suffixed...

Answered: acer 32373

September 09 2022

1 1

The names generated by solve come in only a few flavours, so you might be able to get by without need for consideration of the names in the original expression f. For example,

restart;
assume(M::nonnegative):
f := 10*cos((-1+t)/sqrt(1+M))-10*cos(t/sqrt(1+M)):
g := solve({diff(f, t), t>0 }, t, allsolutions)[1][1]:

indets(rhs(g),And(name,Not(constant),
                  ':-suffixed(:-_)'));

           {_Z2}

indets(rhs(g),And(name,Not(constant),
                  ':-suffixed({:-_B,:-N,:-_Z})'));

           {_Z2}

I don't see a natural reason to attempt a heavy handed approach like as follows. (The original expression might have already had such generated names, in which case a combined approach might be needed.)

convert(indets(f,And(name,Not(constant))),`global`):
convert(indets(rhs(g),And(name,Not(constant))),`global`) minus %;

           {_Z2}

Maple 2015.2 vs Maple 2022.1...

Answered: acer 32373

September 08 2022

0 1

Simpler results can be had in Maple 2022.1, using the method below.

Perhaps someone can come up with something more direct. I programmatically determine the second root, to be less ad hoc about passing a range to Calculus1:-Roots.

The Calculus1:-Roots and maximize commands have a special (peppering) technique for determining values of the _Zxx assumed-integer parameter that solve returns here when passed its allsolutions option.

The solve command has improved in the past few releases, as far as returning a finite number of explicit roots when passed inequalities that denote a finite range. But it still seems not generally as strong as that internal technique of Roots and maximize. For this example solve returns a general formula when supplied its allsolutions option, but seems to return NULL when instead passed various reasonable inequalities, e.g t>0,t<10.

Here it is first with Maple 2015.2,

>

restart;

>	kernelopts(version);

>	with(plots):

>	shock := piecewise(t <0, 0, t < 1, 10, 0): sys := {(M__p+M__a)diff(x(t), t$2)=M__pshock-x(t), x(0)=0, D(x)(0)=0}

sys := {(`#msub(mi("M"),mi("p"))`+`#msub(mi("M"),mi("a"))`)*(diff(x(t), t, t)) = `#msub(mi("M"),mi("p"))`*piecewise(t < 0, 0, t < 1, 10, 0)-x(t), x(0) = 0, (D(x))(0) = 0}

>	sol := unapply(rhs(dsolve(sys)), (M__p,M__a))

proc (M__p, M__a) options operator, arrow; piecewise(t < 0, 0, t < 1, -10*M__p*cos(t/(M__p+M__a)^(1/2))+10*M__p, 1 <= t, 10*M__p*cos((-1+t)/(M__p+M__a)^(1/2))-10*M__p*cos(t/(M__p+M__a)^(1/2))) end proc

>	# Computation of tend for M__p=M__a=2 # # Im expecting to get a value around 6.8

>	gensol := solve({sol(2,2),t>0},t,allsolutions);

${t = -2*arctan((cos(1/2)-1)/sin(1/2))+2*Pi*_Z2}$

>	RootFinding:-NextZero(unapply(sol(2,2),t),0); secondroot := RootFinding:-NextZero(unapply(sol(2,2),t),%);

>	# Maple 2022.1 gives 2*Pi+1/2 directly here, without using # 'expand' with 'Roots', and without even any simplification # done to the result # cands := {Student:-Calculus1:-Roots(expand(sol(2,2)),t=0..secondroot)[]}; tend := [op(cands minus {0})][1]; evalf(tend);

{0, -2*arctan((cos(1/2)-1)/sin(1/2))+2*Pi}

-2*arctan((cos(1/2)-1)/sin(1/2))+2*Pi

>	# Maple 2022.1 gives 0 here, without calling 'expand'. # xtend := simplify(expand(eval(sol(2,2),t=tend)));

>	# Maple 2022.1 gives 40sin(1/4) directly here, without even calling 'simplify' # xmax := maximize(20cos(((t - 1)sqrt(4))/4) - 20cos(sqrt(4)*t/4)); xmax := simplify(xmax); evalf(xmax);

40/(1+cos(1/2)^2/sin(1/2)^2+2*cos(1/2)/sin(1/2)^2+1/sin(1/2)^2)^(1/2)

>	# Maple 2022.1 gives Pi+1/2 here, without even using 'simplify'. # tmax := Student:-Calculus1:-Roots(expand(sol(2,2))=xmax,t=0..10)[1]; tmax := simplify(tmax); evalf(tmax);

-2*arctan((cos(1/2)+1)/sin(1/2))+2*Pi

>	display( pointplot([[tmax,xmax], [tend,xtend]], symbolsize=20, symbol=solidcircle), plot(sol(2,2),t=0..10, size=[400,200]) );

>

Download ToyProblem_ac2015.2.mw

And here that is in Maple 2022.1,

>

restart;

>	kernelopts(version);

>	with(plots):

>	shock := piecewise(t <0, 0, t < 1, 10, 0): sys := {(M__p+M__a)diff(x(t), t$2)=M__pshock-x(t), x(0)=0, D(x)(0)=0}

$sys := {(M__p+M__a)*(diff(x(t), t, t)) = M__p*piecewise(t < 0, 0, t < 1, 10, 0)-x(t), x(0) = 0, (D(x))(0) = 0}$

>	sol := unapply(rhs(dsolve(sys)), (M__p,M__a))

sol := proc (M__p, M__a) options operator, arrow; piecewise(t < 0, 0, t < 1, -10*M__p*cos(t/sqrt(M__p+M__a))+10*M__p, 1 <= t, -10*M__p*cos(t/sqrt(M__p+M__a))+10*M__p*cos((-1+t)/sqrt(M__p+M__a))) end proc

>	# Computation of tend for M__p=M__a=2 # # Im expecting to get a value around 6.8

>	gensol := solve({sol(2,2),t>0},t,allsolutions);

${t = -2*arctan((-1+cos(1/2))/sin(1/2))+2*Pi*_Z2}$

>	RootFinding:-NextZero(unapply(sol(2,2),t),0); secondroot := RootFinding:-NextZero(unapply(sol(2,2),t),%);

>	cands := {Student:-Calculus1:-Roots(sol(2,2),t=0..secondroot)[]}; tend := [op(cands minus {0})][1]; evalf(tend);

>	xtend := simplify(eval(sol(2,2),t=tend));

>	xmax := maximize(20cos(((t - 1)sqrt(4))/4) - 20cos(sqrt(4)t/4)); evalf(xmax);

>	tmax := Student:-Calculus1:-Roots(sol(2,2)=xmax,t=0..10)[1]; evalf(tmax);

>	display( pointplot([[tmax,xmax], [tend,xtend]], symbolsize=20, symbol=solidcircle), plot(sol(2,2),t=0..10, size=[400,200]) );

>

Download ToyProblem_ac2022.1.mw

something like this?...

Answered: acer 32373

September 07 2022

0 0

Is this the kind of thing you were trying to accomplish?

I wasn't sure whether you wanted only the floating-point approximate root where f(n)=g(n), or an exact result for that, or the next highest integer above that where f(n)>g(n).

Download ex1.mw

An exact result is also possible, in terms of the LambertW special-function.

solve({f(n)=g(n),n>5},n):

    { n = -8/ln(2)*LambertW(-1,-1/8*ln(2)) }

evalf(%);

   {n = 43.55926044}

Your Question shows an unevaluated call to intervalsolve, as a returned value. There is no such command in stock Maple. There is such a command in the Gym package (used by Danish schools). If using that command was your intention then you made two mistakes:

1) Your unevaluated return value shows that you failed to refer to its name properly. You could either use the command's fully qualified name (by including its package name in the call), or you could first load the package and then use the command's short name. See below, for examples of both.

2) The calling sequences of arguments that you passed is also incorrect for that Gym package's command. See below.

Here below is that command working in Maple 2021.2, called in a few alternate ways.

(These produce the same results in Maple 2022.1 and Maple 2020.2. You didn't tag your Question with your Maple version number; that can be useful information for us if it's an old version.)

>

restart;

>	kernelopts(version);

>	f := n -> 8*n^2:

>	g := n -> 64nlog[2](n):

>	Gym:-intervalsolve(f(n)=g(n),n=5..infinity);

[-8*LambertW(-1, -(1/8)*ln(2))/ln(2)]

>	Gym:-intervalsolve(f(n)=g(n),n=5..50);

[-8*LambertW(-1, -(1/8)*ln(2))/ln(2)]

>	with(Gym):

>	intervalsolve(f(n)=g(n),n=5..infinity);

[-8*LambertW(-1, -(1/8)*ln(2))/ln(2)]

>	intervalsolve(f(n)=g(n),n=5..50);

[-8*LambertW(-1, -(1/8)*ln(2))/ln(2)]

>

evalf(%);

Download intervalsolve_Gym.mw

`Å`...

Answered: acer 32373

September 06 2022

1 2

plot( x, labels = [ `Å`, y ],
labelfont = [ Times, 14 ] );

Alternatives include the following (the Mapleprimes backend doesn't render them all, but my Maple 2022.1 GUI does):

`&Aring;`
`&#8491;`
`&#x212B;`

There is also,
`#mo("Å");`
which renders in upright Roman, instead of italic. It's possible to get alternate fonts with this upright variant, but it's messier and in that case you might well find Tom's Unit(angstrom) suggestion easier. For the italic, the above seems to work as well with alternate font families.

See also here, if you're more generally interested in fantastic beasts and where to find them.

unapply...

Answered: acer 32373

September 05 2022

0 0

unapply(x*cos(x), x)

or, in your example,

unapply(f,x)

where x is the intended symbol that appears in the expression to which f evaluates.

list of coefficients...

Answered: acer 32373

September 04 2022

1 8

Below, CL is a list of the coefficients of powers of e^(-gamma*tau).

Is that what you're after?

Download How_to_collect_Coefficent_ac.mw

ps. Your expression contain e^(-gamma*tau) which is not the same as exp(-gamma*tau). Is that intentional?

sure...

Answered: acer 32373

September 03 2022

1 10

Sure, you do not have to write out manually all the applications of the individual trig functions to the input values.

There are several ways to do it. You could use the Matrix and Vector commands, with initializers (see Help). You could use map or seq across lists/Vectors of the functions and inputs.

Or you could apply the functions to the inputs using elementwise operations. Below I show an example of how that may construct the Matrix data of the DataFrame. It's not the only way.

The embedded GUI Table doesn't render properly on this site. That's just a Mapleprimes problem. In the actual Maple GUI it looks fine (centered, with appropriate size, etc), like your original. What I changed is how to programmatically construct the DataFrame.

NumericEventHandler('division_by_zero' = proc (f, ops, defVal) NumericStatus('division_by_zero' = false); undefined end proc)

>

F := [sin, cos, tan]; C := [0, (1/6)*Pi, (1/4)*Pi, (1/3)*Pi, (1/2)*Pi]

>


sin(x)
cos(x)
tan(x)

TableauTest_ac.mw

You could also substitute something else for undefined, into the data. Eg,
TableauTest_ac_ne.mw

ps. Your question appears to be about how to programmatically and more easily construct the DataFrame (and its data/Matrix, and its headers). Your question doesn't really seem to be about how to use Tabulate to embed that DataFrame. I've adjusted the Questions's tags accordingly.

Shift-Enter...

Answered: acer 32373

September 02 2022

0 1

You can use the keystrokes Shift-Enter to move the input cursor down to another line within the loop (without causing execution).

In the Maple's documentation (Help) that is called a soft new line here, as one of its keyboard shortcuts.

You'll still need the statement terminators such as a semicolon. (Strictly speaking it's not needed on the last statement before the end do, but I find life easier to always have it.)

I recommend indenting your code involving do-loops, if-then, etc. It helps with legibility. How many spaces to indent is up to you, but being consistent helps.

undefined values...

Answered: acer 32373

September 01 2022

0 9

[edit: I did the following before the OP suggested that the x-axis range he used is incorrect.]

Given the presence of the undefined values in the Matrix, you could do the following,

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

However, it seems that the Float(undefined) values might arise from your relatively poor approach to linear-solving for r, given A.r=B. You are computing the inverse of A and then multiplying. I get the following (Maple2016 the oldest I have at hand right now) if I instead use r:=LinearSolve(A,B).

[edit: It might even be appropropriate to perform a least-squares solve (LeastSquares with QR or SVD method, or LinearSolve with QR method) if the system is underdetermined. I don't know enough about what you're computing to pass judgement on that.]

ps. Every time you ask a Question here you omit the important detail that you are using Maple 18. (I end up adding it to your Question headers, if I remember to look at your .mw attachments -- in the cases that you provide such.) It'd be more helpful to the readership if you could bother to mark your Questions with your particular version number.

>

restart:

>

	Procedure

>

doCalc:= proc( xi , u) #u is the \bar(H): normalize magnetic field magnitude,
                          # where H = bar(H)*H__a
                 # 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*alpha/(4*k__1*q^2/d^2-alpha__3*xi/eta__1 - chi*H^2),
                       n:= 50,
                       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,soln3, Imagroot1, Imagroot2, AreTherePurelyImaginaryRoots;

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

  g:= q-(1-alpha)*tan(q)+ c*tan(q) + w*q:
  f:= subs(q = x+I*y, g):
  b1:= evalc(Re(f)) = 0:
  b2:= y-(1-alpha)*tanh(y) -(alpha__3*xi*alpha/(eta__1*(4*k__1*y^2/d^2+alpha__3*xi/eta__1)))*tanh(y) = 0;
  qstar:= (fsolve(1/c = 0, q = 0 .. infinity)):
  OddAsymptotes:= Vector[row]([seq(evalf(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 complex roots

  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;

AreTherePurelyImaginaryRoots:= type(true, 'truefalse');
  try
# Compute the rest of the Roots
  soln3:= simplify~(fsolve({ b2}, { y = 0 .. infinity}),zero);
  Imagroot1:=subs(soln3, I*y);
  Imagroot2:= -Imagroot1;
  catch:
   AreTherePurelyImaginaryRoots:= type(FAIL, 'truefalse');
  end try;

## odd and all asymptotes
  OddAsymptotes:= Vector[row]([seq(evalf((1/2)*(2*j+1)*Pi), j = 0 .. n)]);
  AllAsymptotes:= sort(Vector[row]([OddAsymptotes, qstar]));

if (xi > 0) then
  if AreThereComplexRoots and not AreTherePurelyImaginaryRoots
  then gg:= [qcomplex1, qcomplex2,op(Roots(g, q = 0.1e-3 .. AllAsymptotes[1])),
              seq(op(Roots(g, q = AllAsymptotes[i] .. AllAsymptotes[i+1])), i = 1 .. n)];
  elif AreThereComplexRoots and AreTherePurelyImaginaryRoots
  then gg:= [qcomplex1, qcomplex2, Imagroot2, op(Roots(g, q = 0.1e-3 .. AllAsymptotes[1])), seq(op      (Roots(g, q = AllAsymptotes[i] .. AllAsymptotes[i+1])), i = 1 .. n)];
  elif not AreThereComplexRoots and AreTherePurelyImaginaryRoots
  then gg:= [Imagroot2, op(Roots(g, q = 0.1e-3 .. AllAsymptotes[1])), seq(op(Roots(g, q =               AllAsymptotes[i] .. AllAsymptotes[i+1])), i = 1 .. n)];
  elif not AreThereComplexRoots and not AreTherePurelyImaginaryRoots
  then gg:= [op(Roots(g, q = 0.1e-3 .. AllAsymptotes[1])), seq(op(Roots(g, q = AllAsymptotes[i] ..      AllAsymptotes[i+1])), i = 1 .. n)];
  end if:
else
if AreThereComplexRoots and not AreTherePurelyImaginaryRoots
  then gg:= [qcomplex1, qcomplex2,op(Roots(g, q = 0.1e-3 .. AllAsymptotes[1])),
              seq(op(Roots(g, q = AllAsymptotes[i] .. AllAsymptotes[i+1])), i = 1 .. n)];
  elif AreThereComplexRoots and AreTherePurelyImaginaryRoots
  then gg:= [qcomplex1, qcomplex2, Imagroot2, op(Roots(g, q = 0.1e-3 .. AllAsymptotes[1])), seq(op      (Roots(g, q = AllAsymptotes[i] .. AllAsymptotes[i+1])), i = 1 .. n)];
  elif not AreThereComplexRoots and AreTherePurelyImaginaryRoots
  then gg:= [op(Roots(g, q = 0.1e-3 .. AllAsymptotes[1])), seq(op(Roots(g, q =                          AllAsymptotes[i] .. AllAsymptotes[i+1])), i = 1 .. n)];
  elif not AreThereComplexRoots and not AreTherePurelyImaginaryRoots
  then gg:= [op(Roots(g, q = 0.1e-3 .. AllAsymptotes[1])), seq(op(Roots(g, q = AllAsymptotes[i] ..      AllAsymptotes[i+1])), i = 1 .. n)];
  end if:
end if:

# Remove the repeated roots if any & Redefine n

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

## Compute the time constants

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

## 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] := int(Efun[i]*Efun[j], z = 0 .. d):
           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] := int(Efun[i]^2, z = 0 .. d):
  end do:

## Calculate integrals for RHS vector

  for i to m do
  F[i] := int(theta_init*Efun[i], z = 0 .. d):
  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 inverse of A and solve A*x=B

  #Ainv := 1/A:
  #r := MatrixVectorMultiply(Ainv, 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):

## compute minimum value of tau
  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,qq;
     end proc:

Run the calculation for supplied value of 'xi

>	ans:=[doCalc(0.1, 2)]: evalf(ans[9]); evalf(ans[10]);

>

Contour plot for pmin (minimum of tau) for different xi and H values

>	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;

>	(gridji,rngj,rngi):=[5,5],0..3,-2.0..-0.0; PP[gridji,4]:=plot3d(max(-1e2,testproc(i,j,4)), j=rngj, i=rngi, grid=gridji): PP[gridji,7]:= plot3d(Im(testproc(i,j,7)), j=rngj, i=rngi, grid=gridji):

>	## Surface plot using 3d PP[gridji,4]:

>	## This is the pmin values data[gridji,4]:=op([1,3],PP[gridji,4]);

>	## Convert Array to Matrix to allow the use listcontplot ConMatrix := data[gridji,4];

>	## Assign the min and max value of pmin to minv and maxv (minv,maxv):=[min[defined],max[defined]](data[gridji,4])[];

>	## GRB color and the contours (color1,color2):="Red","Yellow": conts:=[minv, seq(minv+(i-1)*(maxv-minv)/(NN+2-1),i=2..NN+1), maxv]/2.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))):

>

##

>

## Use the Conmatrix and the conts to plot contours and transform the axes
subsindets(plots:-listcontplot(ConMatrix, contours=conts,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, labels = ["u", "xi"], labelfont = ["TimesNewRoman", 14], axesfont = [14, 14]);

>	## This deals with the Imaginary part

>	## Surface plot using plot3d for Imginary part PP[gridji,7]:

>

##

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

>	## Convert Array to matrix AMatrix := Matrix(data[gridji,7]);

AMatrix := Matrix(5, 5, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (2, 1) = .0, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (2, 5) = .0, (3, 1) = .0, (3, 2) = .0, (3, 3) = -.0, (3, 4) = .0, (3, 5) = .0, (4, 1) = .0, (4, 2) = .0, (4, 3) = .0, (4, 4) = .0, (4, 5) = .0, (5, 1) = .0, (5, 2) = -.0, (5, 3) = .0, (5, 4) = .0, (5, 5) = .0}, datatype = float[8])

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

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

>

colorlist:=ColorTools:-Color~(((c1,c2,N)->[seq(c1+(i-1)*(c2-c1)/(N-1),
                                               i=1..N)])([ColorTools:-Color(color2)[]],
                                                         [ColorTools:-Color(color1)[]],
                                                         nops(conts1))):

>

subsindets(plots:-listcontplot(AMatrix,
                               contours=conts1-~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](conts1[i]),thickness=15,
                 color=colorlist[nops(conts)-i+1]),
            i=1..nops(conts1)),
        size=[500,400],legendstyle=[location=right], axes=boxed, labels = ["u", "xi"], labelfont = ["TimesNewRoman", 14], axesfont = [14, 14]);

>

Download Case3_Ijuptilk_Contourplot_pmin_290822_ac.mw

example...

Answered: acer 32373

September 01 2022

0 1

Is this the kind of thing you want to accomplish?

I did not change the way your formulas work, or which kinds of data structure you use. (Eg. it might be easier to create a as an Array, etc.)

[edit: I figured that you are trying to learn Maple coding, and that a direct movement of your original code into a procedure (which is what you asked) would help in your understanding. Some people learn better by focusing on only a few things at once. So I just put it into a reusable procedure, as is.

Look in particular at how the parameters and the locals of the procedure are all just symbols, and not any indexed names as you attempt had tried.]

>

restart;

>	Digits := 15;

>

T5Scheme:=proc(x1,y1,a,L)
local i,k,nk,N,aa,x,y;
   x[1]:=x1; y[1]:=y1;
   N:=numelems(x[1]):
   for k from 1 to L do
     nk := 3^(k - 1)*(N - 6) + 6;
     for i from 3 to nk - 2 do
       x[k + 1][3*i - 8] := a[-6]*x[k][i + 2] + a[-3]*x[k][i + 1] + a[0]*x[k][i]
                            + a[3]*x[k][i - 1] + a[6]*x[k][i - 2];
       y[k + 1][3*i - 8] := a[-6]*y[k][i + 2] + a[-3]*y[k][i + 1] + a[0]*y[k][i]
                            + a[3]*y[k][i - 1] + a[6]*y[k][i - 2];
       x[k + 1][3*i - 7] := a[-5]*x[k][i + 2] + a[-2]*x[k][i + 1] + a[1]*x[k][i]
                            + a[4]*x[k][i - 1] + a[7]*x[k][i - 2];
       y[k + 1][3*i - 7] := a[-5]*y[k][i + 2] + a[-2]*y[k][i + 1] + a[1]*y[k][i]
                            + a[4]*y[k][i - 1] + a[7]*y[k][i - 2];
       x[k + 1][3*i - 6] := a[-4]*x[k][i + 2] + a[-1]*x[k][i + 1] + a[2]*x[k][i]
                            + a[5]*x[k][i - 1];
       y[k + 1][3*i - 6] := a[-4]*y[k][i + 2] + a[-1]*y[k][i + 1] + a[2]*y[k][i]
                            + a[5]*y[k][i - 1];
     end do;
   end do;
   return evalf(convert(x[L+1],list)), evalf(convert(y[L+1],list));
end proc:

>	## Initial polygon x[1] := [-1, 1/2, 1, 0, -sqrt(2)/2, -1, 1/2, 1, 0, -sqrt(2)/2, -1, 1/2, 1, 0, -sqrt(2)/2];

>	y[1] := [0, sqrt(3)/2, 0, -1, -sqrt(2)/2, 0, sqrt(3)/2, 0, -1, -sqrt(2)/2, 0, sqrt(3)/2, 0, -1, -sqrt(2)/2];

>	s2 := [-1, 1/2, 1, 0, -sqrt(2)/2, -1];

>	t2 := [0, sqrt(3)/2, 0, -1, -sqrt(2)/2, 0];

>	assign(a[-6] = 13/1296, a[-5] = -11/648, a[-4] = -1/16, a[-3] = -107/1296, a[-2] = 179/1296, a[-1] = 9/16, a[0] = 137/144, a[1] = 137/144, a[2] = 9/16, a[3] = 179/1296, a[4] = -107/1296, a[5] = -1/16, a[6] = -11/648, a[7] = 13/1296);

>	s1,t1 := T5Scheme(x[1],y[1],a,3):

>	plot([<<s1> \| <t1>>, <<s2> \| <t2>>], linestyle = [1, 3], color = [black, red]);

Download zz123_ac.mw

file: protocol (URL)...

Answered: acer 32373

August 31 2022

1 2

It works for me, in the following sense of referencing a file on my local machine.

In the Hyperlink properties I set the Type to URL and I set the Target to,

file:///home/acer/mapleprimes/collectexample.html

and clicking on that link in my Maple GUI document makes my web-browser open that local file. (Files with extension .html extension are associated with my web-browser application.)

The key thing above was to use the file: protocol in the URL.

It's not entirely clear to me what you want to happen by linking to a geogebra file. Something with your web-browser, or do you expect some 3rd party app to launch it, or...?

defn of f...

Answered: acer 32373

August 31 2022

1 2

Look at the output of your assignment of an operator to the name f.

Clearly something has gone wrong in your 2D Input, since the A*sin(omega*t + varphi)*I term is missing.

I deleted and retyped that term,

Download Mapleprimes_Book_2_Question_3.6_ac.mw

something...

Answered: acer 32373

August 31 2022

2 15

For the purpose of plotting you may not need high accuracy in the numeric integration.

Also, if you forcibly take only the (significant) real part of calls to N, while also "disabling" overflow exceptions under evalhf and the _d01amc quadrature method, then the float evaluation of calls to J appears to speed up.

(Also, check that might simplified revision of L is correct...)

I didn't experiment with using a finite range of integration. (Maybe that could be faster, while still reasonably accurate? I didn't look...)

Perhaps you could tell us how long your original plot took to compute?

>

restart;

>	with(plots):

>	F:=kappa->kappa;

>	f:=(alpha,delta)->exp(-abs(F(kappa))^2(1+delta^2)/2-abs(F(kappa))alpha)/abs(F(kappa));

proc (alpha, delta) options operator, arrow; exp(-(1/2)*abs(F(kappa))^2*(1+delta^2)-abs(F(kappa))*alpha)/abs(F(kappa)) end proc

>	L:=(alpha,delta,Lambda) -> (lambda^2exp(-alpha^2/2)/4)(2*Int(f(alpha,delta),kappa= -infinity..-Lambda));

proc (alpha, delta, Lambda) options operator, arrow; (1/2)*lambda^2*exp(-(1/2)*alpha^2)*(Int(f(alpha, delta), kappa = -infinity .. -Lambda)) end proc

>	#0.0008209373770*lambda^2 forget(evalf); CodeTools:-Usage( evalf(L(4,1,0.001)) );

memory used=1.57MiB, alloc change=0 bytes, cpu time=29.00ms, real time=30.00ms, gc time=0ns

>	g:=(beta,delta)->exp(-Ikappabeta-abs(F(kappa))^2*(1+delta^2)/2)/abs(F(kappa));

proc (beta, delta) options operator, arrow; exp(-I*kappa*beta-(1/2)*abs(F(kappa))^2*(1+delta^2))/abs(F(kappa)) end proc

>	E:=(omega,gg)->exp(Iomegagg)(1-erf((gg+Iomega)/sqrt(2)));

proc (omega, gg) options operator, arrow; exp(I*omega*gg)*(1-erf((gg+I*omega)/sqrt(2))) end proc

>	MyHandler := proc(operator,operands,default_value) NumericStatus( overflow = false ); return 10^10; end proc: ig1template := proc(kappa) NumericEventHandler(overflow = MyHandler); evalhf(__dummy); end proc:

>

igdum:='Re'(simplify(exp(-alpha^2/2)/8*(g(beta,delta)*(E(abs(F(kappa)),gg)+E(abs(F(kappa)),-gg))))):
J := subs(__igdum=igdum, proc(alpha,delta,Lambda,beta,gg) local ig;
  ig := subs(__dummy= __igdum, eval(ig1template));
  evalf(lambda^2*abs(Int(ig,-infinity..-Lambda,epsilon=1e-4,method=_d01amc)
                     +Int(ig,Lambda..infinity,epsilon=1e-4,method=_d01amc)));
end proc):

>	forget(evalf); CodeTools:-Usage( J(4,1,0.001,8,3) );

memory used=49.92MiB, alloc change=0 bytes, cpu time=319.00ms, real time=319.00ms, gc time=0ns

>	N := (beta,alpha) -> (J(alpha,1,0.001,beta,3)-L(alpha,1,0.001))/lambda^2;

proc (beta, alpha) options operator, arrow; (J(alpha, 1, 0.1e-2, beta, 3)-L(alpha, 1, 0.1e-2))/lambda^2 end proc

>	forget(evalf); CodeTools:-Usage(contourplot(N, 0..10, 0..10, grid=[25,25]));

memory used=20.12GiB, alloc change=76.00MiB, cpu time=2.50m, real time=2.32m, gc time=19.72s

>

Download Negativity_v1_acc.mw

It might even be possible to improve the timing even further by parallelizing computation over a 2D Array (over beta and alpha values), and then using listcontplot.

browser...

Answered: acer 32373

August 31 2022

1 1

@nm You might try it with a different web browser.

A month ago I was unable to submit postings with Chromium, but Firefox worked. That problem went away after a while, as mysteriously as it arrived.

E-Mail Address:
Password:
Remember Me:	Automatically sign in on future visits

E-Mail Address:
Password:
Remember Me:	Automatically sign in on future visits

Ask a Question

Create a Post

32373 Reputation

29 Badges

Social Networks and Content at Maplesoft.com

MaplePrimes Activity

These are answers submitted by acer

with assuming...

type suffixed...

Maple 2015.2 vs Maple 2022.1...

something like this?...

`Å`...

unapply...

list of coefficients...

sure...

Shift-Enter...

undefined values...

example...

file: protocol (URL)...

defn of f...

something...

browser...

Save this setting as your default sorting preference?

Ask a Question

Create a Post

Generating PDF…

Save this setting as your default sorting preference?
Note: You can change your preference any time in your account settings
Don't show this again