エルミート多項式

エルミート多項式

エルミート多項式（Hermite polynomial）は $g(t,x)=\exp(2tx-t^2)$ を指数型母関数として生成される直交多項式です。すなわち、
 \[g(t,x)=\sum_{n=0}^{\infty}H_n(x)\frac{t^n}{n!}\]
と展開したときの係数をエルミート多項式 $H_n(x)$ と定義します。エルミート多項式は以下の漸化式を満たします。
 \[\begin{align*}&H_{n+1}(x)=2xH_n(x)-2nH_{n-1}(x)\\[6pt]&H’_n(x)=2nH_{n-1}(x)\\[6pt]&H’_n(x)=2xH_n(x)-H_{n+1}(x)\end{align*}\]
この漸化式を使って $H_n(x)$ の具体的な表式を求めると
 \[\begin{align*}H_0(x)&=1\\[6pt]H_1(x)&=2x\\[6pt]H_2(x)&=4x^2-2\\[6pt]H_3(x)&=8x^3-12x\\[6pt]H_4(x)&=16x^4-48x^2+12\end{align*}\]
のようになります。

scipy.special.eval_hermite()

scipy.special.eval_hermite() は任意の点におけるエルミート多項式の値を返します。

scipy.special.eval_hermite(n, x, out=None)

n には多項式の次数、x には配列もしくは配列に相当するオブジェクトを渡します。

# SCIPY_EVAL_HERMITE

# In[1]

import numpy as np
from scipy.special import eval_hermite

# 3次のエルミート多項式のx=0, 1, 2, 3, 4における値を計算
x = np.arange(5)
h = eval_hermite(3, x)

print(h)
[ -0.  -4.  40. 180. 464.]

$n=0,\ 1,\ 2,\ 3$ のエルミート多項式のグラフをプロットしてみます。

# In[2]

import matplotlib.pyplot as plt

# x座標データ
x = np.linspace(-4, 4, 129)

# グラフ描画領域を設定
fig = plt.figure(figsize = (6, 6))
ax = fig.add_subplot(111)
ax.set_title("Hermite polynomials", size = 15)
ax.grid()
ax.set_xlim(-4, 4)
ax.set_ylim(-10,10)
ax.set_xlabel("x", size = 14, labelpad = 8)
ax.set_ylabel("Hn(x)", size = 14, labelpad = 8)

# 色のリスト
c = ["red", "blue", "green", "darkorange"]

# H0,H1,H2,H3のグラフを描画
for i in range(4):
    ax.plot(x, eval_hermite(i, x),
            label = "n = {}".format(i), color = c[i-1])

# 凡例を表示
ax.legend()

plt.show()

scipy.special.hermite()

scipy.special.hermite() は、次数 n を指定して直交多項式オブジェクトを生成する関数です。

scipy.special.hermite(n, monic=False)

4 次のエルミート多項式を作ってみます。

# SCIPY_HERMITE

# In[1]

import numpy as np
from scipy.special import hermite

# 4次のエルミート多項式
h4 = hermite(4)

print(h4)
#    4      2
# 16 x - 48 x - 8.882e-16 x + 12

生成されたインスタンスの正式なクラス名を確認しておきます。

# In[2]

# h4のクラスを表示
print(type(h4))

# <class 'scipy.special.orthogonal.orthopoly1d'>

生成したオブジェクトから、任意の点における値を取得することができます。

# In[3]

# x=[0 1 2 3 4 5]
x = np.arange(0, 6)

# xの各要素におけるh4の値を表示
print(h4(x))
# [12.  -20.   76.  876. 3340. 8812.]

各項の係数はインスタンス変数 c に格納されています。

# In[4]

# h4の各項の係数を表示
print(h4.c)
# [1.6000000e+01 0.0000000e+00 -4.8000000e+01 -8.8817842e-16
#  1.2000000e+01]

多項式の根 (方程式 h4=0 の解) を得ることもできます。

# In[5]

# h4の根を表示
print(h4.r)

# [-1.65068012  1.65068012 -0.52464762  0.52464762]

正規化されたエルミート多項式

エルミート多項式に重み関数 $e^{-x^2/2}$ を乗じて、$\sqrt{\sqrt{\pi}2^nn!}$ で割った関数
 \[u_n(x)=\frac{1}{\sqrt{\sqrt{\pi}2^nn!}}H_n(x)\exp\left(-\frac{x^2}{2}\right)\]
は正規直交系をなします。すなわち $u_n$ を全区間で積分すると
 \[\int_{-\infty}^{\infty}u_m(x)u_n(x)=\delta_{mn}\]
となります ($\delta_{mn}$ は、$m=n$ のときに $1$、それ以外は $0$ の値をとるクロネッカーのデルタ記号です)。$u_n(x)$ は量子力学に登場する１次元調和振動子の固有関数です。

scipy.special.hermitenorm()

scipy.special.hermitenorm() は正規化されたエルミート多項式 $u_n$ を返します。

scipy.special.hermitenorm(n, monic=0)

5 次の正規化されたエルミート多項式を生成してみます。

# SCIPY_HERMITE_NORM

# In[1]

import numpy as np
from scipy.special import hermitenorm

# 正規化された5次のエルミート多項式
u5 = hermitenorm(5)

print(u5)
#    5      3
# 1 x - 10 x + 15 x

u5 の次数、係数、根を調べてみます。

# In[2]

# u5の次数
print("u5の次数\n{}{}".format(u5.order))

# u5の各項の係数
print("\nu5の係数\n{}".format(u5.c))

# u5=0の根
print("\nu5の根\n{}".format(u5.r))

# u5の次数
# 5

# u5の係数
# [1.   0. -10.   0.  15.   0.]

# u5の根
# [-2.85697001  2.85697001 -1.35562618  1.35562618  0.]

正規化されたエルミート多項式のグラフを描いてみます。

# In[3]

import matplotlib.pyplot as plt

# x座標データ
x = np.linspace(-4, 4, 129)

# グラフ描画領域を設定
fig = plt.figure(figsize = (6, 6))
ax = fig.add_subplot(111)
ax.set_title("Normalized Hermite Polynomial", size = 14)
ax.grid()
ax.set_xlim(-4, 4)
ax.set_ylim(-10,10)
ax.set_xlabel("x", size = 14, labelpad = 8)
ax.set_ylabel("un(x)", size = 14, labelpad = 8)

# 色のリスト
c = ["red", "blue", "green", "darkorange"]

# u0,u1,u2,u3のグラフを描画
for i in range(4):
    u = hermitenorm(i)
    ax.plot(x, u(x),
            label = "n = {}".format(i), color = c[i-1])

# 凡例を表示
ax.legend()

plt.show()