2
$\begingroup$

I want a contour plot of the complex roots of the following fifth degree dispersion equation in wave velocity v in Mathematica:

[1+(c2^2/c1^2) v^2][1+k t1 c2 v+(c2 d/ k) v(1+ k t2 c2 v+(1/2) k^2 t2^2 c2^2 v^2)]+(e c2 d/ k) v[1+k t2 c2 v+(1/2) k^2 t2^2 c2^2 v^2] = 0;

The parameters are given as:

c1=4631.0; c2=2280.1; d=8066.8; e=0.0168; t1=0.1; t2=0.1; k=1;

Please someone help me. I am learning Mathematica but I need the above help for my research work.

Improve this question
$\endgroup$
4
  • $\begingroup$ Initialize the parameters first. Then correct your expression: Change brackets to parentheses. Brackets are used for functions and arrays in Mathematica, next, simpler to remove the =0, then just assign your expresson to a name such as theFuntion=. So will have theFuntion=expression here. Then finally use NSolve as in theRoots=v/.NSolve[theFunction==0,v] to obtain the roots. Don't know what you mean abut contour plot of the complex roots though. $\endgroup$
    – josh
    Commented Oct 2, 2023 at 12:17
  • $\begingroup$ Thank you for your suggestion. Here v is complex in general. I want Contour plot of complex roots v. $\endgroup$
    – Nantu Sarkar
    Commented Oct 2, 2023 at 12:21
  • $\begingroup$ @Nantu: Suggest you make the changes I suggested, run the code in a notebook and post both the code and roots back into your thread by copy/paste into your post then highlight the code and select {} in the menu above the edit window to format it as code then better explain what you wish to do with them. $\endgroup$
    – josh
    Commented Oct 2, 2023 at 12:54
  • $\begingroup$ Do you mean plot $v(k)$ or just to show 5 roots at k=1? $\endgroup$
    – Alex Trounev
    Commented Oct 2, 2023 at 13:27

2 Answers 2

Reset to default
3
$\begingroup$

This answers your other question that was deleted as a duplicate of this one.

$Version

(* "13.3.1 for Mac OS X ARM (64-bit) (July 24, 2023)" *)

Clear["Global`*"]

eqn = (1 + (c2^2/c1^2)*v^2)*(1 + 
       k*t1*c2*v + (c2*d/k)*v*(1 + k*t2*c2*v + (1/2)*k^2*t2^2*c2^2*v^2)) + (e*
       c2*d/k)*v*(1 + k*t2*c2*v + (1/2)*k^2*t2^2*c2^2*v^2) == 0;

c1 = 4631;
c2 = 22801/10;
d = 80668/10;
e = 168*10^-4;
k = 1;

roots = Solve[eqn, v];

n = Length@roots;

EDIT:

Using ScalingFunctions

Manipulate[
 Plot3D[
  Evaluate[func[v /. roots[[k]]]],
  {t1, 10^-10, 10^-1}, {t2, 10^-10, 10^-1},
  AxesLabel -> (Style[#, 14] & /@ {t1, t2, func}),
  ScalingFunctions -> {"Log", "Log", None},
  PlotRange -> All,
  WorkingPrecision -> 20,
  PlotPoints -> 50,
  MaxRecursion -> 5],
 {{func, Re}, {Re, Im, Abs},
  ControlType -> RadioButtonBar},
 {{k, 1, "root"}, Range[n],
  ControlType -> RadioButtonBar},
 SynchronousUpdating -> False,
 TrackedSymbols :> {func, k}]

enter image description here

With manual scaling

Manipulate[
 Plot3D[
  Evaluate[func[v /. roots[[k]] /.
     {t1 -> 10^logt1, t2 -> 10^logt2}]],
  {logt1, -10, -1}, {logt2, -10, -1},
  AxesLabel -> (Style[#, 14] & /@
     {Log[t1], Log[t2], func}),
  PlotRange -> All,
  WorkingPrecision -> 20,
  PlotPoints -> 50,
  MaxRecursion -> 5],
 {{func, Re}, {Re, Im, Abs},
  ControlType -> RadioButtonBar},
 {{k, 1, "root"}, Range[n],
  ControlType -> RadioButtonBar},
 SynchronousUpdating -> False,
 TrackedSymbols :> {func, k}]

enter image description here

Improve this answer
$\endgroup$
4
  • $\begingroup$ This is nice answer (+1}. $\endgroup$
    – Alex Trounev
    Commented Oct 3, 2023 at 3:15
  • $\begingroup$ Dear Bob Hanlon, Sir thank you very much for your code. It works on my PC but have some problem for different t1 and t2. I need another help from you sir - Actually, according to my model, the ranges of t1 and t2 are [10^-10, 10^-1]. I want to plot the real and imaginary parts of each root on the logarithmic scales of t1 and t2 and indicate log(t1), log(t2) and Re/Im parts of the particular root along the x-axis, y-axis and z-axis, respectively in the 3D plot. It is my kind request to you please help me sir. $\endgroup$
    – Nantu Sarkar
    Commented Oct 4, 2023 at 11:02
  • $\begingroup$ @BobHanlon: Sir, thank you very much for your kind help. The real part of any one of the roots should be like this posted in the link mathematica.stackexchange.com/questions/291185/… I can't get exact figure! But i must appreciate your effort in helping me. If please go through the link pasted above once if possible! $\endgroup$
    – Nantu Sarkar
    Commented Oct 5, 2023 at 8:22
  • $\begingroup$ The target plot has the axes reversed and the plot range restricted (Automatic vs All). $\endgroup$
    – Bob Hanlon
    Commented Oct 5, 2023 at 22:27
4
$\begingroup$

There are 5 roots v=x+Iy which could be defined as follows

c1 = 4631.0; c2 = 2280.1; d = 8066.8; e = 0.0168; t1 = 0.1; t2 = 0.1; k = 1;

 f = (1 + (c2^2/c1^2) v^2) (1 + 
       k t1 c2 v + (c2 d/k) v (1 + 
          k t2 c2 v + (1/2) k^2 t2^2 c2^2 v^2)) + (e c2 d/k) v (1 + 
       k t2 c2 v + (1/2) k^2 t2^2 c2^2 v^2) /. v -> x + I y;
r1 = FindRoot[{Re[f] == 0, Im[f] == 0}, {x, .1}, {y, .1}]

{x -> -5.34699*10^-8, y -> 4.12598*10^-26}

r2 = FindRoot[{Re[f] == 0, Im[f] == 0}, {x, -.1}, {y, -.1}]

{x -> -0.00438575, y -> -0.0043858}

r3 = FindRoot[{Re[f] == 0, Im[f] == 0}, {x, -.1}, {y, .1}]

{x -> -0.00438575, y -> 0.0043858}
 r4 = 
 FindRoot[{Re[f] == 0, Im[f] == 0}, {x, 1/10}, {y, -3}, 
  AccuracyGoal -> 4]

{x -> -4.12454*10^-15, y -> -2.04804}

r5 = 
 FindRoot[{Re[f] == 0, Im[f] == 0}, {x, 1/10}, {y, 3}, 
  AccuracyGoal -> 4]

{x -> -4.12437*10^-15, y -> 2.04804}

Visualization

Show[ContourPlot[{Re[f] == 0, Im[f] == 0}, {x, -1/50, 
   1/50}, {y, -1/50, 1/50}], 
 Graphics[{{Red, PointSize[.015], Point[{x, y}] /. r1}, {Red, 
    PointSize[.015], Point[{x, y}] /. r2}, {Red, PointSize[.015], 
    Point[{x, y}] /. r3}}]]

Figure 1

Show[ContourPlot[{Re[f] == 0, Im[f] == 0}, {x, -5, 5}, {y, -5, 5}], 
 Graphics[{{Red, PointSize[.015], Point[{x, y}] /. r4}, {Red, 
    PointSize[.015], Point[{x, y}] /. r5}}]]

Figure 2

If v depends on t1, t2 then see answer @Bob Hanlon.

Improve this answer
$\endgroup$
6
  • $\begingroup$ Sir thank you very much for your code. It works nice on PC. I need another help from you sir - As you told me in your previous comment, is it possible to contour/density plot of the five roots v(k) as a function of k? Range of k is (0,1]. If possible, how to indicate the extrema of v(k) with mentioning max/min of v and the corresponding point? If so, please help me. $\endgroup$
    – Nantu Sarkar
    Commented Oct 4, 2023 at 10:54
  • $\begingroup$ Do you mean max/min Abs[v[k]]? $\endgroup$
    – Alex Trounev
    Commented Oct 4, 2023 at 12:07
  • $\begingroup$ Sir, I actually mean max/min of Re[v[k]] and Im[v[k]] assuming v is dependent on k as you suggested. I think it would be a nice idea to show where the wave velocity and attenuation attain max/min in density plot. Please help sir. $\endgroup$
    – Nantu Sarkar
    Commented Oct 4, 2023 at 14:10
  • 1
    $\begingroup$ @NantuSarkar in general, it is better to ask these separate & new questions as separate & new questions. You can do this by making new question posts that link to this one, rather than by making additional requests within the comments of answers to your questions. $\endgroup$
    – CA Trevillian
    Commented Oct 5, 2023 at 2:42
  • $\begingroup$ @NantuSarkar Could you please start a new topic with your requests? $\endgroup$
    – Alex Trounev
    Commented Oct 5, 2023 at 4:37

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.