今更聞けないモンテカルロ法(円周率編)

カジノを潰したという都市伝説を持つモンテカルロ法で、円周率を求めます。
C#ですけれど、ソースコードはこちらにあります。
https://github.com/ma3100/AlgorithmPrimer/blob/master/FirstAlgorithmForCSharp/First/MonteCalro.cs
※随時追加予定

モンテカルロ法とは?

モンテカルロ法は、問題を数値計算で解くのではなく確率を利用して解く事です。

例えば今回テーマにしている円周率の場合、**「円周の長さと直径の比」**から一生懸命計算する訳ですが、モンテカルロ法は上記の計算を必要としません。

そのため、今回学ぶ方法で求める円周率は正確な値ではないのですが、導入としては適しているのでその辺りは割り切ってください。

モンテカルロ法で円周率を求める際の考え方

※半径が1の円で説明します。

1. 円周率を求めたい円を四分割し、直径と同じ辺を持つ正方形で1/4の円を囲みます。
2. 1で出来た範囲の中に一様実数乱数を発生させます。3つしか点を打っていませんが、以下の様な状態になります。

1. 一様実数乱数であれば、「1/4の円の中に収まる乱数と全ての乱数」の比は、「1/4の円の面積と正方形の面積の比」と同じになります。

1/4の円の面積 : 正方形の面積 = 1/4の円の中に収まる乱数の数 : 全ての乱数の数
1/4の円の中に収まる乱数の数をa,円の外に出た乱数をbと置くと、
円周率/4 : 1 = a : a + b
円周率 = 4a / (a + b)
= 4a / n(発生した乱数の総数)

実装(C#)

pymonte.cs

        const int RANDOMNUM = 10000;

        void GetMonteCalro()
        {
            double x, y, montePai;
            var random = new Random();
            var inCircle = 0;
            foreach(var i in Enumerable.Range(0,RANDOMNUM))
            {
                x = random.NextDouble();
                y = random.NextDouble();

                if(x * x + y * y <= 1)
                {
                    inCircle++;
                }
            }

            montePai = (double)4 * inCircle / RANDOMNUM;
            Console.WriteLine($"Answer:{montePai}");
        }

※x^2 + y ^2 が半径を超えないのであれば、円の中にあると見なせます。

例えばRANDOMNUMの値を100などにすると、かなり怪しい結果になるのですがこの値を大きくしていくと、だんだんと3.14に近しい数になっていきます。

今回は以上です、次はモンテカルロ法を用いて円の面積を求めていきます。