モンテカルロ法で円周率を求めてみた

学校の授業でモンテカルロ法で円周率を出す方法(ただし近似値)を

習ったのでpythonで組んでみることにしました。

ダウンロード
PI.zip
zip ( 圧縮 ) ファイル 444 Bytes

#coding: UTF-8
import math,time,random
def MonteCarlo(loop,sample):
    C = 0
    for r in range (loop):
        for i in range (sample):
            E = (random.random()**2+random.random()**2)
            if E <= 1:
                C += 1.0
    return C

random.seed()
    
while True:
    Loop = input("Loop >>>")
    Sample = input("Sample >>>")
    
    Start = time.time()

    result = MonteCarlo(Loop,Sample)
    PI = (result/(Sample*Loop))*4.0
    
    print str(PI)+"\n"+str(math.fabs((PI-math.pi)/math.pi*100))+"%\n"+str(time.time() - Start)
        
    if raw_input("Again? Y/N > ").upper() != "Y":
        break


トータルで1000万回ぐらいやると、3.14が安定して出てきますが、そこから先の桁がばらばらになってしまいますね、やっぱり...。

コメント: 0