マチンの公式を使った円周率計算

これの続き。マチンの公式を用いると円周率を計算することができる。マチンの公式は

\begin{align}
\frac{\pi}{4}=4 \tan^{-1} \frac{1}{5} – \tan^{-1} \frac{1}{239}
\end{align}

で表される。これをグレゴリー・ライプニッツ級数で使った結果を使えば

\begin{align}
\pi= \sum_{n=0}^{\infty} \left \{ \frac{ 16 \cdot (-1)^{n}}{2n+1} \left( \frac{1}{5} \right)^{2n+1}-\frac{4 \cdot (-1)^{n}}{2n+1} \left( \frac{1}{239} \right)^{2n+1} \right \}
\end{align}

となるので円周率と等しくなる。収束は早い。

import numpy as np
from matplotlib import pyplot as plt

N = 30
x = np.zeros(N)
x1 = np.zeros(N)
x2 = np.zeros(N)

for n in range(N):
    alpha = 2 * n + 1
    x1[n] = x1[n - 1] + (16 * (1 / 5) ** alpha * (-1) ** n) / alpha
    x2[n] = x2[n - 1] + (4 * (1 / 239) ** alpha * (-1) ** n) / alpha
    x[n] = x1[n] - x2[n]
print(x)