[SymPy] 級数展開

[SymPy] 級数展開

ベキ級数展開

関数の級数展開

 ある関数 $f(x)$ を $x$ のベキ級数(整級数)
 
\[a_0+a_1(x-x_0)^1+a_2(x-x_0)^2+\cdots\]
で表すことを 級数展開 といいます。無限に項をとることによって、理論上は関数と級数展開は厳密に一致しますが、実用のために有限項で打ち切っても、適切な数の項をとることによって、ベキ級数は $x=x_0$ の近くで関数の良い近似を与えます (ただし、$x$ が $x_0$ から離れるほどその精度は下がっていきます)。たとえば、$\sin x$ のマクローリン級数
 
\[\sin x=x-\frac{x^3}{3!}+\frac{x^5}{5!}-\frac{x^7}{7!}+\: \cdots\]
などは級数展開の代表的な例の1つです (マクローリン級数の詳細については、こちらのサイトを参照してください)。
 

sympy.series()

 sympy.series() を使うと任意の関数をベキ級数展開できます。
 以下のコードは Jupyter notebook の使用を前提に記述されています。
 リストA-1 で必要な関数をまとめてインポートしておきます。

# リストA-1
from sympy import init_printing, var, series, poly, cos, sin

# 数式をLaTeX表示
init_printing()

# シンボルxを定義
var('x')

 $\sin x$ を $x=0$ の近くで級数展開してみましょう。

# リストA-2

# sinxをx=0の周りで級数展開
series(sin(x), x)
\[x - \frac{x^{3}}{6} + \frac{x^{5}}{120} + O\left(x^{6}\right)\]

 デフォルト設定では、$x^6$ 以降の項はすべて O(x**6) にまとめられています。オプション引数 n でこの設定を変更することができます。たとえば、O(x**7) まで級数を表示したい場合は次のように記述します。

# リストA-3

# sinxをx=0の周りで級数展開
series(sin(x), x, n = 8)
\[x - \frac{x^{3}}{6} + \frac{x^{5}}{120} - \frac{x^{7}}{5040} + O\left(x^{8}\right)\]

 オプション引数 x0 で展開の中心点を指定できます (デフォルトは 0)。$\sin x$ を $x=1$ のまわりで展開してみます。

# リストA-4

# sinxをx=1の周りで級数展開
series(sin(x), x, x0 = 1)
\[\begin{align*}
\sin{\left (1 \right )} &+ \left(x - 1\right) \cos{\left (1 \right )} - \frac{\left(x - 1\right)^{2} \sin{\left (1 \right )}}{2} - \frac{\left(x - 1\right)^{3} \cos{\left (1 \right )}}{6}\\[6pt]
&+ \frac{\left(x - 1\right)^{4} \sin{\left (1 \right )}}{24} + \frac{\left(x - 1\right)^{5} \cos{\left (1 \right )}}{120} + O\left(\left(x - 1\right)^{6}; x\rightarrow 1\right)\end{align*}\]

 sympy.series() で得た式を sympy.poly() で多項式クラス (sympy.polys.polytools.Polyクラス) のオブジェクトに変換しておくと、色々なメソッドが使えるようになって便利です。$\cos x$ を級数展開して、coeffs()メソッドで $x^n$ の係数をすべて取り出してみます。

# リストA-5

# cosxをx=0の周りで級数展開
f = series(cos(x), x, n = 10)

# fを表示
display(f)

# fを多項式オブジェクトに変換
g = poly(f)

# 係数をすべて取得
c = g.coeffs()
c.reverse()

# 係数を表示
display(c)
\[\begin{align*}
&1 - \frac{x^{2}}{2} + \frac{x^{4}}{24} - \frac{x^{6}}{720} + \frac{x^{8}}{40320} + O\left(x^{10}\right)\\[6pt]
&\left [ 1, \quad 1, \quad - \frac{1}{2}, \quad \frac{1}{24}, \quad - \frac{1}{720}, \quad \frac{1}{40320}\right ]
\end{align*}\]

 

リーダブルコード ―より良いコードを書くためのシンプルで実践的なテクニック (Theory in practice)

新品価格
¥2,592から
(2019/8/6 17:55時点)

母関数の級数展開

 級数 $a_n$ が関数 $p(x)$ を級数展開したときの $x^n$ の係数となっているとき、$p(x)$ は級数 $a_n$ の母関数であるといいます(「級数を生み出す関数」という意味です)。たとえば、前の $2$ つの項を足して次の項をつくるフィボナッチ数列
 
\[1,\ 1,\ 2,\ 3,\ 5,\ 8,\ 13,\ 21,\ 34,\ \cdots\]
の母関数は
 
\[p(x)=\frac{1}{1-x-x^2}\]
で与えられることが知られています。$p(x)$ を sympy.series() で展開して確認してみましょう。

# リストA-6

# フィボナッチ数列の母関数
p = 1/(1 - x - x**2)

# フィボナッチ数列の母関数を級数展開
f = series(p, x, n=10)

display(f)
\[1 + x + 2 x^{2} + 3 x^{3} + 5 x^{4} + 8 x^{5} + 13 x^{6} + 21 x^{7} + 34 x^{8} + 55 x^{9} + O\left(x^{10}\right)\]

 係数は確かにフィボナッチ数列となっていることがわかります。

# リストA-7

# fを多項式オブジェクトに変換
g = poly(f)

# すべての係数を取得
c = g.coeffs()
c.reverse()

display(c)
\[\left [ 1, \quad 1, \quad 1, \quad 2, \quad 3, \quad 5, \quad 8, \quad 13, \quad 21, \quad 34, \quad 55\right ]\]