多角形の面積で円周率を求める
はじめに
半径 \(r\) の円の面積は \(\pi r^2\) ですから、半径1の円の面積は \(\pi\) です。
従って、半径1の円に内接する4角形, 8角形, 16角形 … の面積は円周率 \(\pi\) に近づいてゆくはずです。
図1のように、赤い4角形各辺の垂直二等分線が円に交わる点を新たな頂点として追加すれば青い8角形が得られ, 同様に各辺の垂直二等分線によって緑の16角形が得られ、以下同様に 32角形, 64角形 … を作図することができます。
4,8,16… 角形の面積の増分(図1の赤,青,緑…部分)を求め、それらを足し合わせてゆくことで \(\pi\) を計算します。
図1 : 4角形, 8角形, 16角形… の作図法と面積の増分(赤,青,緑…)
計算法
図2の濃い赤,青,緑の三角形に注目し、それぞれの面積を \(s_0, s_1, s_2\) としましょう。
図2 : 多角形を三角形に分解する
このように大きな三角形の斜辺と円の隙間に小さな三角形を 2 つずつ詰め込んでゆけば、円周率 \(\pi\) は次のように表せるはずです。 \begin{eqnarray} \pi &=& 4s_0+8s_1+16s_2+\cdots \end{eqnarray}
各部の長さを図3 のように定義します。 また、図には描かれていませんが \(y_0=0\) とします (長さ 0 なので描けない)。
図3 : 各部の長さの定義
(\(x_0\) は円の半径1に等しい)
\(x_1,y_1\) を斜辺とする直角三角形と \(x_2,y_2\) を斜辺とする直角三角形は、どちらも底辺の長さ \(x_0\) が円の半径=1 に等しいので、ピタゴラスの定理から次式が成り立ちます。 \begin{eqnarray} x_1^2+y_1^2 &=& 1 \\ \label{x12y121} x_2^2+y_2^2 &=& 1 \end{eqnarray}
図4 : \(x_1,y_1\) を斜辺とする直角三角形(左)と \(x_2,y_2\) を斜辺とする直角三角形(右)
小さな三角形を詰め込んでゆくやり方が同じなら \(n=2\) 以降でも同様に成り立ち、次式のように書けるでしょう。 \begin{eqnarray} x_n^2+y_n^2 &=& 1,\ n=0,1,2\cdots \label{x2y21} \end{eqnarray} 図5の黄色い三角形に注目すると、次式が成り立つことがわかります。 \begin{eqnarray} (2x_2)^2 &=& x_1^2+(1-y_1)^2 \\ &=& x_1^2+1-2y_1+y_1^2 \\ && 式(\ref{x12y121})より x_1^2+y_1^2=1 だから \nonumber\\ 4 x_2^2 &=& 2-2y_1 \\ && 両辺を 4 で割って \nonumber\\ x_2^2 &=& \frac{1-y_1}{2} \end{eqnarray}
図5 : \(x_1,y_1\) と \(x_2\) の関係
小さな三角形を詰め込んでゆくやり方が同じなら \(n=2\) 以降でも同様に成り立ち、次式のように書けるでしょう。 \begin{eqnarray} x_{n+1}^2 &=& \frac{1-y_n}{2},\ n=0,1,2\cdots \label{x2} \end{eqnarray} \(x_{n+1}^2\) が求まれば、\(y_{n+1}^2\) は式(\ref{x2y21})と式(\ref{x2})より、次のように計算できます。 \begin{eqnarray} y_{n+1}^2 &=& 1-x_{n+1}^2 \\ &=& 1-\frac{1-y_n}{2} \\ &=& \frac{1+y_n}{2} \label{y2} \end{eqnarray} 図6の色付き三角形の面積は次式で表せます。 \begin{eqnarray} s_n &=& \frac{1}{2} x_n (1-y_n),\ n=0,1,2\cdots \label{sn} \end{eqnarray}
図6 : 三角形の面積 \(s_0,s_1,s_2,\cdots\)
\(x_0=1,\ y_0=0\) から出発して、式(\ref{x2}),式(\ref{y2})を使って \(x_0,\ y_0\to x_1,\ y_1\to x_2,\ y_2\to \cdots\) と計算してゆけば、円周率 \(\pi\) は次のように計算できます。 \begin{eqnarray} \pi &=& 4s_0+8s_1+16s_2+\cdots \\ && 式(\ref{sn})より s_n=\frac{1}{2} x_n (1-y_n) だから \nonumber\\ &=& 4\cdot\frac{1}{2} x_0 (1-y_0) + 8\cdot\frac{1}{2} x_1 (1-y_1) + 16\cdot\frac{1}{2} x_2 (1-y_2) + \cdots \\ &=& 2 x_0 (1-y_0) + 4 x_1 (1-y_1) + 8 x_2 (1-y_2) + \cdots \\ &=& \sum_{n=0}^\infty 2^{n+1} x_n (1-y_n) \end{eqnarray}
コード例
C++ で書くと次のようになります。
cout << fixed << setprecision(15);
double x=1, y=0, t=2, pi=0;
for (int n=0; n<26; n++)
{
pi+= t*x*(1-y); t*=2;
x = sqrt((1-y)/2);
y = sqrt((1+y)/2);
}
cout << "pi=" << pi << endl;
変数 t は for ループ内で \(t=2^{n+1}\) に等価です。
計算結果
for ループ内の pi+=... の下で表示した x, y, pi の各値を以下に示します (pi の赤字部分は正しくありません)。
n | x | y | pi |
0 | 1.000000000000000 | 0.000000000000000 | 2.000000000000000 |
1 | 0.707106781186548 | 0.707106781186548 | 2.828427124746190 |
2 | 0.382683432365090 | 0.923879532511287 | 3.061467458920718 |
3 | 0.195090322016128 | 0.980785280403230 | 3.121445152258052 |
4 | 0.098017140329561 | 0.995184726672197 | 3.136548490545939 |
5 | 0.049067674327418 | 0.998795456205172 | 3.140331156954753 |
6 | 0.024541228522912 | 0.999698818696204 | 3.141277250932772 |
7 | 0.012271538285719 | 0.999924701839145 | 3.141513801144301 |
8 | 0.006135884649156 | 0.999981175282601 | 3.141572940367091 |
9 | 0.003067956762969 | 0.999995293809576 | 3.141587725277160 |
10 | 0.001533980186282 | 0.999998823451702 | 3.141591421511200 |
11 | 0.000766990318753 | 0.999999705862882 | 3.141592345570118 |
12 | 0.000383495187567 | 0.999999926465718 | 3.141592576584872 |
13 | 0.000191747597257 | 0.999999981616429 | 3.141592634338563 |
14 | 0.000095873799280 | 0.999999995404107 | 3.141592648776986 |
15 | 0.000047936899495 | 0.999999998851027 | 3.141592652386591 |
16 | 0.000023968449458 | 0.999999999712757 | 3.141592653288993 |
17 | 0.000011984225887 | 0.999999999928189 | 3.141592653514593 |
18 | 0.000005992115260 | 0.999999999982047 | 3.141592653570993 |
19 | 0.000002996059946 | 0.999999999995512 | 3.141592653585094 |
20 | 0.000001498029973 | 0.999999999998878 | 3.141592653588619 |
21 | 0.000000749033514 | 0.999999999999719 | 3.141592653589500 |
22 | 0.000000374535284 | 0.999999999999930 | 3.141592653589721 |
23 | 0.000000187304692 | 0.999999999999982 | 3.141592653589776 |
24 | 0.000000093652346 | 0.999999999999996 | 3.141592653589790 |
25 | 0.000000047121609 | 0.999999999999999 | 3.141592653589793 |
26 回反復して得た \(2^{27}\)=1億3421万7728角形の面積 3.141592653589793 は、円周率 \(\pi\) に小数点以下 15 桁まで一致しています。
関連項目
付記
本方式と等価な結果を 1995 年に Kirby Urner さんという方が sci.math に公表されていたらしいのですが、投稿が見当たらず導出方法を確認できませんでした。
【情報元】http://www.pluto.ai.kyutech.ac.jp/~matumoto/dvi/pi.pdf の p14