I try to solve $\int_0^T e^{-x}\frac{I_n(\alpha x)}{x}dx$ where $I_n(x)$ is the modified Bessel function of the first kind and $0<\alpha<1$.

My first approach was to turn this integral into an infinite sum to fit a hypergeometric series:

- Using the infinite series representation of the Bessel function, I got incomplete gamma functions in the sum, which does not sound promising.
- The [multiplication theorem][1] yields an infinite series of integrals where we get rid of the $\alpha$:
$$\int_0^T e^{-x}\frac{I_n(\alpha x)}{x}dx=\alpha^n\sum_{m=0}^{\infty}\frac {\big(\frac {\alpha ^{2}-1}{2}\big)^m}{m!}\int_0^T e^{-x}x^{m-1}I_{n+m}(x)dx$$

According to a table of integrals, the new integrals are:
\begin{align}
\int_0^T e^{-x}x^{m-1}I_{n+m}(x)dx&=\frac{T^{2m+n}}{2^{m+n}}\frac{\Gamma(2m+n)}{\Gamma(m+n+1)\Gamma(2m+n+1)}\\
&\times{}_2F_2[\{m+n+\frac{1}{2},2m+n\};\{2m+2n+1,2m+n+1\};-2T]
\end{align}

Expanding ${}_2F_2$ (let's call $k$ the summation index), we get a double infinite series which might fit the definition of a hypergeometric function of 2 variables. However, I get several Pochhammer symbols with coupled summations indices like this: $\frac{(n)_{2m+k}}{(n+1)_{2m+k}(2n)_{2m+k}}$ which does not fit any hypergeometric function definition.

Another approach could be to get inspiration from the limit $T\rightarrow \infty$ which is the Laplace transform of $\frac{I_n(x)}{x}$ (up to a constant) and it has a closed form (according to a table):
$$\int_0^\infty e^{-x}\frac{I_n(\alpha x)}{x}dx=\frac{\big(\frac{\alpha}{1+\sqrt{1-\alpha^2}}\big)^n}{n}$$

However, I don't find any reference on the way to compute this. 

Do you have any suggestion about how to compute it?


  [1]: https://en.wikipedia.org/wiki/Bessel_function#Multiplication_theorem