# How to get details about how Mathematica did a definite Integral

Assuming[{Element[S, Reals],S>0},Integrate[Exp[-I*S*w]/(w^2 + 1)^(3/2)
,{w, 0, Infinity}]]


gets

 1/2 S(I π BesselI[1, S] + 2 BesselK[1, S] - I π StruveL[-1, S])


However Mathematica cannot do the indefinte integral with a variable lower or upper limit so differentiation cannot be used to see how it did this integral. Is there any way to find out or person who can figure it out, please. Change of variable to Sw would give leading S as a divisor not a factor. Not clear where the 3 parts come from.

Thanks JM

          LaplaceTransform[1/(w^2 + 1)^(3/2), w, s]
-(1/2) π s (BesselY[1, s] + StruveH[-1, s])


is down to 2 terms. But still unclear where they come from

• Have you tried looking at the result of LaplaceTransform[1/(w^2 + 1)^(3/2), w, s]? Commented Mar 27, 2019 at 13:57
• Related: mathematica.stackexchange.com/a/89190/4999 (see "Hints" section). Commented Mar 28, 2019 at 0:20
• This might be helpful: mathematica.stackexchange.com/questions/148/… Commented Mar 29, 2019 at 4:03
• @theorist, since special functions are involved, none of the stuff in the thread you linked to would be of much use for now. Commented Mar 30, 2019 at 6:32

I'm not sure how Integrate found the answer, but here's a way it could have conceivably used.

If we can express your integrand as the product of two MeijerG functions, there's a rich list of integral identities we can use. In particular we can treat your integral as a Mellin transform at s == 1 and use the general identity here:

In our case:

mei = MeijerGReduce[Exp[-I*S*w]/(w^2 + 1)^(3/2), w]


int = MellinTransform[mei, w, 1]

MeijerG[{{1/2}, {}}, {{0, 1/2, 1}, {}}, -S^2/4]/π

FullSimplify[FunctionExpand[int], S > 0]

1/2 I S (-2 + π BesselI[1, S] - 2 I BesselK[1, S] - π StruveL[1, S])


We can see int (defined above) appears in the TracePrint of the Integrate call and so maybe this is the method it chose:

TracePrint[
Integrate[Exp[-I*S*w]/(w^2 + 1)^(3/2), {w, 0, Infinity}],
_MeijerG,
TraceInternal -> True
]

MeijerG[{{1/2}, {}}, {{0, 1/2, 1}, {}}, 16/S^4, -2]
...
...
long spew omitted
...

• Indeed, you can recast the integral as a Mellin convolution: Assuming[S > 0, MellinConvolve[Exp[-I S w]/w^2, (1 + w^2)^(-3/2), w, 1]] Commented Mar 27, 2019 at 16:40
• Neither MeijerGReduce or MellinConvolve copied and pasted seem to work in v7, but maybe has them internally to do the integral . i tried TracePrint on the original integral fortunately in a new notebook but it has spewed for 10 minutes at least. Had to abort with TaskManager! Commented Mar 28, 2019 at 2:51
• @simon, the Mellin stuff are new in version 11, so you can't use them in version 7 (and you made no mention that you were using an old version). Commented Mar 28, 2019 at 4:56
• The output is (2*Inactive[MeijerG][{{-1/2}, {}}, {{0}, {}}, w^2]*Inactive[MeijerG][{{}, {}}, {{0, 1/2}, {}}, I*S*w/2, 1/2])/Pi. The reason I posted it as an image is the presence of Inactive, which doesn't exist in V7. You can just pattern match it away, but be aware that these particular MeijerG expressions auto evaluate. Commented Mar 29, 2019 at 14:43
• That terms has a series expansion 2 π I Sum[(-S/2)^n/(n Gamma[n/2]^2), {n, 1, ∞}]. Maybe it's possible to get a more concise expression from that. Commented Mar 30, 2019 at 14:16

Thanks to everyone for their input. The real and imaginary parts of this integral are in fact given in Abramovitz and Stegun, eqns. 9.6.25 and 12.2.3 resp. sections on $$K_1$$ and $$L_1$$, though that still doesn't explain how they were worked out originally.

$$I_{-\nu}(x)-\mathbf{L}_\nu(x)=\frac{2(x/2)^\nu}{\sqrt{\pi}\,\Gamma(\nu+\tfrac12)}\int_0^\infty\sin(tx)(1+t^2)^{\nu-\tfrac12}\mathrm dt\quad \Re\nu<\frac12,\,x>0$$

There is a as much info there on the imaginary part difference function BesselI[1, S] - StruveL[-1, S]) as StruveL individually including the asymptotic series so their small difference has been noted before. But the Taylor series that Chip found is not there, so thanks.

Unfortunately Mathematica does not recoginise this difference function and does not use the differnce asymptotic series so its computation broke down at S=35 for me.