# A closed form expression for $\int_0^{\infty} (t^2+t^4)^n e^{-t^2-t^4}\,dt$

I was doing some computations for research purposes, which led me to this integral:

$$I(n) = \int_0^{\infty} (t^2+t^4)^n e^{-t^2-t^4}\,dt.$$

This is very suggestively written so as to employ a parametric differentiation technique as so:

$$\left(\frac{\partial^n}{\partial \alpha^n}\right)\int_0^{\infty}e^{-\alpha(t^2+t^4)}\,dt.$$

This integral has a nice, closed form expression:

$$\int_0^{\infty} e^{-\alpha(t^2+t^4)}\,dt = \frac{1}{4} e^{\frac{\alpha}{8}}K_{\frac{1}{4}}\left(\frac{\alpha}{8}\right),$$

where $K_{\nu}$ is the modified Bessel function of the second kind. From here, I would have to employ $n$ differentiations which would be pretty messy to work out due to product rule. $n$ applications of product rule does have a nice combinatorial expression but it is far from explicit. Moreover, Bessel functions can get pretty complicated after differentiating so this seems like a bad approach.

Instead I decided to run some examples on Mathematica and computed the first 22 of these and noticed a very surprising pattern. In what follows $I_{\nu}$ is the modified Bessel function of the first kind.

$$I(0) = \frac{1}{4} e^{\frac{1}{8}}K_{\frac{1}{4}}\left(\frac{1}{8}\right)$$

$$I(1) = \frac{1}{32} e^{\frac{1}{8}}\left(K_{\frac{1}{4}}\left(\frac{1}{8}\right) + K_{\frac{3}{4}}\left(\frac{1}{8}\right)\right)$$

$$I(2) = \frac{3}{128\sqrt{2}} e^{\frac{1}{8}}\pi \left(3 I_{-\frac{1}{4}} \left(\frac{1}{8}\right) + I_{\frac{1}{4}}\left(\frac{1}{8}\right) - I_{\frac{3}{4}}\left(\frac{1}{8}\right) + I_{\frac{5}{4}}\left(\frac{1}{8}\right)\right)$$

$$I(3) = \frac{1}{256\sqrt{2}} e^{\frac{1}{8}}\pi \left(39 I_{-\frac{1}{4}} \left(\frac{1}{8}\right) + 17 I_{\frac{1}{4}}\left(\frac{1}{8}\right) - 14 I_{\frac{3}{4}}\left(\frac{1}{8}\right) + 14 I_{\frac{5}{4}}\left(\frac{1}{8}\right)\right)$$

$$I(4) = \frac{1}{2048\sqrt{2}} e^{\frac{1}{8}} \pi \left(1029 I_{-\frac{1}{4}} \left(\frac{1}{8}\right) + 367 I_{\frac{1}{4}}\left(\frac{1}{8}\right) - 349 I_{\frac{3}{4}}\left(\frac{1}{8}\right) + 349 I_{\frac{5}{4}}\left(\frac{1}{8}\right)\right)$$

$$I(5) = \frac{9}{8192\sqrt{2}} e^{\frac{1}{8}} \pi \left(1953 I_{-\frac{1}{4}} \left(\frac{1}{8}\right) + 619 I_{\frac{1}{4}}\left(\frac{1}{8}\right) - 643 I_{\frac{3}{4}}\left(\frac{1}{8}\right) + 643 I_{\frac{5}{4}}\left(\frac{1}{8}\right)\right)$$

$$I(6) = \frac{1}{16384\sqrt{2}} e^{\frac{1}{8}} \pi \left(185157 I_{-\frac{1}{4}} \left(\frac{1}{8}\right) + 53131 I_{\frac{1}{4}}\left(\frac{1}{8}\right) - 59572 I_{\frac{3}{4}}\left(\frac{1}{8}\right) + 59572 I_{\frac{5}{4}}\left(\frac{1}{8}\right)\right)$$

Repeat ad nauseum. Each of the terms in the denominator seems to be a power of $2$, the third and fourth terms seem to have the same coefficient (modulo a sign) and the signs are $+$, $+$, $-$, $+$. The "nice" output seems to suggest to me that there is a closed-form expression for $I(n)$ in general but I haven't the slightest clue as to how to come up with it. Can anyone shed some light on the matter?

A PDF with more expressions can be found here. (Mathematica output.)

-
Can it be done by substituting $k=t^2+t^4$, partial integration and finally using the gamma function? –  REr Jul 13 '14 at 9:14
No that's a bad way to do it because of some horrible non linearities that crop up. I tried and quickly dropped that idea. –  Cameron Williams Jul 13 '14 at 9:56

Let us make the change of variables $t=\sinh \frac{x}{4}$. Since $$t^2+t^4=\frac{\sinh^2\frac{x}{2}}{4}=\frac{\cosh x-1}{8}$$ and $dt=\frac14 \cosh\frac{x}{4}dx$, the initial integral can be rewritten as \begin{align} I(n)&=\frac{e^{\frac18}}{4^{n+1}}\int_0^{\infty}\cosh\frac{x}{4}\sinh^{2n}\frac{x}{2}e^{-\frac18\cosh x}dx=\\ &=\frac{e^{\frac18}}{2^{4n+3}}\int_{-\infty}^{\infty}e^{\frac{x}{4}}\left(e^{\frac{x}{2}}-e^{-\frac{x}{2}}\right)^{2n}e^{-\frac18\cosh x}dx=\\ &=\frac{e^{\frac18}}{2^{4n+3}}\sum_{k=0}^{2n} \left(-1\right)^k {2n \choose k}\int_{-\infty}^{\infty}e^{(n-k+\frac{1}{4})x-\frac18\cosh x}dx=\\ &=\frac{e^{\frac18}}{2^{4n+2}}\sum_{k=0}^{2n} \left(-1\right)^k {2n \choose k} K_{n-k+\frac14}\left(\frac18\right). \end{align} At the last step, I use the integral representation $\displaystyle K_{\nu}(r)=\frac12\int_{-\infty}^{\infty}e^{-r\cosh x+\nu x}dx$.

Added: Using the recurrence relation $K_{\nu+1}(r)-K_{\nu-1}(r)=\frac{2\nu}{r}K_{\nu}(r)$, the last sum (1) can be expressed in terms of rational multiples of only two Macdonald functions, for example, $K_{\frac14}\left(\frac18\right)$ and $K_{\frac54}\left(\frac18\right)$. However it seems unlikely that the corresponding coefficients can be written in a nice way.

-
This was a very neat solution and the form looks more or less like the answer you would conceivably get from parametric differentiation, without all of the fuss that comes with it. Thanks for such an enlightening approach. It never would have occurred to me to use hyperbolic trig substitution. –  Cameron Williams Jul 13 '14 at 9:52
Very nice result! –  mike Jul 13 '14 at 12:58
It seems the way to simplify this is by making use of the identity $$K_{\nu}(z) = \frac{\pi}{2}\frac{I_{-\nu}(z)-I_{\nu}(z)}{\sin(\nu\pi)}.$$ –  Cameron Williams Jul 13 '14 at 19:16

Not a full answer, but a partial one:
Notice that the powers of two in the denominator are increasing in periodic repetitions of 3, 2, 1. That is, in this case the exponents on the 2 are 2, 5, 7, 8, 11, 13, 14, 17, 19, 20.... This is sequence A047268 on the OEIS website, defined as the sequence of numbers congruent to {1, 3, 5} mod(6).

I fiddled around with the other coefficients, however nothing useful turned out (no known sequences with the numbers themselves or their differences, products, sums...). If you absolutely need a closed form expression for it, then you could just define three new sequences $S_1$, $S_2$ and $S_3$ as the coefficients of the first three modified Bessel functions of the first kind in your expression. Then the closed form expression would be:

$I(n) = 2^{-\delta_{n + 2}}e^{\frac18}(K_{\frac14}(\frac18)+nK_{\frac34}(\frac18))$ $n = 0,1$

$I(n) = 2^{-\frac{\delta_{n + 2}}{2}}e^{\frac18}\pi(S_1I_{-\frac{1}{4}}(\frac18) + S_2I_{\frac{1}{4}}(\frac18) - S_3I_{\frac{3}{4}}(\frac18) + S_3I_{\frac{5}{4}}(\frac18))$ $n \geq 2$.

Where,
$\delta_n$ is sequence A047268,
$S_1, S_2, S_3$ are explicitly given by you. Note that the coefficients are not simplified by being factored out, as is for example in I(2) where 3 is factored out. You need to multiply them back in to get the sequences so there will always be a coefficient of 1 outside the numerator.

Finally, if you really want to get compact, you can write the second expression as such:
$I(n) = 2^{-\frac{\delta_{n + 2}}{2}}e^{\frac18}\pi\sum\limits_{k = 1}^4 S_kI_{\frac{2k - 3}{2}}(\frac18)$.

Where $\delta_n$ and $S_1, S_2, S_3$ are defined as above, and $S_4 = -S_3$ (it appears implicitly in the summation).

Hope this helps!

-

I would think that the parametric differentiation technique is probably easier to carry through:

$$\frac{\partial}{\partial r}K_{\nu}(r)= -(1/2)(K_{\nu-1}(r) + K_{\nu+1}(r))......(1)$$

EDIT: From the answer by O.L., we know that:

$$K_{\nu+1}(r)-K_{\nu-1}(r)=\frac{2\nu}{r}K_{\nu}(r)......(2)$$

Thus we can get rid of $K_{\nu-1}(r)$ and obtain: $$\frac{\partial}{\partial r}K_{\nu}(r)= -K_{\nu+1}(r) +(1/2)\frac{2\nu}{r}K_{\nu}(r)......(3)$$

Therefore, as O.L. pointed out, the final result can be expressed as a linear combination of $K_{1/4}(1/8)$ and $K_{5/4}(1/8)$.

-
The issue is that I don't have to do just one derivative but $n$ of them. It's going to be a nightmare of an expression when all is said and done when done this way. The relatively nice output seems to suggest that there's a better way to do it (instead of having $n+1$ terms like parametric differentiation would give). –  Cameron Williams Jul 13 '14 at 2:47