Bessel functions of order zero

Hello everyone. I try to plot a figure of a journal article. I gave the equations, the expected figure and my incorrect figure in the attached image file. I wrote the code given below to obtain the figure on the left but it is wrong. May you help me to find my error? Thanks a lot!

phi = linspace(0.001, 1000, 1000000);
v = 2;
lRp = sqrt((1./phi).*((v+1)^2-1));
I1 = besselj(1,lRp);
I0 = besselj(0,lRp);
h = 1 - 2./lRp*(I1/I0);
plot(phi,h)
set(gca, 'XScale', 'log')

[Walter Roberson]

I doubt that you want to do a least-squared fitting between I1 and I0, so I suspect

h = 1 - 2./lRp.*(I1./I0);

However, this looks nearly flat. It is not actually flat: there are discontinuities at the places where I1 crosses 0, but it is difficult to see them unless you zoom in.

[Hasan Celik]

Thanks! I realized I forgot the dot between I1 and I0.

Answer 1

[David Goodmanson]

Hello Hasan,

Per the reference, try the same code with

besseli(n,lrp)

[modified bessel function ] instead of besselj. These are basically regular bessel functions of imaginary argument,

I_0(x) = J_0(i*x)
I_1(x) = -i*J_1(i*x)

[Hasan Celik]

Thanks!

Changing (I1/I0) to (I1./I0) worked fine.