円周率を求めよう!
- (by K, 2004.01.30)
- 中学生にも理解できる範囲で、円周率の計算方法を説明できないかと言われて、考えた話。
- 読み物としてちょっと面白いだけで、学術的価値はありません。
円周率 パイ
数学的な計算方法
- まず1辺が1メートルの正方形を考えます。
- この左下の頂点を原点として、下の辺をx軸に見立てて0.1ずつ目盛を付けるとします。同じく左の辺をy軸に見立てて0.1ずつ目盛を付けることにします。さて、こうすると、0.1*0.1な正方形が100個できたことになります。
- 次に、原点を中心として、(1.0, 0.0)から(0.0, 1.0)を円弧にするような円盤を考えます。つまり半径1の円を1/4だけ切り取ったような形です。・・・さてこの円盤の中は、一体いくつの0.1*0.1のミニ正方形が含まれているでしょうか?
- 「仮に」78個だとしたら、この1/4円の面積は0.78平方メートルであると近似できることになります。
- そして半径1の円の面積はパイであって、この円盤はその1/4なので、
- と推定できることになります。これはなかなか良い感じです。
- 次に、ミニ正方形がいくつ含まれているのかを判定する方法を考えます。
- たとえば4つの頂点が
A(0.4, 0.3) B(0.5, 0.3)
D(0.4, 0.2) C(0.5, 0.2)
- のミニ正方形ABCDは円盤の中に含まれているでしょうかいないでしょうか?
- これは三平方の定理を使うと判定できるのです。 0.4*0.4+0.3*0.3=0.25 なので、A点と原点の距離は0.5ということになります。これは円の半径よりも小さいので、A点は円盤の中です。
- 一般に、点(x, y)が円盤の中かどうかは、x*x+y*yを計算して、それが1以下かどうかを見ればいいのです。
- これで4点とも円盤の中であれば、そのミニ正方形は間違いなく円盤の中にあります。また4点とも円盤の外であれば、そのミニ正方形は間違いなく円盤の外にあります。
- 4点のうち、一部だけが円盤の中に入っている場合は、そのミニ正方形が円弧によって分割されていることを意味します。・・・この分割された面積を0%とみなすか100%とみなすかで、円周率の推定の下限と上限を求めることができます。
- この計算方法のポイントは、円周率の推定にあたって、定規などで「測る」という作業をなくしてしまったことです。
- だから0.1ずつの目盛を0.000000001ずつの目盛とかにすれば、相応に高い精度で円周率を求められます。
- それどころか、目盛を十分に細かくすれば、どこまででも好きな精度で円周率を求められます。
プログラムで実験
- では、実際に計算したらどうなるかやってみましょう。
#include <stdio.h>
#define N 10
int main()
{
int i, j, k, c[5];
double x0, y0, x1, y1;
for (k = 0; k < 5; k++)
c[k] = 0;
for (j = 0; j < N; j++) {
y0 = (double) j / N;
y1 = (double) (j + 1) / N;
for (i = 0; i < N; i++) {
x0 = (double) i / N;
x1 = (double) (i + 1) / N;
k = 0;
if (x0 * x0 + y0 * y0 <= 1.0) k++;
if (x1 * x1 + y0 * y0 <= 1.0) k++;
if (x0 * x0 + y1 * y1 <= 1.0) k++;
if (x1 * x1 + y1 * y1 <= 1.0) k++;
c[k]++;
}
}
for (k = 0; k < 5; k++)
printf("c[%d]=%d\n", k, c[k]);
return 0;
}
- とても単純なプログラムです。<math.h>すら使っていません。ミニ正方形の4つの頂点が円盤の中にあるかどうかを判定して、分類して集計しているだけです。
- 実行するとこんな結果になります。
c[0]=12, c[1]=5, c[2]=8, c[3]=6, c[4]=69
- これはつまり、完全に円盤の中に入っていたものが69個で、頂点が3つだけ円盤の中に入ってたのが6個で、2つ入っていたのが8個で、1つ入っていたのが5つで、完全に円盤の外にあったものが12個あった、ということです。
- c[1〜3]の面積を0とみなしてパイを推定すると、2.76にしかなりません。これはかなり3.14から遠いです。がっかりですね。・・・逆にc[1〜3]を100%のミニ正方形だと仮定してパイを推定すると、3.52となり、これも3.14から遠いです。
- でもとりあえず上限と下限は分かりました。Nを大きくすれば、より誤差を狭めることが出来そうです。
- ここでN=100へ進んでもいいのですが、c[1]を25%、c[2]を50%、c[3]を75%として計算すると、
69+5*0.25+8*0.50+6*0.75=78.75
- となって、円盤の中にはミニ正方形が78.75個相当あると推定することもできます。この推定値からパイを計算すればなんと3.15となり、N=10の計算だけでも3.14にかなり近づくことができます!
- ではN=100で計算してみます。
c[0]=2047, c[1]=58, c[2]=82, c[3]=59, c[4]=7754
- これを使うとパイの値は、3.1016〜3.1812だということが分かります。
- さらに先の25%〜75%の方法を使えば、パイの値の近似値は3.1415になります。
- もっと先に進んでN=1000にしてみます。
c[0]=213614, c[1]=585, c[2]=828, c[3]=586, c[4]=784387
上限下限:3.137548〜3.145544, (25%〜75%で近似)=3.141547
- さらにN=1万をやってみます。
c[0]=21450238, c[1]=5857, c[2]=8284, c[3]=5858, c[4]=78529763
上限下限:3.14119052〜3.14199048, (25%〜75%で近似)=3.14159051
考察
- 上限下限方式は数学的な意味が一番高いと思います。いくつより大きくていくつより小さいかはっきりするからです。たとえばN=1万の結果を見れば 3.14119052〜3.14199048 なので、3.141までは確定させることができます。しかしミニ正方形を1億個も作って比較して分類したのに、たったの4桁しか分からないというのは、なかなかさみしい気もします。
- これを改善する方法としては、まず上限値と下限値を足して2で割ってみる、という単純な方法があります。もう一つは、上でやったように、25%〜75%で近似するという方法もあります。僕の直感では、単純に足して2で割るよりは、25%〜75%で近似するほうが正解に近いかなと思っているのですが、さてどうでしょうか?
- 推定したパイが3.14159265からどのくらいずれているのか?
| N=10 | N=100 | N=1000 | N=1万 | 単純に足して2で割る | -0.00159265 | -0.00019265 | -0.00004665 | -0.00000215 | 25%〜75%で近似 | +0.00840735 | -0.00009265 | -0.00004565 | -0.00000214 |
- こうして見ると、25%〜75%の近似が特に優れているとかそういうことはなくて、適当にやった上限と下限の単純平均でもそんなに悪くなさそうです。
リンク
雑談
- 小学校のゆとり教育が始まったときに、円周率を3に近似するのは円を円に内接する正六角形に近似するのに相当するという批判があった。実は愚かにも僕もそれをすっかり信じた。しかしこれは少々おかしい気がする。正六角形の面積は、正三角形x6に等しいが、斜辺が1の正三角形の面積はsqrt(3)/4であって、これを6倍しても2.6にすら届かない(半径1の円の面積は円周率に等しい)。
- いやおかしくなかった。円周で考えると、正六角形の周の長さは6だから、直径の2で割れば3なのか。なるほど。
- となると日能研で紹介されていた東京大学の入試問題の、円周率が3.05以上であることを証明せよというのは、正12角形の周の長さを出せばいいのかな。それでも足りないのかな。足りなければ24角形だなあ。・・・いつか暇になったら計算しよう。
こめんと欄
- ルール違反でページを編集する荒らしが何度もきたので、当分の間このページは凍結します。コメントを書きたい場合はimpressionsへ。 -- K 2008-04-18 (金) 17:39:07
- ちなみに凍結中でもこのようにコメントは書けますが、書き損じた場合に修正ができません。書き損じない自信があればここに書いてもいいです。 -- K 2008-04-18 (金) 17:42:59
- 書き始めから数えてもう10年以上になるのですが、ずっと書き途中になっていました。さすがに恥ずかしくなったので、予定していた内容を最後まで書きました。 -- K 2015-11-10 (火) 18:41:34
|