最近在做一个数值逼近的算法,里面用到了埃尔米特多项式。所以就花了些时间推导了一遍,推导笔记放在这里算是给自己做个备忘。
埃尔米特多项式 (Hermite Polynomials)简介
埃尔米特多项式是一组正交的多项式。就如许多其他的以人名命名的数学公式一样,埃尔米特多项式其实也并不是埃尔米特第一个提出的。 Laplace 在 1810 年一篇论文中就给出了埃尔米特多项式的系数,Chebyshev 则在 1859 年的一篇论文中详细的讨论了埃尔米特多项式的各种性质。可惜 Chebyshev 的这篇论文并没有引起学术圈应由的重视。 Charles Hermite 在 1864 年的一篇文章中才提到埃尔米特多项式,这已经比 Laplace 最初的研究成果晚了 54 年。
由于物理学家的小圈子与数学家的小圈子相对独立,现在有两种定义略有不同的埃尔米特多项式。一种被称为“统计学家的埃尔米特多项式”,另一种的被称为“物理学家的埃尔米特多项式”。这两种埃尔米特多项式的系数是有联系的,这里我只介绍“物理学家的埃尔米特多项式”。
前 10 个 埃尔米特多项式如下:
H0(x)H1(x)H2(x)H3(x)H4(x)H5(x)H6(x)H7(x)H8(x)H9(x)H10(x)===========12x4x2−28x3−12x16x4−48x2+1232x5−160x3+120x64x6−480x4+720x2−120128x7−1344x5+3360x3−1680x256x8−3584x6+13440x4−13440x2+1680512x9−9216x7+48384x5−80640x3+30240x1024x10−23040x8+161280x6−403200x4+302400x2−30240H0(x)=1H1(x)=2xH2(x)=4x2−2H3(x)=8x3−12xH4(x)=16x4−48x2+12H5(x)=32x5−160x3+120xH6(x)=64x6−480x4+720x2−120H7(x)=128x7−1344x5+3360x3−1680xH8(x)=256x8−3584x6+13440x4−13440x2+1680H9(x)=512x9−9216x7+48384x5−80640x3+30240xH10(x)=1024x10−23040x8+161280x6−403200x4+302400x2−30240
由微分方程引入埃尔米特多项式
常微分方程:
y′′−2xy′+2ny=0, (n∈N+)y″−2xy′+2ny=0, (n∈N+)
被称为 埃尔米特 方程。
x=0x=0 是埃尔米特方程的常点,所以这个微分方程的解可以在 x=0x=0 的邻域表示为泰勒级数:
y=∑k=0∞akxky=∑k=0∞akxk
对这个级数求导可以得到:
y′=∑k=1∞kakxk−1y′′=∑k=2∞k(k−1)akxk−2y′=∑k=1∞kakxk−1y″=∑k=2∞k(k−1)akxk−2
带入埃尔米特方程得到:
∑k=2∞k(k−1)akxk−2−2∑k=1∞kakxk+2n∑k=0∞akxk=0∑k=2∞k(k−1)akxk−2−2∑k=1∞kakxk+2n∑k=0∞akxk=0
考察 x0x0 的系数,有:
2a2x0+2na0x0=02a2x0+2na0x0=0
所以:
a2=−na0a2=−na0
考察 x1x1 的系数,有:
3⋅2⋅a3x1−2a1x1+2na1x1=03⋅2⋅a3x1−2a1x1+2na1x1=0
所以
a3=2−2n6a1a3=2−2n6a1
考察 xmxm 的系数,有:
(m+2)(m+1)am+2xm−2mamxm+2namxm=0(m+2)(m+1)am+2xm−2mamxm+2namxm=0
所以:
am+2=2m−2n(m+2)(m+1)amam+2=2m−2n(m+2)(m+1)am
不难看出,对于偶次项,有:
a2m=∏m−1p=0(4p−2n)(2m)!a0a2m=∏p=0m−1(4p−2n)(2m)!a0
对于奇次项,有:
a2m+1=∏mp=1(4p−2−2n)(2m+1)!a1a2m+1=∏p=1m(4p−2−2n)(2m+1)!a1
可以看出奇数和偶数是两套独立的递推系数,因此可以这样写:
y(x)=a0y0(x)+a1y1(x)y(x)=a0y0(x)+a1y1(x)
y0(x)y0(x) 只含有 xx 的偶次项,y1(x)y1(x)只含有 xx 的奇次项。
当 nn 是偶数时,y0(x)y0(x) 只有有限项,只要我们将 a1a1 设为 00,那么 y(x)y(x) 就退化为多项式了。
同理,当 nn 是奇数时,y1(x)y1(x) 只有有限项,只要我们将 a0a0 设为 00,那么 y(x)y(x) 就退化为多项式了。
按照这个思路,选择合适的 a0a0 和 a1a1 ,可以使得多项式的最高次为 (2x)n(2x)n:
当 n=0n=0 时, y(0)=a0=1y(0)=a0=1
当 n=1n=1 时, y(0)=a1x=2xy(0)=a1x=2x
当 n=2n=2 时, y(0)=a0(1−2x2)=4x2−2y(0)=a0(1−2x2)=4x2−2
当 n=3n=3 时, y(0)=a1(x−23x3)=8x3−12xy(0)=a1(x−23x3)=8x3−12x
当 n=4n=4 时, y(0)=a0(1−4x2+43x4)=16x4−48x2+12y(0)=a0(1−4x2+43x4)=16x4−48x2+12
这正是我们的埃尔米特多项式。
埃尔米特多项式的母函数和递推公式
上面的方式虽然推导出了埃尔米特多项式,但是用起来并不方便。实际上,函数 Ψ(t,x)=exp(2tx−t2)Ψ(t,x)=exp(2tx−t2) 在 t=0t=0 处的泰勒展开可以写为:
Ψ(t,x)=exp(2tx−t2)=∑n=0∞Hn(x)tnn!Ψ(t,x)=exp(2tx−t2)=∑n=0∞Hn(x)tnn!
其中的系数 Hn(x)Hn(x) 就是埃尔米特多项式。下面我们就来证明这个结论。
首先,容易验证:
∂Ψ(t,x)∂x=2tΨ(t,x)∂Ψ(t,x)∂t+2(t−x)Ψ(t,x)=0∂Ψ(t,x)∂x=2tΨ(t,x)∂Ψ(t,x)∂t+2(t−x)Ψ(t,x)=0
把 Ψ(t,x)Ψ(t,x) 的泰勒展开带入上面的两个式子:
∑n=0∞H′n(x)tnn!=∑n=0∞2Hn(x)tn+1n!∑n=0∞Hn(x)tn−1(n−1)!+2(t−x)∑n=0∞Hn(x)tnn!=0∑n=0∞Hn′(x)tnn!=∑n=0∞2Hn(x)tn+1n!∑n=0∞Hn(x)tn−1(n−1)!+2(t−x)∑n=0∞Hn(x)tnn!=0
比较两边 tt 的同幂项,得到:
H′0(x)=0H′n(x)=2nHn−1(x)Hn+1(x)+2nHn−1(x)−2xHn(x)=0H0′(x)=0Hn′(x)=2nHn−1(x)Hn+1(x)+2nHn−1(x)−2xHn(x)=0
所以 :
H0(x)=c0H1(x)=2x+c1⋯H0(x)=c0H1(x)=2x+c1⋯
可知, Hn(x)Hn(x) 确实是多项式。继续推导可知:
Hn(x)−2xHn−1(x)+2(n−1)Hn−2(x)=0Hn(x)−xnH′n(x)+H′n−1(x)=0Hn(x)−xnH′n(x)+12nH′′n−1(x)=0H′′n−1(x)−2xH′n(x)+2nHn(x)=0Hn(x)−2xHn−1(x)+2(n−1)Hn−2(x)=0Hn(x)−xnHn′(x)+Hn−1′(x)=0Hn(x)−xnHn′(x)+12nHn−1″(x)=0Hn−1″(x)−2xHn′(x)+2nHn(x)=0
最后这个式子就是埃尔米特方程,而埃尔米特方程的多项式解只能是埃尔米特多项式。
上面的推导过程中同时还给出了埃尔米特多项式的递推公式:
Hn(x)−2xHn−1(x)+2(n−1)Hn−2(x)=0Hn(x)−2xHn−1(x)+2(n−1)Hn−2(x)=0
有了这个递推公式,再计算埃尔米特多项式就容易多了。