スポンサーリンク
はじめに
円周率といえば、理系の人間ならだれでも(?)憧れを抱く無理数ですね。今回はそんな円周率をプログラムによる計算で求める方法を実践していきたいと思います。
環境
OS:Windows 7
IDE:Visual Studio 2015 C++
CPU:Intel Core i5-2520M 2.50GHz
RAM:4.00GB
ライプニッツの公式について
ライプニッツの公式は、単純な分数の和と差だけで円周率を表現してしまう公式です。
$$π/4 = 1-1/3+1/5-1/7+\cdots$$
これをプログラムに起こして、実際に計算させてみます。
この公式の証明を知りたい方は、Wikiへどうぞ。。。
プログラムにする
とりあえず、long double型でやってみることにします。
#include <iostream>
#include <math.h>
#include <stdio.h>
using namespace std;
int main()
{
int i = 1;
int n = 0;
long double temp = 0.0;
cin >> n;
for (i = 1; i < n; i++) {
if ((i % 2) == 1) {
temp += (1.0 / (2 * i - 1));
}
else {
temp -= (1.0 / (2 * i - 1));
}
}
long double ans = temp * 4.0;
printf("%15.14lf", ans);
return 0;
}
足したり引いたりする分数の数(精度)を数字で入力すると、計算が始まります。並列化の手法を用いれば、さらに高速化できると思いますが、今回は単純な計算手順の説明にとどめさせていただきます。
公式をそのままプログラム化してやるだけでコンピュータが勝手に計算していってくれるというのは、初めてのプログラミングではとても感動したものです。
スポンサーリンク収束
この手法ではlong double型レベルの精度しか出ない上、お遊び程度の収束速度ですが、まあ、いいでしょう。
分数の個数 | 10 | 100 | 1,000 | 10,000 | 100,000 | 1,000,000 | 10,000,000 | ∞理論値 |
---|---|---|---|---|---|---|---|---|
値 | 3.25236593471888 | 3.15169340607112 | 3.14259365434004 | 3.14169266359053 | 3.14160265368972 | 3.14159365359077 | 3.14159275358980 | 3.14159265358979 |
以上のような結果となりました。10,000,000コの分数を使って計算すると、3.141592までが合致していて、そのあとも5ケタ合致していることがわかります。ただ、これ以上計算量を増やすと、時間がかかるようになってしまい、実用的ではありません。先ほども述べましたが、「お遊び」程度に考えて、やってみてください。