この記事では VBA でモンテカルロ法(Monte Carlo method)を実装して、円周率の近似値や定積分を求める方法を解説します。
【VBA】モンテカルロ法による円周率の計算
下の図のように半径
そして、この正方形の中に無作為に点を
なので、円周率は
と計算できるはずです。
'[VBA] モンテカルロ法による円周率の近似計算 Sub Monte_Carlo_Pi() Dim x As Double, y As Double Dim r As Double, p As Double Dim n As Long, nc As Long Randomize For n = 1 To 100000 '無作為な点(x,y)を1つ生成 x = RND y = RND '原点からの距離の2乗を計算 r = x ^ 2 + y ^ 2 '原点からの距離が 1 以内ならば '四分円に入ったカウント数に1を加える If r <= 1 Then nc = nc + 1 End If Next n '円周率の近似値を計算 p = 4 * nc / n Debug.Print p End Sub
Monte_Carlo_Pi() を実行するとイミディエイトウィンドウに円周率の近似値が表示されます。乱数を使っているので、実行する度に違った値となります。
Monte_Carlo_Pi() の構造について簡単に解説します。
For n = 1 To 100000
としてあるので、上限値を 100,000 としてありますが、より高い精度を求めるのであれば、この上限値を書き換えてください。x と y に 0 から 1 の乱数を入れていますが、これが発生した点の座標 (x, y) となっています。
r <= x ^ 2 + y ^ 2
という記述で (x, y) の原点からの距離の2乗を変数 r に入れます。
If r <= 1 Then
は r が 1 以内であること、つまり点が四分円の中にあるという条件で、それが満たされた場合に限って
nc = nc + 1
として四分円に入ったカウント数 nc に 1 を加えます。最後に nc と n の比率を使って
p = 4 * nc / n
として円周率の近似値を p に入れています。
【VBA】モンテカルロ法による定積分の計算
関数
をモンテカルロ法で計算してみます。
図にある面積
に入った回数が
というように近似計算できます。モンテカルロ法による定積分の計算プログラムの一例を載せておきます。
'[VBA] モンテカルロ法による定積分の計算 Sub MonteCarlo_Integral() Dim x As Double, y As Double Dim a As Double, itg As Double Dim n As Long, k As Long, ct As Long '積分の上限値 a = 3 Randomize '試行回数の設定 n = 1000 For k = 1 To n '無作為点(x,y)を生成 x = a * RND y = a ^ 2 * RND 'y < f(x) ならば If y < x ^ 2 Then 'カウント数に 1 を加える ct = ct + 1 End If Next k '積分近似値を計算します itg = ct * a ^ 3 / n Debug.Print itg End Sub
真値 (9) と比較しやすいように、上限を
エクセルや数学に関するコメントをお寄せください