Integration of products of spherical Hankel functions yield different results

Question

Background

I am trying to solve the following integral $J_{nm}^{(j)}(x) =\int_1^\infty y^{(-j)}h_n^{(1)}(x y)j_m(x^* y)\mathrm{d}y$

where $x$ is complex valued, the asterisk denotes the complex conjugate, $h_n^{(1)}$ is the spherical Hankel function of the first kind, and $j_m$ is the spherical Bessel function of the first kind.

In the end, I will need numerical values for the integrals but in order to find an appropriate quadrature I was playing around with these integrals in Mathematica.

For certain values for $n$, $m$, $j$, and $x$ Mathematica will be able to evaluate the integrals analytically using the Mathematica's Integrate-function. The problem is that if I use the numerical quadrature NIntegrate the result will be very different. This might not be surprising, since Bessel functions are known to be notoriously difficult to numerically integrate. However, the integration also yields incorrect results on very small subintervals where integrand of $H_{nm}^{(j)}$ appears to be linear. So the reason does not seem to be the usual problems with the oscillations of the integrand.

An example:

The problem appears in particular if the imaginary part of $x$ is small.

Integrating on the subinterval $[1, 1.1]$ for $x = (1 + i) \times 10^{-8}, n = 0, m = 1, j = 0$

yields for the analytic integration $\int_1^{1.1} y^{(-j)}h_n^{(1)}(x y)j_m(x^* y)\mathrm{d}y = 6.42699\times 10^7 + 6.42699\times10^7 i$ while the numerical integration yields $\int_1^{1.1} y^{(-j)}h_n^{(1)}(x y)j_m(x^* y)\mathrm{d}y = -0.0333333 - 3.5\times10^{-10} i.$

Looking at the graph of the kernel, the result from the numerical quadrature looks much more reasonable.

Question

I looked the problem up in the literature and there are many articles discussing these types of integrals. However, the specifics never quite match my case. Now I have some questions:

• Where is the difference coming from and which result can I trust? Why?
• Is my assumption correct that integration of a complex integrand along the real axis in this case means that I do not have to deal with complex analysis (calculus of residues)?
• Why does the problem arise only for certain values of $j$, $n$, $m$, and $x$. In particular, why does it appear especially if $\mathrm{Im(x)} \ll 1$?

Minimal Working Example

Here is the code:

J[y_, x_, n_, m_, j_] := y^(-j) * SphericalHankelH1[n, x * y] * SphericalBesselJ[m, Conjugate[x] * y];
Integrate[J[z, (1+I) * 10^(-8), 0, 1, 0], {z, 1, 1.1}]
NIntegrate[J[z, (1+I) * 10^(-8), 0, 1, 0], {z, 1, 1.1}]

Thank you!

Answer 1

Integrate is a symbolic solver and sometimes messes up numerically when inexact coefficients are used.

Using a high-enough precision to avoid round-off error:

Integrate[jj[z, (1 + I)*10^(-8), 0, 1, 0], {z, 1, 1.1`50}]
(*  -0.0333333329833333333333333 - 3.499999970577778*10^-10 I  *)

Using N[] on an exact result; note that machine-precision is insufficient to evaluate the result accurately:

res = Integrate[jj[z, (1 + I)*10^(-8), 0, 1, 0], {z, 1, 11/10}];
N[res]
N[res, 16]
(*
  1.32231*10^7 + 1.32231*10^7 I
  -0.03333333298333333 - 3.5000000*10^-10 I
*)