理系の憧れ！円周率の求め方　～ランダム編～

Tweet

はじめに

円周率は理系が崇拝する数値で、さまざまな求め方が数学者たちによって研究されてきました。今回は、その中でも、最も単純な求め方を実践していきたいと思います。

環境

OS:Windows 7

IDE:Visual Studio 2015 C++

CPU:Intel Core i5-2520M 2.50GHz

RAM:4.00GB

モンテカルロ法による円周率計算

モンテカルロ法とは、乱数によって数値を計算する手法です。

ここに紙が１枚、鉛筆が１本あるとします。最初に、紙に1cm×1cmの正方形を描きます。その中に半径1cmの四分円を描きます。そして、必ず正方形の中に鉛筆が落ちるように、ランダムに鉛筆を落としていきます。

これを何回、何万回と繰り返します。
正方形の面積は$$1cm^2$$四分円の面積は$$(π/4) cm^2$$つまり、正方形に入っている四分円の中に鉛筆が落ちる確率は、$$(四分円の面積)/(正方形の面積) = π/4$$となります。この性質を利用して、計算していきます。

プログラムにする

最初に打つ点の数を入力します。四分円の中に入っているかどうかは、ランダムな点と四分円の中心との距離が１以上、１未満で判定することとしました。double精度でやることにします。

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <iostream>
#include <math.h>
using namespace std;
double frand() {
	double num;
	num = (double)rand() / 32767.0;
	return num;
}
int main()
{
	srand(time(NULL));
	double x, y;
	int countall, countin;
	countall = 0;
	countin = 0;
	x = 0;
	y = 0;
	int n=0;
	cin >> n;
	for (int i = 0; i < n; i++) {
		x = frand();
		y = frand();
		if (pow(x,2.0) + pow(y,2.0) <= 1) {
			countin++;
		}
		countall++;
	}
	double pi = ((double)countin / countall)*4;
	cout << pi;
    return 0;
}

とても単純ですね。

収束

この手法はPCのリソースをかなり使うにもかかわらず、収束が大変遅いです。原理が単純だということだけが取り柄です。

以上のような結果となりました。乱数を使用しているため、実行するたびに値が変わります。ただし、回数を増やせば増やすほど、正確な値に近づいていきます。とはいえ、時間的にせいぜい１０億くらいが限界ではないでしょうか…これだけの回数を重ねても、３ケタまでしか合っていません。やはり、他のやり方で、再現性のあるプログラムを作るべきだと思います。

というわけで、次回はライプニッツの公式編を書きたいと思います。お楽しみに。