Calculations with Hermite Polynomials

Questions
, ,
1

I want to evaluate moments of the Harmonic oscillator wave functions. I try the command

Integrate[x^p Exp[-x^2/2] HermiteH[n,x],{x,-Infinity,Infinity}, Assumptions->{Element[p,Integers,Element[n,Integers],n>=0,p>=0]

Mathematica replies that the integral doesn’t converge, which is clearly false.

Mathematica does correctly calculate the intgral if n and p are specified:

Integrate[x^2 Exp[-x^2/2] HermiteH[4,x],{x,-Infinity,Infinity}]
Sqrt(2 pi) 108

A related problem:
I want to express expressions like
x^p HermiteH[n,x]
as a sum of Hermite polynomials of different order. That is, I want to find the coefficients a[k] in the expansion
x^p HermiteH[n,x]=Sum[a[k] HermiteH[k,x],{k,0,n+p}].

This can be computed by repeated application of the identity:

x HermiteH[n,x]=HermiteH[n+1,x]/2 -n HermiteH[n-1,x]

How do I get Matematica to do this?

2

The result that the harmonic oscillator energies are quantized is often derived from the necessity of the wave function not diverging at infinity, so I’m skeptical of your statement that Mathematica’s statement about convergence is false.

3

Since HermiteH[n,x] is a polynomial of degree n, the coefficient of Exp[-x^2/2] in the integrand is a polynomial of degree n+p.
For large x, this behaves as x^(n+p) Exp[-x^2/2], which clearly goes to 0 faster than 1/x, so the integral converges..

4

As I remember, it’s of degree n if n is an integer.

5

In my Mathematica Integrate command I included the restriction that n and p both be positive integers.

6

It’s hard to tell if you correctly specified the restriction since you didn’t put your code in a code sample window.

7

Moments involve the square of the wave function. For example, the m^th moment is given by Integrate[x^m psi_n(x)^2, {x, -inf, inf}].

Many results involving Hermite polynomials can be found by making use of the generating function exp(2xt-t^2) = sum …

8

Mathematica can handle the integral for specific integer values of the first parameter n, but not for general n.

Integrate[
 x^p Exp[-x^2/2] HermiteH[3, x], {x, -Infinity, Infinity}, 
 Assumptions -> {p >= 0}]

(* Out[246]= -2^(2 + p/2) (-1 + (-1)^p) (1 + 2 p) Gamma[1 + p/2] *)

This has to do with the interplay of how convergence testing is done, and how assumptions are handled-- in effect integrality is replaced by realness so as to avoid nongeneric oversimplifications in case of discrete-valued parameters. For general real n>=0 the convergence test fails at infinity.

Even if the convergence were correctly assessed, I’m not sure a symbolic result would be produced. Below I use GenerateConditions in a way that, I think, bypasses the convergence assessment. It comes back unevaluated.
Integrate[
x^p Exp[-x^2/2] HermiteH[n, x], {x, -Infinity, Infinity},
Assumptions → {p >= 0, n >= 0}, GenerateConditions → False]

Out[247]= Integrate[
E^(-(x^2/2)) x^p HermiteH[n, x], {x, -[Infinity], [Infinity]},
Assumptions → {p >= 0, n >= 0}, GenerateConditions → False]
Integrate[
x^p Exp[-x^2/2] HermiteH[n, x], {x, -Infinity, Infinity},
Assumptions → {p >= 0, n >= 0}, GenerateConditions → False]

Out[247]= Integrate[
 E^(-(x^2/2)) x^p HermiteH[n, x], {x, -\[Infinity], \[Infinity]}, 
 Assumptions -> {p >= 0, n >= 0}, GenerateConditions -> False]
9

I’m not interested in calculating the quantum mechanical expectation value, which is the moment of the square of the wave function. I’m using the QM Harmonic oscillator wave functions as a basis set in which to expand a function I want to compute the moments of, so I want the moments of the HO wave functions, not its square.

10

Any suggestions for the second problem?

11

Are you looking for a symbolic Sum result for unspecified parameter values of p, or for code that does the decomposition into an explicit Plus expression when given specific integer values for p?

If the latter, could use a replacement rule, applied repeatedly.

rul = (x_^p_.* HermiteH[n_, x_] /; IntegerQ[p] && p >= 1) :> 
   x^(p - 1)*HermiteH[n + 1, x]/2 - n *x^(p - 1)*HermiteH[n - 1, x];

Example:

In[187]:= Expand[x^3*HermiteH[n, x] //. rul]

(* Out[187]= -2 n HermiteH[-3 + n, x] + 3 n^2 HermiteH[-3 + n, x] - 
 n^3 HermiteH[-3 + n, x] + 3/2 n^2 HermiteH[-1 + n, x] - 
 3/4 HermiteH[1 + n, x] - 3/4 n HermiteH[1 + n, x] + 
 1/8 HermiteH[3 + n, x] *)
12

Thank you for your help. I’d prefer a general formula for any n and p, but at least I can now compute specific values of p.

13

Here are a few things to consider.

(1) Actually your formula had a sign error (which i found by comparing HermiteH[4,x] behavior with HermiteH[n,x] and then substituting n->4).

So the corrected version of my replacement rule code would be

rul = (x_^p_.* HermiteH[n_, x_] /; IntegerQ[p] && p >= 1) :> 
   x^(p - 1)*HermiteH[n + 1, x]/2 + n *x^(p - 1)*HermiteH[n - 1, x];

(2) I do not know how to set up and solve an appropriate partial difference equation to obtain a general form. I suspect there is no closed form since it would involve a sum that likely has no closed form.

(3) An alternative to rule replacement is as follows (works for specific p but general n).

Clear[f]
f[0, n_] := HermiteH[n, x]
f[p_, n_] := f[p, n] = Expand[f[p - 1, n + 1]/2 + n*f[p - 1, n - 1]]

Example:

f[3, n]

(* Out[48]= 2 n HermiteH[-3 + n, x] - 3 n^2 HermiteH[-3 + n, x] + 
 n^3 HermiteH[-3 + n, x] + 3/2 n^2 HermiteH[-1 + n, x] + 
 3/4 HermiteH[1 + n, x] + 3/4 n HermiteH[1 + n, x] + 
 1/8 HermiteH[3 + n, x] *)