Advanced Search

ijuptilk

100 Reputation

3 Badges

2 years, 264 days
Contact ijuptilk

MaplePrimes Activity


These are replies submitted by ijuptilk

plotting the imaginary against the obtai...

ijuptilk 100
April 10 2022
0 0

@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...

ijuptilk 100
April 09 2022
0 0

@tomleslie 

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

ijuptilk 100
April 09 2022
0 0

@tomleslie 

Thank you very much. I really appreciate this.

Question...

ijuptilk 100
April 09 2022
0 0

@tomleslie 

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

ijuptilk 100
April 09 2022
0 0

@tomleslie

Hi Tom,

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

Question...

ijuptilk 100
April 08 2022
0 0

@tomleslie

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

Thank you so much for your assistance. 

@tomleslieHi,In the procedure, you ...

ijuptilk 100
April 08 2022
0 0

@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...

ijuptilk 100
April 07 2022
0 0

@tomleslie 

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

ijuptilk 100
April 07 2022
0 0

@tomleslie 

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

ijuptilk 100
April 07 2022
0 0

@tomleslie 

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

ijuptilk 100
April 06 2022
0 0

@tomleslie 

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

ijuptilk 100
April 06 2022
0 0

@tomleslie 

Thank you so muc for your assistance

finding minimum...

ijuptilk 100
April 06 2022
0 0

@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):

[-.8, -.4, 0., .4, .8, 1.2, 1.6, 2.0]

(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 q = x+I*y 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 c from the definition of "q."NULL

> 

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

Vector[row](%id = 18446744073892603774)

(4.1.1)

NULL

Compute Odd asymptote

 

First, Since tan*q = sin*q*(1/(cos*q)), then an asymptote occurs at cos*q = 0. In general, we have
"q= ((2 k+1)Pi)/(2). "
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);

OddAsymptotes := Vector[row](6, {(1) = 1.570796327, (2) = 4.712388981, (3) = 7.853981635, (4) = 10.99557429, (5) = 14.13716694, (6) = 17.27875960})

 

Vector[row](%id = 18446744073892604134), Vector[row](%id = 18446744073892604374), Vector[row](%id = 18446744073892604614), Vector[row](%id = 18446744073892604854), Vector[row](%id = 18446744073892605094)

  

0.22421550e-1, .825360260, 1.570796327, 1.570796327, 1.570796327

 

3, 2, 1, 1, 1

 

7.853981635, 4.712388981, 1.570796327, 1.570796327, 1.570796327

(4.2.1)

Compute x and y

 

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

> 

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;
 

true

 

false

(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](6, {(1) = 1.570796327, (2) = 4.712388981, (3) = 7.853981635, (4) = 10.99557429, (5) = 14.13716694, (6) = 17.27875960})

 

Vector[row](%id = 18446744073892608710), Vector[row](%id = 18446744073892608950), Vector[row](%id = 18446744073892609190), Vector[row](%id = 18446744073892609430), Vector[row](%id = 18446744073892609670)

  

Error, (in Student:-Calculus1:-Roots) unexpected option(s): q = 0.1e-3 .. (Vector[row](7, {(1) = 1.570796327, (2) = 4.712388981, (3) = 7.831560085, (4) = 7.853981635, (5) = 10.99557429, (6) = 14.13716694, (7) = 17.27875960}))

 
> 

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

NULL

Remove the repeated roots if any

 
> 

#qq := MakeUnique(gg):

Redefine n

 
> 

m := numelems(qq);

6

(4.6.1)

Compute the `τ_n`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

 

NULL

> 

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:

NULL

NULL

NULL

``

Download Error_student_calculus_package.mw

finding minimum...

ijuptilk 100
April 05 2022
0 0

@tomleslie 

Hi Tomleslie,

Please, how can I get various pp for different xi, say xi = -2..5, and then plot pp versus xi?

finding minimum...

ijuptilk 100
April 05 2022
0 0

@tomleslie 

Great, thanks

4 5 6 7 8 Page 6 of 8