ウォリスの公式は円周率の近似などにあたって役立つ公式で、ウォリス積分(Wallis integral)の式から導出を行うことができます。当記事ではウォリス積分に基づくウォリスの公式の導出とPythonを用いた円周率の近似計算例について取り扱いました。
作成にあたっては「チャート式シリーズ 大学教養 微分積分」の第$4$章「積分($1$変数)」を主に参考にしました。
・数学まとめ
https://www.hello-statisticians.com/math_basic
ウォリスの公式の導出
ウォリス積分の公式
ウォリス積分の公式は下記のように表される。
$$
\large
\begin{align}
\int_{0}^{\frac{\pi}{2}} \cos^{2m}{x} dx &= \int_{0}^{\frac{\pi}{2}} \sin^{2m}{x} dx = \frac{(2m-1)!!}{2m!!} \cdot \frac{\pi}{2} \quad (1) \\
\int_{0}^{\frac{\pi}{2}} \cos^{2m+1}{x} dx &= \int_{0}^{\frac{\pi}{2}} \sin^{2m+1}{x} dx = \frac{2m!!}{(2m+1)!!} \quad (2) \\
m & \in \mathbb{N}
\end{align}
$$
ウォリスの公式の導出
$\displaystyle x \in \left[ 0, \frac{\pi}{2} \right]$の$x$について、$0 \leq \sin{x} \leq 1$であるので下記が成立する。
$$
\large
\begin{align}
\sin^{2n+2}{x} \leq \sin^{2n+1}{x} \leq \sin^{2n}{x}
\end{align}
$$
上記は$\displaystyle x \in \left[ 0, \frac{\pi}{2} \right]$の任意の$x$について成立するので、下記が成立する。
$$
\large
\begin{align}
\int_{0}^{\frac{\pi}{2}} \sin^{2n+2}{x} dx \leq \int_{0}^{\frac{\pi}{2}} \sin^{2n+1}{x} dx \leq \int_{0}^{\frac{\pi}{2}} \sin^{2n}{x} dx
\end{align}
$$
上記に対し、ウォリス積分をそれぞれ適用することで下記が得られる。
$$
\large
\begin{align}
\int_{0}^{\frac{\pi}{2}} \sin^{2n+2}{x} dx \leq & \int_{0}^{\frac{\pi}{2}} \sin^{2n+1}{x} dx \leq \int_{0}^{\frac{\pi}{2}} \sin^{2n}{x} dx \\
\frac{(2n+1)!!}{(2n+2)!!} \cdot \frac{\pi}{2} \leq & \frac{(2n)!!}{(2n+1)!!} \leq \frac{(2n-1)!!}{(2n)!!} \cdot \frac{\pi}{2} \\
\frac{2n+1}{2n+2} \cdot \frac{\pi}{2} \leq & \frac{(2n)!!}{(2n+1)!!} \cdot \frac{(2n)!!}{(2n-1)!!} \leq \frac{\pi}{2}
\end{align}
$$
ここで$\displaystyle \lim_{n \to \infty} \frac{2n+1}{2n+2} = 1$なので、下記が成立する。
$$
\large
\begin{align}
\lim_{n \to \infty} \frac{(2n)!!}{(2n+1)!!} \cdot \frac{(2n)!!}{(2n-1)!!} &= \frac{2^{2} \cdot 4^{2} \cdot 6^{2} \cdot \cdots}{(1 \cdot 3)(3 \cdot 5)(5 \cdot 7) \cdots} \\
&= \prod_{n=1}^{\infty} \frac{(2n)^{2}}{(2n-1)(2n+1)} \\
&= \frac{\pi}{2}
\end{align}
$$
円周率の近似
下記のような計算を行うことでウォリスの公式に基づいて円周率の近似を行うことができる。
import numpy as np
pi = 2.
for i in range(1000):
n = i+1
pi = pi * (2*n)**2 / ((2*n-1)*(2*n+1))
if n%200 == 0:
print("n: {}, estimated_pi: {:.3f}".format(n, pi)) ・実行結果
n: 200, estimated_pi: 3.138
n: 400, estimated_pi: 3.140
n: 600, estimated_pi: 3.140
n: 800, estimated_pi: 3.141
n: 1000, estimated_pi: 3.141