均勻分佈

Numerical Recipes, The Art of Scientific Computing, Third Edition, Chapter 7 Random Numbers:

Ran 程式

typedef unsigned int Uint;             // 32-bit unsinged integer
typedef unsigned long long int Ullong; // 64-bit unsinged integer
typedef double Doub;                   // 64-bit floating point

注意，由於 Ullong 和 Doub 都是 64-bit 的，所以把 Ullong 轉換為 Doub 時有可能喪失精度。即，上述程式中的 int64() 的不同的值有可能對應 doub() 的相同的值。

1.0 / (1.0 + ULLONG_MAX) == 1.0 / ULLONG_MAX == 5.42101086242752E-20
0.0 <= doub() <= 1.0  // 注意: doub() 有可能等於 1.0

NrRandom 類

使用上述 Ran 程式中的演算法，可以得到以下 C# 程式：

namespace Skyiv {
  using System;

  public sealed class NrRandom : Random
  {
    ulong v1, v8, v16, v32, u, v = 4101842887655102017, w = 1;
    int c1 = 0, c8 = 0, c16 = 0, c32 = 0;

    protected override double Sample() {
      return 1 / (1.0 + uint.MaxValue) * NextUInt32(); }
    public override double NextDouble() { return Sample(); }
    public override int Next() { return Next(int.MaxValue); }
    public NrRandom() : this((ulong)DateTime.Now.Ticks) { }

    public NrRandom(ulong seed)
    {
      u = v ^ seed; NextUInt64();
      v = u;        NextUInt64();
      w = v;        NextUInt64();
    }

    public ulong NextUInt64()
    {
      u = u * 2862933555777941757 + 7046029254386353087;
      v ^= v >> 17; v ^= v << 31; v ^= v >> 8;
      w = 4294957665U * (w & 0xffffffff) + (w >> 32);
      ulong x = u ^ (u << 21); x ^= x >> 35; x ^= x << 4;
      return (x + v) ^ w;
    }

    public override int Next(int max)
    {
      if (max < 0) throw new ArgumentOutOfRangeException(nameof(max));
      return (int)(Sample() * max);
    }

    public override int Next(int min, int max)
    {
      if (min > max) throw new ArgumentOutOfRangeException(nameof(min));
      return (int)(min + (long)(Sample() * (max - (long)min)));
    }

    public override void NextBytes(byte [] buffer)
    {
      if (buffer == null) throw new ArgumentNullException(nameof(buffer));
      for (int i = 0; i < buffer.Length; i++) buffer[i] = NextByte();
    }

    public uint NextUInt32()
    {
      if (c32-- != 0) return (uint)(v32 >>= 32);
      c32 = 1; return (uint)(v32 = NextUInt64());
    }

    public ushort NextUInt16()
    {
      if (c16-- != 0) return (ushort)(v16 >>= 16);
      c16 = 3; return (ushort)(v16 = NextUInt64());
    }

    public byte NextByte()
    {
      if (c8-- != 0) return (byte)(v8 >>= 8);
      c8 = 7; return (byte)(v8 = NextUInt64());
    }

    public bool NextBoolean()
    {
      if (c1-- != 0) return ((v1 >>= 1) & 1) != 0;
      c1 = 63; return ((v1 = NextUInt64()) & 1) != 0;
    }
  }
}

簡要說明：

• NrRandom(ulong seed) 對應 Ran(Ullong j)。
• ulong NextUInt64() 對應 Ullong int64()。
• 上述兩個方法是本演算法的核心部分，其他方法都基於這兩者。
• double NextDouble() 要求在 [0.0, 1.0) 內，所以和 Doub doub() 有所不同。

我想，在需要高質量的偽隨機數生成器的場合，這個 NrRandom 類可以替代 .NET Framework 中的 Random 類。而且，它還有以下的方法：

• ulong NextUInt64()
• uint NextUInt32()
• ushort NextUInt16()
• byte NextByte()
• bool NextBoolean()

參考資料

1. MSDN: Random Class (System)
2. 豆瓣：Numerical Recipes, The Art of Scientific Computing, 3rd Edition

---