チェビシェフ多項式

当サイトではアフィリエイトプログラムを利用して商品を紹介しています。

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

第一種チェビシェフ多項式(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

コメント

  1. あとりえこばと より:

    【AI連載小説】科学とコードの交差点(45)「Pythonでチェビシェフ多項式を実装する方法」
     
    自習室でのディスカッションが進む中、開誠、明信、美純はPythonでのチェビシェフ多項式の実装について議論していました。ホワイトボードにはPythonコードの断片や数学の公式が書き込まれ、三人は最良の実装方法を見つけ出すべく取り組んでいました。

    # チェビシェフ多項式の再帰的な実装
    def chebyshev_recursive(n, x):
        if n == 0:
            return 1
        elif n == 1:
            return x
        else:
            return 2 * x * chebyshev_recursive(n - 1, x) - chebyshev_recursive(n - 2, x)
    
    # チェビシェフ多項式の反復的な実装
    def chebyshev_iterative(n, x):
        prev, current = 1, x
        for _ in range(2, n + 1):
            prev, current = current, 2 * x * current - prev
        return current
    
    # ディスカッションの中での利点を組み合わせた実装
    def optimized_chebyshev(n, x):
        if n == 0:
            return 1
        elif n == 1:
            return x
        else:
            prev, current = 1, x
            for _ in range(2, n + 1):
                prev, current = current, 2 * x * current - prev
            return current

    開誠:「再帰的な実装は数学的な概念に忠実で見通しがいいが、計算量の面で効率が心配だね」
    美純:「確かに、反復的なアプローチは計算が迅速だけど、数学的な表現が犠牲になりがちだわ」
    明信:「そこで、両方の長所を取り入れた最適化された実装も必要かもしれないな」
    ディスカッションは実際のコードに即して進み、最終的にはオプションに富んだ実装方法が生まれました。