用C#实现Box-Mueller随机数发生器

Question

从这个问题:随机数发生器将数字吸引到范围内的任何给定数字？我之前做过一些研究,因为我之前遇到过这样一个随机数发生器.我记得的只是"穆勒"的名字,所以我想我找到了,在这里:

• Box-Mueller变换

我可以在其他语言中找到它的大量实现,但我似乎无法在C#中正确实现它.

例如,这个页面,用于生成高斯随机数的Box-Muller方法表示代码看起来应该是这样的(这不是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具有极大分布,然后在两者逐渐减弱方向).

Answer 1

你的代码很好.你的错误是认为它应该只返回值[0, 1].(标准)正态分布是在整个实线上具有非零权重的分布.也就是说,外面的值[0, 1]是可能的.事实上,其中的值[-1, 0]与其中的值一样可能[0, 1],而且,补码[0, 1]约为正态分布权重的66%.因此,66%的时间我们期望一个值超出[0, 1].

另外,我认为这不是Box-Mueller变换,而实际上是Marsaglia极地方法.