2 years, 32 days

## plotting the imaginary against the obtai...

@tomleslie

Hi Tom,

You are right, I want to plot the imaginary parts of the two complex numbers (- 23.3886061400*I and 23.3886061400*I) against the xi in the list.

[-10, -9.6000000000, -9.2000000000, -8.8000000000, -8.0000000000, -7.6000000000, -7.2000000000, -6.8000000000, -6.0000000000, -5.6000000000, -5.2000000000, -4.8000000000, -4.4000000000, -4.0000000000, -3.6000000000, -2.8000000000, -2.4000000000, -0.4000000000]

The idea is that once we know the xi that gives us the minimum value of Re(p_n), then we can plot the Im(p_n) for that xi values.

Note that imaginary for the other xi values are zeros except for the two above.  The reason for doing this is that I want to find out if the Im(p_n) is ever non-zero for these values of xi listed above.

## Question...

Sorry again, how do I plot the Im(p_n) for the values of xi we got from this: nstart := [ seq( if( member( complex(extended_numeric), whattype~({entries(ans[j,5], nolist)})), ans[j,1], NULL), j=1..op([2,1,2], ans))];

Because I added the seq( Im(gamma1*alpha/(4*k__1*qq[i]^2/d^2-alpha3*xi/eta__1)), i=1..m)  to the return list, and I got 5 entries in the ans.

And if I want to plot nstart vs convert(ans[..,5], they don't have the same length.

xi = [-10, -9.6, -9.2, -8.8, -8.0, -7.6, -7.2, -6.8, -6.0, -5.6, -5.2, -4.8, -4.4, -4.0, -3.6, -2.8, -2.4, -.4] ( i.e the xi  for real p_n).

## Question...

Thank you very much. I really appreciate this.

## Question...

Hi Tom,

Please could help me to explain what each of these codes does? Sorry, I like said I'm using maple for such for the first time.

[ seq( if( member( complex(extended_numeric), whattype~({entries(ans[j,5], nolist)})), ans[j,1], NULL), j=1..op([2,1,2], ans))];

member( complex(extended_numeric)
whattype~({entries(ans[j,5], nolist)}))
ans[j,1], NULL)
op[2,12]

Also, how do I get the correspond imaginary p_n for the list the xi-values which return complex entries?

## Question...

@tomleslie

Hi Tom,

Please, how can extract all the xi that give complex values of p_n?

## Question...

@tomleslie

Kindly disregard the last question I asked. I understand now.

Thank you so much for your assistance.

## @tomleslieHi,In the procedure, you ...

@tomleslie

Hi,

In the procedure, you returned 3 outputs:

min(Re), theta_int , vol_sol;

On plotting min(Re) vs xi, you used

ans:= Array([seq( [j, doCalc(j)], j=-10..10, 0.4)]):

plot(convert(ans[..,1],list), convert(ans[..,2],list))

If doCalc(j) corresponds to convert(ans[..,2]), and convert(ans[..,2]) corresponds to min(Re).
Why did you use index 2 to access the returned value instead of 1, since min(Re) is the first entry in the return list?

## How to include a plot in a procedure...

Hi Tomleslie,

Thank you so much. This is exactly what I wanted Please, how can learn, maple at the advanced level?

## How to include a plot in a procedure...

Hi,

I want it evaluated for several values of 't' and plotted over a range of z-values, thus producing mulitple curves in a single 2D plot. If so what are the t-values? You plot for t =1..10 with size of 0.1.

## How to include a plot in a procedure...

Hi,

The plot theta_int is perfect because it doesn't depend on \xi. My major concern is actually plotting theta_sol and vol_sol.  The theta_sol and v_sol are functions of (z,t). Don't worry about the complex you see, because is the complex and its conjugate, which is real. I will appreciate it if you could plot theta_sol or v_sol for a particular value of \xi. tried to run the code after the changed theta_int to theta_sol.

## How to include a plot in a procedure...

Hi Tomleslie,

Thank I really understand how the procedure works now.

But,  please how do include a plot in the procedure? I have added another code to the earlier codes. I want to plot theta_sol and v_sol  in the code below. I also what the procedure to return Real and Imaginary instead of the only real part.

New_cutPaste_(1).mw

Sorry for bugging you.

## finding minimum...

Thank you so muc for your assistance

## finding minimum...

@tomleslie

Thank you for your reply. This is my first to using maple extensively, that's why I have to take it bit by bit.  I tried to use

x :=[(seq(0.1*i, i=-10..12))];

pp_min:= [seq (min( seq( Re(gamma1*alpha/(4*k[1]*qq[i]^2/d^2-alpha3*x[j]/eta[1])), i=1..m) ), j= 1..m)];
display(plot(x, pp_min));

to sort the problem.

However, my supervisor said that is not correct, because I only use one value of xi in qq[i]. So I had to result to define xi using sequence from the beginning and use the result to compute qq[i]. However, I encounter a challenge. I received an error on studentcalculus ().

Find attached the problem:

Error_student_calculus_package.mw

Import packages

 > restart: with(ArrayTools): with(Student:-Calculus1): with(LinearAlgebra): with(ListTools):with(RootFinding):with(ListTools):

Parameters

 > n:= 5: gamma1 := .1093: alpha3 := -0.1104e-2: k[1] := 6*10^(-12): d:= 0.2e-3: x:=[(seq(0.4*i, i=-2..n))]; theta0:= 0.1e-3: eta[1]:= 0.240e-1: alpha:= 1-alpha3^2/(gamma1*eta[1]): c:= seq(alpha3*x[i]*alpha/(eta[1]*(4*k[1]*q^2/d^2-alpha3*x[i]/eta[1])),i=1..n): theta_init:= theta0*sin(Pi*z/d):
 (2.1)

Assign g for q and plot

 > g := seq(q-(1-alpha)*tan(q)-c[i]*tan(q), i=1..n): plot(g[1], q = 0 .. 3*Pi, view = [DEFAULT, -30.. 30]):

Set q as a complex

Assume  and subsitute the result into g and equate the real and complex part to zero, and solve for x and y.

 > f := seq(subs(q= x+I*y, g[i]), i=1..n ): b1 := seq(evalc(Re(f[i])) = 0, i=1..n): b2 := seq(evalc(Im(f[i])) = 0, i=1..n):

Compute the Special Asymptotes

This asymptote is coming from the  from the definition of

 > qstar := Vector[row]([seq(fsolve(1/c[i] = 0, q = 0 .. infinity),i=1..n), Vector[row](n,shape=zero)]);
 (4.1.1)

Compute Odd asymptote

First, Since , then an asymptote occurs at  In general, we have

Next, we compute the entry of the Oddasymptotes that is close to qstar (special asymptote) as assign it to
ModifiedOaddAsym, and then find the minimum of the ModifiedOaddAsym. Searchall Function returns

the index of an entry in a list.

 > OddAsymptotes := Vector[row]([seq(evalf((1/2*(2*j+1))*Pi), j = 0 .. n)]); ModifiedOddAsym := seq(abs(~[-](OddAsymptotes, qstar[j])),j=1..n); qstarTemporary := seq(min(ModifiedOddAsym[i]),i=1..n); indexOfqstar2 := seq(SearchAll(qstarTemporary[i], ModifiedOddAsym[i]),i=1..n); qstar2 := seq(OddAsymptotes(indexOfqstar2[i]), i=1..n);
 (4.2.1)

Compute x and y

Here, we solve for and  within the min. and max. of qstar2 and qstar, and substitute the results into .

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

Compute the rest of the Roots

In this section we compute the roots between each asymptotes.

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

Remove the repeated roots if any

 > #qq := MakeUnique(gg):

Redefine n

 > m := numelems(qq);
 (4.6.1)

Compute the time constants

 > for i to m do p[i] := gamma1*alpha/(4*k[1]*qq[i]^2/d^2-alpha3*xi/eta[1]); 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: