MENU

モンテカルロ法

モンテカルロ法を用いてπを求める

ja.wikipedia.org

f:id:Manao55:20200816235645j:plain

ある問題を確率（乱数）を用いて解くことをモンテカルロ法といいます。円周率πをこの方法で求めてみます。

0 ~ 1の一葉実数乱数を2つ発生させ、それらを x、y とします。こうした乱数の組をいくつか発生させると、1×1 の正方形の中に、(x,y) で示される点は均一にばらまかれると考えられます。

したがって、正方形の面積と 1/4 円の面積の比は、ばらまかれた乱数の数に比例するはずです。

1/4 円の中にばらまかれた乱数の数を a、円外にばらまかれた乱数の数を b とすると、次の関係が成立します。

{ \displaystyle \frac{\pi}{4}:1 = a:a+b }
{ \displaystyle \pi=\frac{4a}{a+b}=\frac{4a}{n} }

プログラム

#include <iostream>
#include <iomanip>
#include <math.h>
#include <stdlib.h>
using namespace std;

#define NUM 330000

double rnd(void);

int main(){
    double x,y,pai;
    int i,in=0;

    for (i=0;i<NUM;i++){
        x=rnd();
        y=rnd();
        if (x*x+y*y<=1)
            in++;
    }
    pai=(double)4*in/NUM;
    cout << "Pi=" << pai << endl;

    return 0;
}
double rnd(void)
{
    return ((double)rand()/32767.1);
}

実行結果

Pi=3.14156