How to make integral of Hankel function at infinite?

Hi guys,

I have a question on numerical integral of Hankel functions, need help from all you. Thank you in advance

How do we treat with integral of a function contains Hankel function from 0 to infinite? for example, function=x.besselh(1,0,x)?

statement "quad" not allow at infinite.

Answer 1

[Mike Hosea]

Use INTEGRAL. If you don't have R2012a or later, use QUADGK. If you don't have QUADGK, it's time to upgrade.

[tran]

Quadgk(f,0,inf) with f=@(x) x.*bessel(0,1,x) gives us warning:

Warning: Reached the limit on the maximum number of intervals in use. Approximate bound on error is 1.1e+045. The integral may not exist, or it may be difficult to approximate numerically. Increase MaxIntervalCount to 738 to enable QUADGK to continue for another iteration. > In quadgk>vadapt at 317 In quadgk at 216

ans =

3.5928e+044 -1.0851e+045i

How can we trust this result?

[Mike Hosea]

I didn't even think about the function you were integrating. I assumed it was integrable. I'm not saying QUADGK can handle any integrable problem, because there are some that it can't, but if you're going to try to integrate a function that oscillates with increasing amplitude as x increases, I don't think it matters what method you use.

Answer 2

[Matt Fig]

This one is easy. Because:

besselh(1,0,x) % Zero for all x.

we know the integral from 0 to inf is 0.

[tran]

Sorry, I made a mistake. I should be besselh(0,1,x) and thus it is not easy as you think.

[Walter Roberson]

By converting to the equivalent BesselJ and BesselK functions and evaluating those, you get 1 + 0i

[tran]

Did you get any warning from matlab? I did not get the same result. I use quadgk with version R2009b. Note that function is x*besselh(0,1,x), not only besselh

[Matt Fig]

This is a difficult integral to evaluate numerically, as Mike noted. What you can do is approximate it because you know it converges. This takes a while to run, so be warned. For greater precision you will need to adjust the conditional on the WHILE loop:

f = @(x) besselh(0,1,x);
I = 0;
A = 1;
ii = 1;

while abs(A)>1e-3
   A = quadgk(f,(ii-1)*20,ii*20); % integrate intervals and sum
   I = I + A;
   ii = ii + 1;
end

I end up with: I = 0.9991 - 0.0001i Which tends to confirm Walter's solution. You could also break the besselh into besselj and bessely then run two loops. This might be actually faster.

[tran]

This way does not give us good result and be limited. Indeed if just besselh, we can break into 2 functions j and y to do. But I want to discuss about function x.besselh(0,1,x).

[Matt Fig]

tran, I understand that the approach I showed above is limited. Are you listening to Mike's and my advice about the difficulty of integrating this particular function? Just because you can write something down does not mean that it is easy to calculate! I have offered a way to get something towards a solution of your problem. I assume something is better than nothing, which is what quadgk(f,0,inf) gives...

If you want to evaluate this integral numerically, feel free to find another approach. I don't know if another approach will work faster and/or be more accurate. If you find one, please let us know how you did it!

[tran]

Hi Matt. Thank you very much for your advice. Discussion is a very good way to gain our understanding, right? Yes, indeed s_omething is better than nothing_ . I like to hear advice from all you

[Walter Roberson]

Looks to me like the two component parts of besselh both oscillate infinitely often towards x=infinity, and it appears that although the oscillations decay that they do so much more slowly than x increases; this leads to the indeterminate (infinity times 0) + I * (infinity times 0) as the limit at infinity, which is undefined.