Integrating products of Hermite polynomials

Question

Given Hermite polynomials $\Phi_0 = 1, \Phi_1 = \xi, \Phi_2 = \xi^2 - 1, \ldots$, I want to calculate

$$ \int_{-\infty}^\infty W(\xi) \prod_{n \in N} \Phi_n(\xi) \:\mathrm{d} \xi $$ where $N$ is any list of polynomial degrees and the weighting function $W(\xi)$ is $$ W(\xi) = \frac{1}{\sqrt{2\pi}} e^{-\xi^2/2} $$ At the moment, I'm only interested in integrals of up to three terms. Is there a straightforward method for calculating such integrals exactly? And has the method been implemented in code already?

Wolfram Alpha is able to calculate such integrals exactly and, for those products I've tried, the results are always integers.

Answer 1

I realised that I can calculate the product exactly using Gauss-Hermite quadrature using $N+1$. This seems the most straightforward method to implement and works for any $N$.

Answer 2

Lets call the three Hermite polynomials $A,B,C=\Phi_k$. Then as the first $J+1$ Hermite polynomials form a basis of the polynomials of degree $J$, we can express

$$ AB = \sum_{j=0}^{\deg(AB)}a_j \Phi_j$$ By orthogonality, the answer is then $$\int_{\mathbb R} W(\xi) ABC \ d\xi = \langle AB, \Phi_k\rangle = a_k\|\Phi_k\|^2 $$ where the inner product is $$ \langle f,g\rangle = \int_\mathbb R f(\xi)g(\xi) W(\xi) d\xi$$ and the norm is the natural one induced by the inner product. Depending on your convention of what the Hermite polynomial is, i.e. what $W$ is, $\Phi$ may not have unit norm so you will have to check that. I don't know what the optimal choice of $C$ is, maybe you should take $C$ to be the largest degree polynomial? A nice special case is if $\deg C > \deg(AB)$ then the answer is automatically $0$.

Answer 3

I'm working with the physicist defn. of Hermite polynomials in the following. The most efficient method for dealing with $$ S_{k,m,n}:=\frac{1}{\sqrt{\pi}} \int_{-\infty}^\infty dx\, e^{-x^2}H_k(x)\,H_m(x)\,H_n(x)$$ is to combine the product formula $$ H_m(x)\,H_n(x)=2^n\,n! \sum_{j=0}^n\binom{m}{n-j}\frac{H_{m-n+2n}(x)}{2^jj!} \quad, \quad m \ge n$$ and the orthogonality relationship $$ \frac{1}{\sqrt{\pi}} \int_{-\infty}^\infty dx\, e^{-x^2} H_p(x)\,H_q(x) = 2^p\,p! \,\delta_{p,q}$$ where the Kronecker delta is 1 if $p=q$ and 0 otherwise. Your answer becomes, $$S_{k,m,n}=\frac{2^k\,k!}{2^j\,j!} 2^n\,n!\binom{m}{n-j} \quad, \quad j=(k-m+n)/2$$ with the caveats that the expression is 0 if $j$ is a half integer, and also $m\ge n.$