0
$\begingroup$

I have a complicated function that I would like to solve for the roots. The function comprises Bessel and Hankel functions:

freq = 5*^9;
Omg = 2*Pi*freq;
rho = 15 ;
vco = 5907.5 (*core velocity*);
vcl = 5944 (*clad velocity*);
m = 0 ;
U = Sqrt[Omg^2/vco^2 - Omg^2/v];
W = Sqrt[Omg^2/vcl^2 - Omg^2/v];

NSolve[U*BesselJ[m + 1, U*rho]*HankelH1[m, W*rho] - W*HankelH1[m + 1, W*rho]*BesselJ[m, U*rho]==0, v]

However, NSolve does not evaluate. I have also tried using FindRoot but I have no idea what starting value to give. Plotting the function (with all variables equal to 1 for simplicity) is of no help for guessing starting values.

enter image description here

$\endgroup$
1
  • $\begingroup$ I'm getting the feeling that there are no roots. Is that possible? $\endgroup$ – Michael E2 Nov 22 '20 at 17:48
1
$\begingroup$
$Version

(* "12.1.1 for Mac OS X x86 (64-bit) (June 19, 2020)" *)

Clear["Global`*"]

freq = 5*^9;
Omg = 2*Pi*freq;
rho = 15;
vco = 5907.5 (*core velocity*);
vcl = 5944 (*clad velocity*);
m = 0;
U = Sqrt[Omg^2/vco^2 - Omg^2/v];
W = Sqrt[Omg^2/vcl^2 - Omg^2/v];

EDIT 2: Replacing Rationalize with Rationalize[#, 0]& to ensure conversion to rationals.

expr = U*BesselJ[m + 1, U*rho]*HankelH1[m, W*rho] - 
     W*HankelH1[m + 1, W*rho]*BesselJ[m, U*rho] // Rationalize[#, 0] & // 
   FullSimplify;

(sol = NSolve[{expr == 0, 10^7 < v < 10^8}, v, WorkingPrecision -> 50]) // N

(* {{v -> 3.533*10^7}} *)

expr /. sol // Quiet

(* {0.*10^-198964 + 0.*10^-198964 I} *)

In the definitions for U and W should the v be v^2? That provides a more reasonable value

Sqrt[v] /. sol // N

(* {5943.9} *)

EDIT: Using v^2

U = Sqrt[Omg^2/vco^2 - Omg^2/v^2];
W = Sqrt[Omg^2/vcl^2 - Omg^2/v^2];

expr = U*BesselJ[m + 1, U*rho]*HankelH1[m, W*rho] - 
     W*HankelH1[m + 1, W*rho]*BesselJ[m, U*rho] // 
    Rationalize[#, 0] & // FullSimplify;

(sol = NSolve[{expr == 0, 5000 < v < 6000}, v, 
    WorkingPrecision -> 50]) // N

(* {{v -> 5943.99}} *)
$\endgroup$
4
  • $\begingroup$ Yes it should be v^2. However, I thought introducing a square root would complicate things so I just called the whole thing v. Well, two of the solution are just the parameters vco and vcl. But there is a new third solution that I didn't have before. I don't see why // Rationalize // FullSimplify would help here. It seems the real helper was the 10^7 < v < 10^8 you added. $\endgroup$ – pmac Nov 21 '20 at 22:19
  • $\begingroup$ And I tried with v^2 and it gives different solutions, so I think using v instead really helps the numerical solver $\endgroup$ – pmac Nov 21 '20 at 22:21
  • $\begingroup$ Also, how do you allow for complex solutions? $\endgroup$ – pmac Nov 21 '20 at 22:35
  • $\begingroup$ Already upvoted for the method, but I think it's a harder problem than NSolve can handle. Consider Plot[RealExponent@expr, {v, 5800, 6000}, WorkingPrecision -> 100]. The gap in the graph is due to a singularity at expr /. v -> 5944. $\endgroup$ – Michael E2 Nov 22 '20 at 15:49

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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