# A closed expression involving Hermite polynomials

I stumbled upon the following integral:

$$\int e^{-\left(\frac{\gamma}{2} + ip \right)^2} H_n \left( q-\frac{\gamma}{2} \right) H_n \left( q+\frac{\gamma}{2} \right) \, d\gamma \, ,$$

and I really need a closed form of it. I'm pretty sure it exists, and probably in terms of $\Gamma$ functions, since

$$\int e^\frac{-\gamma^2}{2} H_n(\gamma) d\gamma = 4^n \sqrt{2} \,\Gamma \left( n + \frac{1}{2} \right) \, .$$

The problem is that Mathematica can't compute this in closed form, although it can do it term by term and hint that a closed expression probably exists:

Table[Integrate[
HermiteH[i, q + 1/2 \[Gamma]] HermiteH[i,
q - 1/2 \[Gamma]] Exp[-(\[Gamma]/2 + I p)^2], {\[Gamma], -Infinity,
Infinity}],{i,0,2}]


gives

$$\left\{2 \sqrt{\pi },4 \sqrt{\pi } \left(2 p^2+2 q^2-1\right),16 \sqrt{\pi } \left(2 \left(p^4+2 p^2 \left(q^2-1\right)+q^4\right)-4 q^2+1\right)\right\} \, .$$

Can someone give me a tip?

• Use a generating function for Hermit polynomials $\exp(2xt-t^2)=\sum_n H_n(x) t^n/n!$ – yarchik Sep 5 '16 at 21:50
• How would that help, since I don't have a summation? – QuantumBrick Sep 5 '16 at 21:57
• But you know how to extract each term using series. – yarchik Sep 5 '16 at 21:58
• My suggestion is not to use Mathematica, but to consult a table of integrals like Abramowitz and Stegun. – QuantumDot Sep 5 '16 at 22:00
• Can we make any assumptions about the parameters - Real, Integer, Positive...? – mikado Sep 5 '16 at 22:04

This is an integral I happen to have worked with before. The answer is a Laguerre polynomial, as can be verified like this for the first three indices as in the question:

ll =
2 Table[LaguerreL[i, 2 (p^2 + q^2)] (-2)^i Sqrt[Pi] i!, {i, 0, 2}];

l1 = Table[
Integrate[
HermiteH[i, q + 1/2 γ] HermiteH[i, q - 1/2 γ] Exp[-(γ/2 + I p)^2],
{γ, -Infinity, Infinity}], {i, 0, 2}];

Simplify[ll/l1]

(* ==> {1, 1, 1} *)


Here, ll is my closed form, and l1 is the list copied from the question. I verify that they give the same answers by forming the ratios.

The formula originally comes from Gradshtein and Ryzhik, equation 7.377, and I used it in my very first paper, Conductance quantization and backscattering below equation (13) to expand a shifted harmonic oscillator function in terms of the unshifted one. That's why I remembered where to look this up...

• I can't thank you enough. Indeed, I just remembered the seminal work by Groenewald, "On the principles of quantum mechanics", uses thais formula to calculate the Wigner function for a harmonic oscillator. Thank you again! – QuantumBrick Sep 6 '16 at 11:13

It is interesting to solve this integral without relying on Gradshteyn and Ryzhik. As I suggested in the comments, one can use the generating function approach. The solution consists of 7 steps. Human input is needed on step 5. All other steps are pretty much automatic.

1. Substitute generating functions $\exp(2xt-t^2)=\sum_n\frac{1}{n!}H_n(x)t^n$ and integrate

g[1] = Integrate[Exp[-(Y/2 + I p)^2] Exp[2(q-Y/2)t-t^2] Exp[2(q+Y/2)s-s^2],
{Y, -Infinity, Infinity}]


We have for g[1]=$2 \sqrt{\pi } e^{-2 i p (s-t)+2 q (s+t)-2 s t}=\sum_{mn}\frac{1}{m!n!}I_{m,n}(p,q)t^ms^n$, where we are interested in $I_{n,n}(p,q)$.

2. Next we need to extract from this generating function coefficients in front of $t^ns^n$, i.e., diagonal elements of the expansion. Let us introduce the substitution {t -> x/y, s -> y} and search for coefficients in front of $x^n y^0$.

g[2] = g[1] /. {t -> x/y, s -> y};

3. It is the exponential function with the $2\sqrt\pi$ prefactor and power:

pow = Cases[g[2], Exp[x_] -> x] // First;


i.e., g[2]=2 Sqrt[Pi] Exp[pow]

4. We simplify the argument of the exponent little bit:

g[3] = Collect[pow, {y, x}]


$2(q+ip)\,x/y+2(q-ip)\,y -2 x$

5. Now, we can easily extract the coefficient in front of $y^0$ in the series expansion of g[2] (this is the only nontrivial step in this derivation). For the summation, we rely on the amazing capabilities of Mathematica:

g[4] = 2 Sqrt[Pi] E^(-2 x) Sum[1/(2 n)! Binomial[2n,n] (2 I p + 2 q)^n x^n
(-2 I p + 2 q)^n,{n, 0, Infinity}]
(*2 E^(-2 x) Sqrt[Pi] Hypergeometric0F1[1, 4 (p^2 + q^2) x]*)

6. At this point we have the generating function, i.e., $2\sqrt{\pi} e^{-2x}{\,}_0F_1\left[;1,(p^2 + q^2) x\right]=\sum_n \frac{1}{(n!)^2}I_{n,n}x^n$. Finally, we extract corresponding series coefficients.

g[5] = (n!)^2 SeriesCoefficient[g[4], {x, 0, n}];

7. The solution in terms of the generalized hypergeometric function can be simplified as follows:

FullSimplify[FunctionExpand[g[5]], Assumptions -> n >= 0 && n \[Element] Integers]


And we are done:

g[6]=(-1)^n 2^(1 + n) n! Sqrt[Pi] LaguerreL[n, 2 (p^2 + q^2)]

• Actually you're still using Gradshteyn & Ryzhik because the generating function from which you start is equation GR 8.957... – Jens Sep 9 '16 at 4:05
• @Jens GR is such a comprehensive book that I am not surprised that the generating function of Hermite polynomials can be found there. I this case, I vaguely remembered that the function has an exponential form which is useful here. To be sure I just check it on wiki. But thank you to putting a reference to this equation in GR, this makes the answer more complete. – yarchik Sep 9 '16 at 6:39