2013-05-08 60 views
0

我試圖生成使用蒙特卡洛方法從一個正態分佈值,按照該網站http://math60082.blogspot.ca/2013/03/c-coding-random-numbers-and-monte-carlo.html沒有足夠的隨機蒙特卡羅

我修改了代碼從原來的一個位,所以它計算出方差和平均對於直接生成的數字來檢查方法是否正在工作,而不是單獨進行測試(同樣的差異真的,但只是一個正面)。

問題不管我做什麼

,方差爲遠高於1,平均不爲零。是否有可能產生的僞隨機數不夠隨機?

代碼

請注意,上面給出的網站的作者是人誰寫的代碼

#include <cstdlib> 
#include <cmath> 
#include <ctime> 
#include <iostream> 
using namespace std; 
// return a uniformly distributed random number 
double uniformRandom() 
{ 
    return ((double)(rand()) + 1.)/((double)(RAND_MAX) + 1.); 
} 

// return a normally distributed random number 
double normalRandom() 
{ 
    double u1=uniformRandom(); 
    double u2=uniformRandom(); 
    return cos(8.*atan(1.)*u2)*sqrt(-2.*log(u1)); 
} 

int main() 
{ 
    double z; 
    int N=1000; 
    double array[N]; 
    double mean=0 ,variance=0; 
    srand(time(NULL)); 

    for(int i=0;i<N;i++) 
    { 
     z=normalRandom(); 
     cout << i << "->"<< z<< endl; 
     mean+=z; 
     array[i]=z; 
    } 

    mean=mean/N ; 
    cout << " mean = " << mean << endl; 

    for(int i=0;i<N;i++) 
    { 
     variance = variance + (mean - array[i])*(mean - array[i]); 
    } 
    variance = variance/N; 
    cout << " variance = " << variance << endl; 

    return 0; 
} 

UPDATE

顯然受到用戶的指出,我搞砸了,因爲一個非常愚蠢的錯誤,程序不工作。

+0

我以前讀過從統一隨機生成正態分佈隨機的正態分佈通常不是特別準確 - 最好使用正態分佈RNG。 (不知道在哪裏可以找到,但是不知道C++。) – 2013-05-08 00:41:28

+0

@HotLicks:我不知道你在哪裏閱讀,但這是一個非常不穩定的說法。生成正態分佈形成一個統一的分佈是相當普遍的(http://en.wikipedia.org/wiki/Marsaglia_polar_method)[thing](http://en.wikipedia.org/wiki/Box%E2%80 %93Muller_transform)。 *然而*,使用一個可憐的統一生成器(比如'rand()')會導致很差的正態分佈。 – 2013-05-08 00:42:52

+0

@HotLicks注意到,顯然是我這樣做的朋友,希望它成爲蒙特卡洛方法。 – 2013-05-08 00:44:19

回答

3

您似乎以錯誤的方式計算了meanmean應該在N之間進行平均,而您只能求和所有數組元素。目前mean實際上是sum

mean = mean /N 
+0

感謝您指出了這一點.B方差仍然沒有接近1. – 2013-05-08 00:38:04

+2

@nerorevenge您可能會嘗試生成超過2000個rand數並查看方差是否接近1?rand()可能不是一個好的隨機數生成器(只是我瘋狂的猜測)。 – taocp 2013-05-08 00:40:38

+1

@nerorevenge:你錯了。平均值的不正確值會影響方差的值。修復平均值後,修復代碼中的所有問題。 – carlosdc 2013-05-08 00:46:26

2

rand()是大多數實現中的質量非常低的隨機數生成器。有些Linux版本會從內核熵池中獲取價值,但不能保證跨平臺(例如,在Windows上)使用Mersenne Twister。 Boost庫實現了一個。

編輯:taocp答案突出顯示了一個編碼問題,但RNG問題仍然適用。

+0

注意到了,代碼正在linux機器上運行,所以你的答案可能是正確的,讓我給你gback。 – 2013-05-08 00:40:43

+0

Gback?無論如何,並不是所有的Linux libc都使用熵池。我無法分辨哪些是做什麼,哪些不做,但我非常確定'rand'對於在PC上播放「猜測我的號碼」非常有用,而不是進行統計分析或密碼學。 – 2013-05-08 00:45:08

+1

如果您使用的是最新版本的gcc,則不需要增加此功能。 ''現在具有所有這些功能。 – Yuushi 2013-05-08 02:20:19