在C中實現Box-Mueller隨機數發生器＃

From this question: Random number generator which gravitates numbers to any given number in range?自從我遇到過這樣一個隨機數發生器之前，我做了一些研究。我只記得是名爲「穆勒」，所以我想我找到了，在這裏：在C中實現Box-Mueller隨機數發生器＃

Box-Mueller transform

我可以找到其他語言的它衆多的實現，但我似乎無法到在C＃中正確實現它。

此頁面，例如，The Box-Muller Method for Generating Gaussian Random Numbers說，代碼應該是這樣的（這不是C＃）：

#include <stdio.h> 
#include <stdlib.h> 
#include <math.h> 
#include <time.h> 

double 
gaussian(void) 
{ 
    static double v, fac; 
    static int phase = 0; 
    double S, Z, U1, U2, u; 

    if (phase) 
     Z = v * fac; 
    else 
    { 
     do 
     { 
     U1 = (double)rand()/RAND_MAX; 
     U2 = (double)rand()/RAND_MAX; 

     u = 2. * U1 - 1.; 
     v = 2. * U2 - 1.; 
     S = u * u + v * v; 
     } while(S >= 1); 

     fac = sqrt (-2. * log(S)/S); 
     Z = u * fac; 
    } 

    phase = 1 - phase; 

    return Z; 
}

現在，這是我在C＃實現以上。請注意，變換會產生2個數字，因此上面的「階段」是個訣竅。我只是丟棄第二個值並返回第一個值。

public static double NextGaussianDouble(this Random r) 
{ 
    double u, v, S; 

    do 
    { 
     u = 2.0 * r.NextDouble() - 1.0; 
     v = 2.0 * r.NextDouble() - 1.0; 
     S = u * u + v * v; 
    } 
    while (S >= 1.0); 

    double fac = Math.Sqrt(-2.0 * Math.Log(S)/S); 
    return u * fac; 
}

我的問題是，具有以下情形，在我的代碼不會在0-1的範圍內返回一個值，我無法理解的原代碼如何既可以。

U = 0.5，V = 0.1
S成爲0.5*0.5 + 0.1*0.1 = 0.26
FAC變得〜3.22
返回值是這樣〜0.5 * 3.22或〜1.6

即內不0 .. 1。

我在做什麼錯誤/不理解？

如果我修改我的代碼，以便代替乘fac與u，我通過S乘，我得到的範圍從0到1的值，但它有錯誤的分佈（好像有周圍具有極大分佈0.7- 0.8，然後在兩個方向上逐漸減小。）

來源

2011-04-28 Lasse Vågsæther Karlsen

請注意，我已經檢查了上述代碼的兩個示例，通常是C或Java，它們看起來都差不多。 – 2011-04-28 11:12:23

你確定C代碼生成你想要的？ – Euphoric 2011-04-28 11:45:13

你的代碼沒問題。你的錯誤是認爲它應該只在[0, 1]內返回值。（標準）正態分佈是在整個實線上具有非零權重的分佈。也就是說，[0, 1]以外的值是可能的。實際上，[-1, 0]內的值與[0, 1]內的值一樣可能，此外，補充的[0, 1]約佔正態分佈權重的66％。因此，66％的時間我們預計[0, 1]以外的值。

此外，我認爲這不是Box-Mueller變換，而是實際上是Marsaglia極座標法。

來源

2011-04-28 16:20:01 jason

嗯，在文獻中他們只是稱之爲Box-Muller極座標形式，參見http://en.wikipedia.org/wiki/Box%E2%80%93Muller_transform#Polar_form – 2013-04-17 06:03:09

我認爲該函數返回極座標。所以你需要兩個值才能得到正確的結果。

此外，高斯分佈不在0 .. 1之間。它可以很容易地結束爲1000，但這種情況發生的可能性非常低。

來源

2011-04-28 11:49:43 Euphoric

我不是數學家，也不是統計學家，但是如果我想到這個問題，我不會指望高斯分佈在確切範圍內返回數字。根據你的實現，平均值爲0，標準偏差爲1，所以我期望在鐘形曲線上分佈的值在0的中心，然後隨着數字偏離0而減小。所以序列肯定會覆蓋兩個+/-數字。

然後，因爲它是統計數據，爲什麼僅僅因爲std.dev是1而難以限制爲-1..1？在統計上可以有一些發揮，並且仍然滿足統計要求。

來源

2011-04-28 11:55:45

均勻隨機變量確實在0..1之內，但高斯隨機變量（這是Box-Muller算法生成的）可以在實線上的任何位置。詳情請參閱wiki/NormalDistribution。

來源

2011-04-28 16:11:38 Sasha

這是一個蒙特卡羅方法，所以你不能限制結果，但你可以做的是忽略樣本。

// return random value in the range [0,1]. 
double gaussian_random() 
{ 
    double sigma = 1.0/8.0; // or whatever works. 
    while (1) { 
     double z = gaussian() * sigma + 0.5; 
     if (z >= 0.0 && z <= 1.0) 
      return z; 
    } 
}

來源

2015-02-26 12:29:44

在C中實現Box-Mueller隨機數發生器＃

回答

相關問題