『Python数値計算ノート』ではアフィリエイトプログラムを利用して商品を紹介しています。

チェビシェフ多項式

第一種チェビシェフ多項式

第一種チェビシェフ多項式(Chebyshev polynomials of the first kind)
 \[T_n(\cos\theta)=\cos(n\theta)\tag{1}\]
によって定義されます。$x=\cos\theta$ とおけば、
 \[T_n(x)=\cos(n \mathrm{Arccos}x)\tag{2}\]
と表すことができます。逆三角関数 $\mathrm{Arccos}x$ が含まれるので、$T_n(x)$ の定義域は $[-1,1]$ に制限されます。第一種チェビシェフ多項式は漸化式
 \[\begin{align*}&T_0(x)=1\\[6pt]&T_1(x)=x\\[6pt]&T_{n+1}(x)=2xT_n(x)-T_{n-1}(x)\end{align*}\tag{3}\]
を満たします。$T_n(x)$ の具体的な表式を並べると
 \[\begin{align*}&T_0(x)=1\\[6pt]&T_1(x)=x\\[6pt]&T_2(x)=2x^2-1\\[6pt]&T_3(x)=4x^3-3x\\[6pt]&T_4(x)=8x^4-8x^2+1\end{align*}\]
となります。

scipy.special.eval_chebyt()

scipy.special.eval_chebyt(n, x) を使うと、n 次の第一種チェビシェフ多項式の点 x における値を評価できます。

# SCIPY_SPECIAL_EVAL_CHEBYT

# In[1]

import numpy as np
from scipy.special import eval_chebyt

x = np.array([-1, -0.5, 0, 0.5, 1])

# 9次の第一種チェビシェフ多項式のxにおける値を計算
t = eval_chebyt(9, x)

print(t)
# [-1.  1.  0. -1.  1.]

第一種チェビシェフ多項式のグラフを描画するコードも載せておきます。

# In[2]

# 第1種チェビシェフ多項式のグラフ

import matplotlib.pyplot as plt

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

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

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

# T0,T1,T2,T3のグラフを描画
for i in range(4):
    ax.plot(x, eval_chebyt(i, x),
            label = "n = {}".format(i), color = c[i-1])

# 凡例を表示
ax.legend()

plt.show()

第1種チェビシェフ多項式 (Chebyshev polynomials of the first kind)
scipy.special.eval_chebyt() は値をとるだけです。係数の取得や微分・積分など、多項式自体を操作する必要があるときは、後述する scipy.special.chebyt() で多項式オブジェクトを作成してください。

scipy.special.chebyt()

scipy.special.chebyt(n) は n 次の第一種チェビシェフ多項式オブジェクトを生成します。

# SCIPY_SPECIAL_CHEBYT

# In[1]

import numpy as np
from scipy.special import chebyt

# 3次の第一種チェビシェフ多項式を生成
t3 = chebyt(3)

print(t3)
#    3
4 x - 3 x

生成されたオブジェクトの引数に数値を渡して、任意の点におけるチェビシェフ多項式の値を取得できます。

# In[2]

# x=[-1.0 -0.5 0.0 0.5 1.0]
x = np.linspace(-1, 1, 5)

print(t3(x))
# [-1.  1. -0. -1.  1.]

多項式オブジェクトは、次数、係数、根、変数などのデータ属性(インスタンス変数)を備えています。

# In[3]

# 次数を取得
print("次数:{}".format(t3.o))

# 係数を取得
print("係数:{}".format(t3.c))

# 根を取得
print("根:{}".format(t3.r))

# 変数を取得
print("変数:{}".format(t3.variable))

# 次数:3
# 係数:[ 4.  0. -3.  0.]
# 根:[ 0.8660254 -0.8660254  0.       ]
# 変数:x

deriv()メソッドで多項式を微分できます。

# In[4]

# t3を微分
print(t3.deriv())
#     2
# 12 x - 3

 integ()メソッドで多項式を積分できます。

# In[5]

# t3を積分
print(t3.integ())
#    4       2
# 1 x - 1.5 x

第二種チェビシェフ多項式

第二種チェビシェフ多項式(Chebyshev polynomials of the second kind)は
 \[U_{n}(\cos t)=\frac{\sin(n+1)t}{\sin(t)}\tag{4}\]
で定義される多項式です。$x=\cos t$ とおくと、
 \[U_n(x)=\frac{\sin((n+1)\mathrm{Arccos}x)}{\sin(\mathrm{Arccos}x)}\tag{5}\]
となります。$U_n(x)$ は漸化式
 \[\begin{align*}&U_0(x)=1\\[6pt]&U_1(x)=2x\\[6pt]&U_{n+1}=2xU_n(x)-U_{n-1}(x)\end{align*}\tag{6}\]
を満たします。この漸化式だけでは区間は制限されませんが、定義式 (4) に逆三角関数 $\mathrm{Arccos}x$ が含まれているので、一般に $U_n(x)$ は区間 [$-1,\ 1$] に制限されます。漸化式を使うと、$U_n(x)$ の表式
 \[\begin{align*}&U_0(x)=1\\[6pt]&U_1(x)=x\\[6pt]&U_2(x)=4x^2-1\\[6pt]&U_3(x)=8x^3-4x\end{align*}\]
を順次得ることができます。

scipy.special.eval_chebyu()

scipy.special.eval_chebyu(n, x) は次数 n の第二種チェビシェフ多項式の x における値を評価します。ただし、n に非整数の値を渡すことも可能で、引数 x に区間制限は設けられていません。

# SCIPY_SPECIAL_EVAL_CHEBYU

# In[1]

import numpy as np
import matplotlib.pyplot as plt
from scipy.special import eval_chebyu

# グラフ描画領域を設定
fig = plt.figure(figsize = (6, 6))
ax = fig.add_subplot(111)
ax.set_title("Chebyshev polynomials of the 2nd kind", size = 14)
ax.grid()
ax.set_xlim(-1, 1)
ax.set_ylim(-3, 3)
ax.set_xlabel("x", size = 14, labelpad = 8)
ax.set_ylabel("Un(x)", size = 14, labelpad = 8)

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

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

# U0,U1,U2,U3のグラフを描画
for i in range(4):
    ax.plot(x, eval_chebyu(i, x),
            label = "n = {}".format(i), color = c[i-1])

# 凡例を表示
ax.legend()

plt.show()

第二種チェビシェフ多項式 (Chebyshev polynomials of the second kind) 2
第二種チェビシェフ多項式を (4) で定義する限り、$n$ は必ずしも整数値である必要はありません。ただし、$n$ が非整数の場合は多項式の直交性は失われ、$U_n(x)$ を漸化式を使って導き出すことはできません。
 
先ほど述べたように、scipy.special.eval_chebyu() は n に非整数の値を渡すこともできます。コードは省略しますが、n に半整数を渡したときの $U_n(x)$ のグラフも掲載しておきます。
 
Python 第二種チェビシェフ多項式 (Chebyshev polynomials of the second kind)

scipy.special.chebyu()

scipy.special.chebyu(n) は n 次の第二種チェビシェフ多項式オブジェクトを生成します。

# SCIPY_SPECIAL_CHEBYU

# In[1]

# 5次の第二種チェビシェフ多項式を生成
u5 = chebyu(5)

print(u5)
#     5      3
# 32 x - 32 x + 6 x

生成された第二種チェビシェフ多項式オブジェクトは次数、係数、根などの属性値をもっています。

# In[2]

# 次数を取得
print("次数:{}".format(u5.o))

# 係数を取得
print("係数:{}".format(u5.c))

# 根を取得
print("根:{}".format(u5.r))

# 変数を取得
print("変数:{}".format(u5.variable))

# 次数:5
# 係数:[ 32.   0. -32.   0.   6.   0.]
# 根:[-0.8660254 -0.5        0.8660254  0.5        0.       ]
# 変数:x

deriv()メソッドで微分、integ()メソッドで積分できます。

# In[3]

# u5を微分
u5_d = u5.deriv()

# u5を積分
u5_i = u5.integ()

print("導関数\n{}\n".format(u5_d))
print("積分\n{}".format(u5_i))

# 導関数
#      4      2
# 160 x - 96 x + 6

# 積分
#        6     4     2
# 5.333 x - 8 x + 3 x

コメント