MaximaでFourier級数展開する

数学の具体的な計算にMaximaを使って、数学もMaximaも同時に学んでしまいましょう。

今回はMaximaを使って関数のFourier（フーリエ）級数展開を見てみたいと思います。

いくつかの関数について、そのFourier級数展開から円周率の級数表式が得られますが、

おまけとしてこれをMaximaで数値的に確かめてみます。

また、正弦関数の無限乗積による恒等式についてもコメントしておきます。

三角関数の直交性
Fourier級数展開
Fourier級数展開の例

三角関数の直交性

関数

$\displaystyle{\frac{1}{\sqrt{2\pi}},~\frac{\cos x}{\sqrt{\pi}},~\frac{\cos 2x}{\sqrt{\pi}},~\frac{\cos 3x}{\sqrt{\pi}},\, \dots,~\frac{\sin x}{\sqrt{\pi}},~\frac{\sin 2x}{\sqrt{\pi}},~\frac{\sin 3x}{\sqrt{\pi}},\,\dots}$

は $-\pi\leq x\leq \pi$ の変域で正規化された直交関数系をなします。実際に計算してみましょう：

$\displaystyle{\int^\pi_{-\pi} dx \frac{\cos mx}{\sqrt{\pi}} \frac{\cos nx}{\sqrt{\pi}}=\delta_{m,n}}$

declare([m,n], integer)$
integrate(cos(m*x)/sqrt(%pi)*cos(n*x)/sqrt(%pi), x, -%pi, %pi);
integrate(cos(m*x)/sqrt(%pi)*cos(m*x)/sqrt(%pi), x, -%pi, %pi);

$\displaystyle{\int^\pi_{-\pi} dx \frac{\sin mx}{\sqrt{\pi}} \frac{\sin nx}{\sqrt{\pi}}=\delta_{m,n}}$

integrate(sin(m*x)/sqrt(%pi)*sin(n*x)/sqrt(%pi), x, -%pi, %pi);
integrate(sin(m*x)/sqrt(%pi)*sin(m*x)/sqrt(%pi), x, -%pi, %pi);

$\displaystyle{\int^\pi_{-\pi} dx \frac{\sin mx}{\sqrt{\pi}} \frac{\cos nx}{\sqrt{\pi}}=0}$

integrate(sin(m*x)/sqrt(%pi)*cos(n*x)/sqrt(%pi), x, -%pi, %pi);
integrate(sin(m*x)/sqrt(%pi)*cos(m*x)/sqrt(%pi), x, -%pi, %pi);

$\displaystyle{\int^\pi_{-\pi} dx \frac{1}{\sqrt{2\pi}} \frac{\cos mx}{\sqrt{\pi}}=0}$

integrate(1/sqrt(2*%pi)*cos(m*x)/sqrt(%pi), x, -%pi, %pi);

$\displaystyle{\int^\pi_{-\pi} dx \frac{1}{\sqrt{2\pi}} \frac{\sin mx}{\sqrt{\pi}}=0}$

integrate(1/sqrt(2*%pi)*sin(m*x)/sqrt(%pi), x, -%pi, %pi);

$\displaystyle{\int^\pi_{-\pi} dx \frac{1}{\sqrt{2\pi}} \frac{1}{\sqrt{2\pi}}=1}$

integrate(1/sqrt(2*%pi)*1/sqrt(2*%pi), x, -%pi, %pi);

以上より正規直交性がわかります。

Fourier級数展開

変域 $-\pi\leq x\leq \pi$ で定義された関数 $f(x)$ の正規化された三角関数による展開

$\displaystyle{f(x)=A_0\frac{1}{\sqrt{2\pi}}+A_1\frac{\cos x}{\sqrt{\pi}}+B_1\frac{\sin x}{\sqrt{\pi}}+A_2\frac{\cos 2x}{\sqrt{\pi}}+B_2\frac{\sin 2x}{\sqrt{\pi}}+\dots\\ \, \ \quad=A_0\frac{1}{\sqrt{2\pi}}+\sum_{n=1}^\infty \left(A_n\frac{\cos nx}{\sqrt{\pi}}+B_n\frac{\sin nx}{\sqrt{\pi}}\right)}$

をFourier級数展開といいます。

関数の正規直交性から、展開係数は

$\displaystyle{A_0=\int^\pi_{-\pi} dx \, \frac{1}{\sqrt{2\pi}} \, f(x),\quad A_n=\int^\pi_{-\pi} dx \, \frac{\cos nx}{\sqrt{\pi}} \,f(x), \quad B_n=\int^\pi_{-\pi} dx \, \frac{\sin nx}{\sqrt{\pi}}\, f(x)}$

によって求めることができます。

Fourier級数展開の例

ここでFourier級数の具体例を見てみましょう。

例は寺澤寛一『自然科学者のための数学概論』（岩波書店）を参照しました。

例１：f(x)=x

$\displaystyle{f(x)=x \quad (-\pi\leq x\leq \pi)}$

f(x):=x$
A0:integrate(1/sqrt(2*%pi)*f(x), x, -%pi, %pi);
A(n):=integrate(cos(n*x)/sqrt(%pi)*f(x), x, -%pi, %pi)$
B(n):=integrate(sin(n*x)/sqrt(%pi)*f(x), x, -%pi, %pi)$
makelist(A(n), n, 1, 10);
makelist(B(n), n, 1, 10);

より

$\displaystyle{A_0=0, \quad A_n=0, \quad B_n=(-1)^{n+1}\sqrt{\pi}\,\frac{2}{n} }$

そのFourier級数展開は

$\displaystyle{x=2\left\{\frac{\sin x}{1}-\frac{\sin 2x}{2}+\frac{\sin 3x}{3}-\frac{\sin 4x}{4}+\dots\right\}}$

となります。

また、Fourier級数の展開係数を求めるコマンドはMaximaに実装されています：

load(fourie)$ /* fourier ではなく fourie であることに注意 */
totalfourier(x, x, %pi);

ただし、このコマンドで得られる展開係数は

$\displaystyle{a_0=\frac{A_0}{\sqrt{2\pi}},\quad a_n=\frac{A_n}{\sqrt{\pi}},\quad b_n=\frac{B_n}{\sqrt{\pi}}}$

として再定義されたものです。よってこのコマンドで得られるFourier級数展開は

$\displaystyle{f(x)=a_0+\sum_{n=1}^\infty \left(a_n\cos nx+b_n\sin nx\right)}$

という形になります。

円周率（Leibniz／Eulerの級数）

さて、上で得られたFourier級数展開で $x=\pi/2$ とおくことで等式

$\displaystyle{\frac{\pi}{4}=1-\frac{1}{3}+\frac{1}{5}-\frac{1}{7}+\dots}$

がわかります。この級数はLeibnizの級数、またはEulerの級数と呼ばれます。

この級数を数値的に検証してみましょう：

P[n]:=4+4*sum((-1)^k/(2*k+1), k, 1, n)$

とすることで、級数を $n$ 項まで足したときの値が数値的に求められます。

例えば１万項を足したときには

bfloat(P[10000]);

より、級数の値は

$3.1416926\dots$

となります。円周率の値は

bfloat(%pi);

より

$\pi=3.14159265358979\dots$

ですから、１万＝ $10^4$ 項の和で４桁程度までは数値的に再現できていることがわかります。

10万項の和だと

bfloat(P[100000]);

より

$3.14160265\dots$

となって精度が１桁ほど上がります。

しかし計算量が多いようで、出力には 20.2 秒ほどかかりました。

ちなみに、１万項の場合で 0.3 秒、20万項の場合では 80.5 秒かかりました。

項の数を $n$ 倍すると時間はだいたい $n^2$ 倍になるようなので、

欲張って100万項の場合を計算しようとすると、2000秒＝30分くらいかかるはずです。

（もちろん、使っているコンピュータの性能によりますが）

例２：f(x)=|x|

$\displaystyle{f(x)=|x| \quad (-\pi\leq x\leq \pi)}$

kill(all)$ load(fourie)$
totalfourier(abs(x), x, %pi);

より展開係数は

$\displaystyle{a_0=\frac{\pi}{2}, \quad a_n=\frac{2((-1)^n-1)}{\pi n^2}, \quad b_n=0 }$

なので、Fourier級数展開

$\displaystyle{|x|=\frac{\pi}{2}-\frac{4}{\pi}\left(\frac{\cos x}{1^2}+\frac{\cos 3x}{3^2}+\frac{\cos 5x}{5^2}+\dots\right)}$

を得ます。とくに $x=0$ とおくことで等式

$\displaystyle{\frac{\pi^2}{8}=\frac{1}{1^2}+\frac{1}{3^2}+\frac{1}{5^2}+\frac{1}{7^2}+\dots}$

がわかります。

級数を数値的に見てみましょう：

kill(all)$
P[n]:=8*sum(1/(2*k+1)^2, k, 0, n)$

から、１万項まで足したものから円周率を求めると、数値的には

bfloat(P[10000])$ sqrt(%);

より

$3.1415608\dots$

です。

例３：f(x)=x^2

$\displaystyle{f(x)=x^2 \quad (-\pi\leq x\leq \pi)}$

kill(all)$ load(fourie)$
totalfourier(x^2, x, %pi);

より展開係数は

$\displaystyle{a_0=\frac{\pi^2}{3}, \quad a_n=\frac{4(-1)^n}{n^2},\quad b_n=0 }$

なので、Fourier級数展開

$\displaystyle{x^2=\frac{\pi^2}{3}-4\left(\frac{\cos x}{1^2}-\frac{\cos 2x}{2^2}+\frac{\cos 3x}{3^2}-\frac{\cos 4x}{4^2}+\dots\right)}$

を得ます。とくに $x=0$ とおくことで

$\displaystyle{\frac{\pi^2}{12}=\frac{1}{1^2}-\frac{1}{2^2}+\frac{1}{3^2}-\frac{1}{4^2}+\frac{1}{5^2}+\dots}$

がわかります。

数値的に見てみましょう：円周率を求めてみると

kill(all)$
P[n]:=12*sum((-1)^k/(k+1)^2, k, 0, n)$
bfloat(P[10000])$ sqrt(%);

から

$3.14159266\dots$

です。

例４：f(x)=-1,1

$\displaystyle{f(x)=\begin{cases} 1 \quad \ \ \ (0< x\leq \pi) \\ -1 \quad (-\pi\leq x\leq 0) \end{cases} }$

ここで、奇関数

$\displaystyle{f(-x)=-f(x)}$

に対して、その性質を使うとFourier級数の展開係数が

$\displaystyle{A_0=\int^\pi_{-\pi} dx \, \frac{1}{\sqrt{2\pi}} \, f(x) \\ \quad =\int^0_{-\pi} dx \, \frac{1}{\sqrt{2\pi}} \, f(x) +\int^\pi_0 dx \, \frac{1}{\sqrt{2\pi}} \, f(x) \\ \quad = \int^0_{\pi} d(-x) \, \frac{1}{\sqrt{2\pi}} \, f(-x) +\int^\pi_0 dx \, \frac{1}{\sqrt{2\pi}} \, f(x) \\ \quad = \int_0^{\pi} dx \, \frac{1}{\sqrt{2\pi}} \, f(-x) +\int^\pi_0 dx \, \frac{1}{\sqrt{2\pi}} \, f(x)=0 }$

また、 $\cos$ と $\sin$ がそれぞれ偶関数、奇関数であることから、同様にして

$\displaystyle{A_n=0 \\ B_n=2\int^\pi_0 dx \, \frac{\sin nx}{\sqrt{\pi}} \, f(x)}$

となります。

（ $f(x)$ が偶関数である場合には $B_n=0$ となり、

$A_n$ は積分範囲を $\displaystyle{0\leq x \leq \pi}$ にすることができます）

これを使うと上記関数のFourier級数展開は

f(x):=1$
B(n):=2*integrate(sin(n*x)/sqrt(%pi)*f(x), x, 0, %pi)$
makelist(B(n), n, 1, 10);

より

$\displaystyle{B_n=\frac{2(1+(-1)^{n+1})}{n\sqrt{\pi}}}$

よって

$\displaystyle{f(x)=\frac{4}{\pi}\left(\frac{\sin x}{1}+\frac{\sin 3x}{3}+\frac{\sin 5x}{5}+\frac{\sin 7x}{7}+\dots\right)}$

がわかります。

例５：f(x)=sin(ax)

$\displaystyle{f(x)=\sin ax }$

ただし $a$ は整数でない実数とします。

kill(all)$ load(fourie)$
totalfourier(sin(a*x), x, %pi);

より、奇関数なので $a_n=0$ そして

$\displaystyle{b_n=\frac{2}{\pi}\left(\frac{\sin(n\pi-a\pi)}{2(n-a)}-\frac{\sin(n\pi+a\pi)}{2(n+a)}\right)\\ \quad =\frac{2}{\pi}\left(\frac{\sin(n\pi)\cos(a\pi)-\cos(n\pi)\sin(a\pi)}{2(n-a)}-\frac{\sin(n\pi)\cos(a\pi)+\cos(n\pi)\sin(a\pi)}{2(n+a)}\right) \\ \quad =\frac{2(-1)^{n+1}\sin(a\pi)}{\pi}\left(\frac{1}{2(n-a)}+\frac{1}{2(n+a)}\right)\\ \quad =\frac{2(-1)^{n}\sin(a\pi)}{\pi}\frac{n}{a^2-n^2} }$

ゆえに

$\displaystyle{\sin ax=-\frac{2\sin(a\pi)}{\pi}\left( \frac{\sin x}{a^2-1^2} -\frac{2 \sin 2x}{a^2-2^2} +\frac{3 \sin 3x}{a^2-3^2} -\frac{4 \sin 4x}{a^2-4^2}+\dots\right) }$

例６：f(x)=cos(ax)

$\displaystyle{f(x)=\cos ax }$

ただし $a$ は整数でない実数とします。

kill(all)$ load(fourie)$
totalfourier(cos(a*x), x, %pi);

より、偶関数なので $b_n=0$ そして

$\displaystyle{a_n=\frac{2}{\pi}\left(\frac{\sin(n\pi-a\pi)}{2(n-a)}+\frac{\sin(n\pi+a\pi)}{2(n+a)}\right)\\ \quad =\frac{2}{\pi}\left(\frac{\sin(n\pi)\cos(a\pi)-\cos(n\pi)\sin(a\pi)}{2(n-a)}+\frac{\sin(n\pi)\cos(a\pi)+\cos(n\pi)\sin(a\pi)}{2(n+a)}\right) \\ \quad =\frac{2(-1)^{n+1}\sin(a\pi)}{\pi}\left(\frac{1}{2(n-a)}-\frac{1}{2(n+a)}\right)\\ \quad =\frac{2(-1)^{n}\sin(a\pi)}{\pi}\frac{a}{a^2-n^2}\\ a_0 =\frac{\sin(a \pi )}{a \pi }}$

ゆえに

$\displaystyle{\cos ax=\frac{2a\sin a\pi}{\pi}\left(\frac{1}{2a^2} -\frac{\cos x}{a^2-1^2} +\frac{\cos 2x}{a^2-2^2} -\frac{\cos 3x}{a^2-3^2} +\dots\right) }$

を得ます。

正弦関数の無限乗積

ここで $x=\pi$ とおいて両辺を $\sin a\pi$ で割ると

$\displaystyle{\frac{\cos a\pi}{\sin a\pi}=\frac{2a}{\pi}\left(\frac{1}{2a^2} +\frac{1}{a^2-1^2} +\frac{1}{a^2-2^2}+\frac{1}{a^2-3^2} +\dots\right) \\ \qquad \ \ =\frac{1}{\pi a} -\frac{2a}{\pi}\left(\frac{1}{1^2-a^2} +\frac{1}{2^2-a^2}+\frac{1}{3^2-a^2} +\dots\right) }$

となります。両辺をさらに $a$ について $0$ から $z$ まで積分すると

$\displaystyle{\frac{1}{\pi}\log (\sin \pi z)-\frac{1}{\pi}\log(\pi z)=-\int^z_0 da\frac{2a}{\pi}\left(\frac{1}{1^2-a^2} +\frac{1}{2^2-a^2}+\frac{1}{3^2-a^2} +\dots\right) }$

ここで

$\displaystyle{\int^z_0 da\frac{2a}{\pi}\,\frac{1}{n^2-a^2}=-\frac{1}{\pi}\left[ \log(n^2-a^2)\right]_{a=0}^z=-\frac{1}{\pi} \log\left(1-\frac{z^2}{n^2}\right)}$

より

$\displaystyle{\log \left(\frac{\sin \pi z}{\pi z}\right)=\sum_{n=1}^\infty \log\left(1-\frac{z^2}{n^2}\right)= \log \prod_{n=1}^\infty\left(1-\frac{z^2}{n^2}\right) }$