# How to find the Eigenvalues for complicated Bessel boundary value problem?

I am trying to find energy eigenvalues in a problem. A function depends on the radius R (of a quantum dot) and the energy. The function was obtained from matching wavefunctions at the boundary of the quantum dot. The problem is complicated. The function contains energy "outside" and "inside" Bessel functions' arguments. I used FindRoot to find eigenvalues. I was able to solve the problem using @Range[#]. I would like to instead find eigenvalues by varying the energy smoothly over a value range of choice. I am also not sure that my approach is correct. Perhaps the problem's solution may have to be completely different.

This is the code I have thus far

1. Define the equation as a function of R and \[Epsilon]:
    f[R_, \[Epsilon]_] =
Sqrt[(\[Epsilon] - 0.1)/(\[Epsilon] - 0.3)]
BesselJ[-(1/2), R Sqrt[(\[Epsilon] + 0.2)^2 - 0.1^2]]/
BesselJ[1/2, R Sqrt[(\[Epsilon] + 0.2)^2 - 0.1^2]] -
Sqrt[(\[Epsilon] + 0.1)/(\[Epsilon] - 0.1)]
HankelH1[-(1/2), R Sqrt[\[Epsilon]^2 - 0.1^2]]/
HankelH1[1/2, R Sqrt[\[Epsilon]^2 - 0.1^2]];

1. Solve for [Epsilon] over a range of R values:
    results =
Table[{R, \[Epsilon] /.
FindRoot[f[R, \[Epsilon]] == 0, {\[Epsilon], #}] & /@
Range[20]} , {R, 10, 20, 1}];

1. Separate the results for plotting (as eigenvalues have real and imaginary parts):
rePts = Flatten[Table[{x[[1]], Re[y]}, {x, results}, {y, x[[2]]}], 1];
imPts = Flatten[Table[{x[[1]], Im[y]}, {x, results}, {y, x[[2]]}], 1];

1. Plot:
Show[ListPlot[rePts, PlotStyle -> Blue, PlotLegends -> {"Real Part"},
AxesLabel -> {"R", "\[Epsilon]"},
PlotLabel -> "Re(\[Epsilon]) Im(\[Epsilon]) vs R"],
ListPlot[imPts, PlotStyle -> Red, PlotLegends -> {"Imaginary Part"}],
PlotRange -> All]

• What do you mean by, "find eigenvalues by varying the energy smoothly over a value range"? Feb 5 at 17:06
• Dear bbgodfrey! Thanks for your answer. I believe the solution you suggested - for which I am thankful - is varying the energy smoothly over a range of values. So, your intuition for what my question was asking was correct! Thanks again! Feb 6 at 14:42

Plotting f is a good starting point for finding its zeroes, here for R = 10.

cplot = ComplexPlot[f[10, x], {x, -3/5 I, 20}, AspectRatio -> 1/3,
ImageSize -> Large, LabelStyle -> {10, Bold, Black},
FrameLabel -> {Re[\[Epsilon]], Im[\[Epsilon]]}]


Zeroes should be expected where multiple colors come together, according to the ComplexPlot documentation; here where red shifts to yellow. There are 63 such points on this plot. An obvious choice for finding these points is

NSolve[f[10, x] == 0 && 0 <= Re[x] <= 20 && -1 <= Im[x] <= 1, x]


but it is slow and misses several zeroes. The most effective approach I have found is that used by the OP but with many more sample points.

Union[Table[x /. Quiet@Check[FindRoot[f[10, x], {x, i - I/100}], x -> Nothing],
{i, 1/10, 20, 1/10}], SameTest -> (Abs[#1 - #2] < 10^-6 &)];
ptplot = ComplexListPlot[%, AspectRatio -> 1/3, PlotRange -> {{0, 20}, {-.6, 0}},
ImageSize -> Large, PlotStyle -> {PointSize[Medium], Black},
LabelStyle -> {10, Bold, Black}, AxesLabel -> {Re[\[Epsilon]], Im[\[Epsilon]]}]


Superimposing the two plots validates that the zeroes lie at the red to yellow transitions.

Show[cplot, ptplot]