モンテカルロシミュレーションで出したπの評価

やること

古典的なモンテカルロシミュレーションの訓練である、円と乱数を使ったシミュレーションでπを出して、評価する。

モンテカルロシミュレーションとは

確率を用いて、値を推定する方法。推定とは言っても、十分に信頼できる値を手軽に出せることから、モンテカルロ法の有用性は高い。現在でも金融商品のリスク計算や統計学における数値計算などで用いられている。

πを出す方法（条件・式）

0~1の乱数を元に、0~1のxy平面にランダムに点を打つ。円内の点の数Pと全体の点の数Nの比が、円S1と平面全体S2の面積比に対応している。

ランダムな点が円の中に入っている様子が分かる。図の単位円は1/4されているので、面積比はS1/S2=π/4である。円内の点の検出には原点からの距離を用いた。距離が1以下の場合は変数pに1を加算する。距離は円の領域x^2+y^2<=1により求まる。プログラム言語はPythonを用いた。

コード

MCsim

import random
def main():
  sta,order=map(int,input("start order").split()) #何個の点を打つか決める。
  for order in range(order):
    N=sta*10**order #orderは指数
    p=0
    for n in range(N):
      x=make_r()
      y=make_r()
      if x**2+y**2<=1:
        p+=1
    pi=4*p/N
    print("pi=",pi,"N=",N)

#0以上1以下の乱数作り関数
def make_r():
  while True:
    r=(random.random())*2
    if r>1:
      continue
    else:
      break
  return r
  
if __name__=="__main__":
  main()

結果

結果は以下のグラフのようになった。

点の数が増えるにつれてpi(π)の値に近づいている。1000万の点を打った時、πは3.1409824になった。

評価

ランダムに打つ点の数が1万を超えたあたりから小数点第二位までの精度で見れるようになった。1000万の点の場合、四捨五入すると小数点第三位まで一致している。

考察

十分なサンプルの点があれば、モンテカルロシミュレーションでは十分に精度高く値を推定することができる。乱数の生成を組み込み関数に頼ってしまったため、次は乱数にもこだわりたい。