A problem with integration of modified Bessel function

[mozam]

Hello,

In my work, I have to solve the following integral: \int {exp(-aX^2)I_0(b\sqrt(cX^2+dX+e))}dX
where I_0() is the modified Bessel function. I did not find the solution in any table of integral.

Any help is appreciated.

Thanks a lot in advance.

 

[jackmell]

How about a formal approach you probably won't like:

Let:

I_0(z)=\sum_{k=0}^{\infty}\frac{\left(1/4 z^2\right)^k}{\left(k!)\right)^2}

then you could re-write the integral as:

\int e^{-ax^2}I_0(b\sqrt{cx^2+dx+e})dx=\int e^{-ax^2}\sum_{k=0}^{\infty}\frac{\left(1/4 b^2(cx^2+dx+e\right)^k}{\left(k!\right)^2}dx

so that ends up being polynomial terms in that sum so if it is legitimate to switch the order of integration and summation, I could re-write it as:

\sum_{k=0}^{\infty} \int e^{-ax^2} P_k(x) dx

Some parts of those terms can be integrated directly I assume while others I would think, have to be expressed in terms of the expintegral. Still though, would be interesting if say five or ten terms of that sum gives decent results for some range of integration but I suppose the rate of convergence would be heavily dependent on the size of a.

 

[mozam]

Thank you Jackmell for your reply.

I tried the way you suggested, but it is quite difficult to have a closed form expression of the integral. The expression of P(x) that I obtain is:
<br /> P_{k}(x)=\frac{\left(b^2(cx^2+dx+e\right)^k}{4\left(k!\right)^2} = \frac{(b^{2k})}{4\left(k!\right)^2}\left(cx^2+dx+e\right)^k<br />
with \left(cx^2+dx+e\right)^k = \sum_{i=0}^{k} \sum_{j=0}^{k-i}\dbinom{k}{i}\dbinom{k-i}{j}e^{k-i-j}d^{j}c^{i}x^{2i+j}.

I still don't see how to approximate this integral.

 

[jackmell]

Ok, didn't think to expand it in a double sum. Now, let me suggest something I'm not sure of but figuring if it's ok is part of the fun of doing math. So we have terms of the form:

e^{-ax^2} x^n

now:

\int e^{-ax^2}x^n dx=-1/2 x^{1+n} (ax^2)^{1/2(-1-n)} \Gamma\left(\frac{1+n}{2},ax^2\right)

Then can we not say:

<br /> \int e^{-ax^2}I_0(b\sqrt{cx^2+dx+e})dx<br />
<br /> =\sum_{k=0}^{\infty}\frac{b^{2k}}{4(k!)^2}\sum_{i=0}^{k}\sum_{j=0}^{k-i}\binom{k}{i}\binom{k-i}{j}e^{k-i-j}d^j c^i \left\{-1/2 x^{1+n}(ax^2)^{1/2(-1-n)}\Gamma\left(\frac{1+n}{2},ax^2\right)\right\},\quad n=2i+j<br />

Now, I'm not saying that's right. I'm just asking is it or is it close and I just made some minor errors? Tell you what, if it was mine, you can bet I'd be plugging that into Mathematica to see if that's "appears" to be converging to numerical results in some reasonable interval for reasonable values of the constants.

 

[mozam]

Sorry, one more precision:
<br /> -\infty &lt; x &lt; \infty<br />
So that,
<br /> \sum_{j=0}^{k-i}(...)\int e^{-ax^2}x^{2i+j} dx = \sum_{p=0}^{k-i}(...)\{\int_{-\infty}^{\infty} e^{-ax^2}x^{2(i+p)} dx + \int_{-\infty}^{\infty} e^{-ax^2}x^{2(i+p)+1} dx\}<br />
given that j is either odd or even,i.e., j=2p or j=2p+1.
Also,
<br /> \int_0^\infty{x^{2(i+p)} e^{-a x^2}\,dx} = \frac{(2(i+p))!}{(i+p)! 2^{2(i+p)+1}} \sqrt{\frac{\pi}{a^{2(i+p)+1}}}<br />
And
<br /> \int_0^\infty{x^{2(i+p)+1} e^{-a x^2}\,dx} = \frac{(i+p)!}{2 a^{i+p+1}}<br />
But the expression obtained does not help me so much because another thing is that the coefficients (except a and b) are function of another variable y. let say for instance c(y), d(y) and e(y). That's why I expect an expression of the integral so that I can later integrate with respect to y.

Thank you for any suggestion

 

[jackmell]

Ok. But just for the record, that formula I gave above should read:

\int e^{-ax^2}I_0(b\sqrt{cx^2+dx+e})dx
=\sum_{k=0}^{\infty}\frac{b^{2k}}{4^k(k!)^2}\sum_{i= 0}^{k}\sum_{j=0}^{k-i}\binom{k}{i}\binom{k-i}{j}e^{k-i-j}d^j c^i \left\{-1/2 x^{1+n}(ax^2)^{1/2(-1-n)}\Gamma\left(\frac{1+n}{2},ax^2\right)\right\},\quad n=2i+j

And I tested it in Mathematica using the first 15 terms cus' well I like math too:

In[218]:=
myf[a_, b_, c_, d_, e_, x_] := Exp[(-a)*x^2]*
    BesselI[0, b*Sqrt[c*x^2 + d*x + e]]; 
a = 1; 
b = 2; 
c = 3; 
d = 4; 
e = 5; 
x = 10; 
high = N[Sum[(b^(2*k)/(4^k*k!^2))*
      Sum[Binomial[k, i]*Binomial[k - i, j]*e^(k - i - j)*d^j*
        c^i*((-2^(-1))*x^(1 + 2*i + j)*(a*x^2)^
          ((1/2)*(-1 - 2*i - j))*Gamma[(1 + 2*i + j)/2, 
          a*x^2]), {i, 0, k}, {j, 0, k - i}], {k, 0, 15}]]; 
x = 1; 
low = N[Sum[(b^(2*k)/(4^k*k!^2))*
      Sum[Binomial[k, i]*Binomial[k - i, j]*e^(k - i - j)*d^j*
        c^i*((-2^(-1))*x^(1 + 2*i + j)*(a*x^2)^
          ((1/2)*(-1 - 2*i - j))*Gamma[(1 + 2*i + j)/2, 
          a*x^2]), {i, 0, k}, {j, 0, k - i}], {k, 0, 15}]]; 
N[high - low]
mynum = NIntegrate[myf[a, b, c, d, e, v], {v, 1, 10}]

Out[228]= 92.2159

Out[229]= 92.216

Which surprised me that it's so close.

Also Mozam, in your particular case, you may have to just integrate it numerically or at least do so numerically for a range then do a curve fit of the data to approximate the symbolic solution.

 

[mozam]

I see what you mean. Lot of work in perspective...

Another idea that I had is the use of classical integration methods and some known integrals like:

<br /> \int_0^{\infty} e^{-ax}J_0(b\sqrt{cx^2+2dx})dx=\frac{1}{\sqrt(a^2+b^2)}exp\left[d\left(a-\sqrt(a^2+b^2)\right)\right]<br />

or

<br /> \int_0^{\infty} e^{-ax^2}I_{\nu}(bx})dx=\frac{\pi}{2\sqrt(a)}exp\left(\frac{b^2}{8a}\right)I_{\frac{1}{2}\nu}\left(\frac{b^2}{8a}\right)<br />

But how to do that?