2013-08-17 48 views
1

我試圖通過使用隨機性作爲概率替換的簡單方法來創建Buffoon的Needle實驗。 pi的值可以從等式 pi = 2 * ln/th中找到,其中l =針的長度,n =針落下的次數,t =線的寬度,h =針穿過線的次數。 我假設l = t,從而將我的等式減少到pi = 2 * n/h。 現在我已經做了兩個代碼。 代碼1:這是關於蟒蛇的隨機選擇

import math, random 
h = 0.0 
n = 0.0 
for i in range(0,10000): 
    a = random.random() 
    if a > 0.64: 
     h = h+1 
    else: 
     n = n+1 
re = 2*(n+h)/n 
print "Value of pi is ", re 
err = (math.pi - re)*100/(math.pi) 
print "Percentage error is ", abs(err) 

現在這一個運行良好,並給了我足夠好的結果。 但是,下面的代碼一遍又一遍地重複相同的答案。 代碼2:

import random, time, math 
h=1.0 
n=1.0 
err = 0.0 
while err < 0.1: 
    a = random.random() 
    if a > 0.64: 
     h = h+1 
    else: 
     n = n+1 
    re = 2*(n+h)/n 
    err = (math.pi - re)*100/(math.pi) 
print "Number of attempts is ", n+h 

有人能告訴我爲什麼嗎?

+1

未之間均勻地產生的值,它回答您的問題,但你的程序,在這兩種情況下,接近於2/0.64,不PI。這與布馮的針沒有關係。 –

回答

0

您的while循環的條件是倒退。它應該是:

while err > 0.1: 

您還應該設置err1.0或東西比0.1高。

最後,計算錯誤時,你應該使用abs

err = abs(math.pi - re)/(math.pi) 
+0

嗨攪拌機和Antti Haapala感謝您的建議。代碼現在工作正常。我有更多的問題,我將很快發佈。再次感謝。 –

0

你想要的誤差爲儘可能大的開始,然後循環,直到它變得足夠小。使用

import sys 

error = sys.float_info.max 
epsilon = 0.1 
while err > epsilon: 
    ... 

另外,通常差異也可能是負數,所以要獲得錯誤,您希望使用abs來僅比較幅度。因此,我們得到一個更地道程序

import random, time, math 
import sys 

h = n = 1.0 
err = sys.float_info.max 
epsilon = 0.001 

while err > epsilon: 
    a = random.random() 
    if a > 0.64: 
     h += 1 
    else: 
     n += 1 

    re = 2 * (n + h)/n 
    err = abs((math.pi - re)/math.pi) 

print "Value of pi is ", re 
print "Number of attempts is ", n + h 
+0

好的。我會試試看。謝謝。 –

+1

不需要'sys.float_info。max'只需使用'float('+ inf')'。 – Bakuriu

+0

...雖然他們不一樣,但在這裏它的作品。 –

0

其他人指出你停止條件的問題,但你的代碼有一個更根本的缺陷。這是不是布豐的針實驗。你已經找到了一個稍微複雜的方法來計算比率2.0/0.64。如果您重複運行程序並對結果取平均值,您會發現平均值收斂到3.125,而不是Pi。布馮的實驗是基於針的中心位置和落下角度是否穿過一組平行線中的一根,這兩者都不在您的程序中考慮。

如果你想要一個相對簡單的蒙特卡洛方法來估計Pi,在以(0,0)爲中心的2x2平方上隨機生成點,並且看看這些點的哪些比例落入單位半徑的內切圓內。這個比例應該與這兩個區域的比例相同,因此,將估計的比例乘以4得到一個體面的Pi估計。在僞碼:

count = 0 
repeat N times 
    x <- U(-1,1) 
    y <- U(-1,1) 
    if x^2 + y^2 <= 1 
    increment count 
pi_est <- 4 * count/N 

U(-1,1)其中指示-1和1