正态分布的随机数发生器 in C#

主要参考《numerical recipes in c++ 2/e》p.292～p.294 和《simulation modeling and analysis 3/e》p.465～p.466。

box 和 muller 在 1958 年给出了由均匀分布的随机变量生成正态分布的随机变量的算法。设 u1, u2 是区间 (0, 1) 上均匀分布的随机变量，且相互独立。令

x1 = sqrt(-2*log(u1)) * cos(2*pi*u2);
x2 = sqrt(-2*log(u1)) * sin(2*pi*u2);

那么 x1, x2 服从 n(0,1) 分布，且相互独立。等于说我们用两个独立的 u(0,1) 随机数得到了两个独立的 n(0,1)随机数。

marsaglia 和 bray 在 1964 年提出了一种改进算法，避免使用三角函数。以下的实现代码用的就是这种改进算法。

//
// gaussian random number generator class
// ref. ``numerical recipes in c++ 2/e'', p.293 ~ p.294
//
  public class gaussianrng
  {
    int iset;
    double gset;
    random r1, r2;
   
    public gaussianrng()
    {
      r1 = new random(unchecked((int)datetime.now.ticks));
      r2 = new random(~unchecked((int)datetime.now.ticks));
      iset = 0;
    }
   
    public double next()
    {
      double fac, rsq, v1, v2;   
      if (iset == 0) {
        do {
          v1 = 2.0 * r1.nextdouble() - 1.0;
          v2 = 2.0 * r2.nextdouble() - 1.0;
          rsq = v1*v1 + v2*v2;
        } while (rsq >= 1.0 || rsq == 0.0);
       
        fac = math.sqrt(-2.0*math.log(rsq)/rsq);
        gset = v1*fac;
        iset = 1;
        return v2*fac;
      } else {
        iset = 0;
        return gset;
      }
    }
  }