随机数生成器：原理、类型与应用

### 随机数生成器：原理、类型与应用 #### 1. 引言 用计算机生成“随机”数看似有些违背常理，毕竟计算机是人类设计出的最精确、最具确定性的机器。任何程序的输出都是可预测的，因此并非真正意义上的“随机”。然而，实用的计算机“随机数生成器”却被广泛使用。 在计算机生成序列的语境中，一个虽不精确但可行的随机性定义是：生成随机序列的确定性程序应与使用该输出的计算机程序不同，并且在所有可测量的方面都统计不相关。也就是说，当两个不同的随机数生成器与特定应用程序结合使用时，它们应该产生统计上相同的结果。若不满足这一点，那么至少其中一个生成器就不是好的生成器。 从实用角度看，随机性取决于使用者（或程序员）的判断。对于一个应用来说足够随机的数，对另一个应用可能就不够了。不过，有一系列统计测试可以帮助找出应用程序可能检测到的相关性。好的随机数生成器应该通过所有这些测试，或者用户至少要了解它们未通过的测试，以便判断这些测试是否与当前情况相关。 #### 2. 均匀偏差 均匀偏差是指位于指定范围内（通常是 0 到 1）的随机数，该范围内的每个数字出现的可能性相同。它是随机数的一种常见形式，也是许多其他类型随机数生成的基础。 ##### 2.1 系统提供的随机数生成器 大多数 C 语言实现都有用于初始化和生成“随机数”的库函数。在 ANSI C 中，相关函数如下： ```c #include <stdlib.h> #define RAND_MAX ... void srand(unsigned seed); int rand(void); ``` 使用 `srand(seed)` 函数并传入一个任意种子来初始化随机数生成器。每个初始值通常会产生不同的随机序列，或者至少是某个极长序列的不同起始点。相同的种子值总是会返回相同的随机序列。通过连续调用 `rand()` 函数可以获取序列中的连续随机数，该函数返回的整数通常在 0 到 `int` 类型可表示的最大正值（包含）之间。若要获取 0.0（包含）到 1.0（不包含）之间的随机浮点数，可以使用以下表达式： ```c x = rand()/(RAND_MAX+1.0); ``` 不过，要对系统提供的类似 `rand()` 函数保持高度怀疑。系统提供的 `rand()` 函数几乎都是线性同余生成器，通过递推关系 `Ij+1 = aIj + c (mod m)` 生成整数序列。虽然这种方法理论上可以提供不错的随机数，但许多 ANSI C 库的实现存在缺陷。 常见问题包括： - **返回值范围小**：ANSI 标准规定 `rand()` 返回 `int` 类型的值，在许多机器上这只是一个两字节的量，导致 `RAND_MAX` 通常不大。这在蒙特卡罗积分等情况下可能会造成严重问题。 - **实现不佳**：ANSI 委员会未规定标准算法，一些实现者试图改进示例代码，却弄巧成拙。例如，一个流行的 32 位 PC 兼容编译器交换了返回值的高低 16 位，破坏了生成器的随机性。 此外，线性同余方法虽然速度快，但存在连续调用的序列相关性问题。如果用 k 个随机数在 k 维空间中绘制点，这些点不会均匀填充 k 维空间，而是会位于 (k - 1) 维的“平面”上。而且，线性同余生成器的低阶（最低有效）位通常不如高阶位随机。 ##### 2.2 便携式随机数生成器 Park 和 Miller 提出了一种“最小标准”生成器，基于简单的乘法同余算法 `Ij+1 = aIj (mod m)`，其中 `a = 16807`，`m = 231 - 1 = 2147483647`。该生成器自 1969 年提出后，通过了所有新的理论测试，并积累了大量成功应用案例。 由于 `a` 和 `m - 1` 的乘积超过了 32 位整数的最大值，不能直接用高级语言实现该算法。Schrage 提出了一种技巧，用于在不使用大于 32 位中间值的情况下，对两个 32 位整数进行模 32 位常数的乘法运算。以下是“最小标准”生成器的实现： ```c #define IA 16807 #define IM 2147483647 #define AM (1.0/IM) #define IQ 127773 #define IR 2836 #define MASK 123459876 float ran0(long *idum) { long k; float ans; *idum ^= MASK; k=(*idum)/IQ; *idum=IA*(*idum-k*IQ)-IR*k; if (*idum < 0) *idum += IM; ans=AM*(*idum); *idum ^= MASK; return ans; } ``` `ran0` 的周期为 `231 - 2 ≈ 2.1 × 109`。不过，该生成器存在一些问题，例如连续随机数之间存在细微的序列相关性。为了消除这些低阶序列相关性，`ran1` 函数在 `ran0` 的基础上进行了洗牌操作： ```c #define IA 16807 #define IM 2147483647 #define AM (1.0/IM) #define IQ 127773 #define IR 2836 #define NTAB 32 #define NDIV (1+(IM-1)/NTAB) #define EPS 1.2e-7 #define RNMX (1.0-EPS) float ran1(long *idum) { int j; long k; static long iy=0; static long iv[NTAB]; float temp; if (*idum <= 0 || !iy) { if (-(*idum) < 1) *idum=1; else *idum = -(*idum); for (j=NTAB+7;j>=0;j--) { k=(*idum)/IQ; *idum=IA*(*idum-k*IQ)-IR*k; if (*idum < 0) *idum += IM; if (j < NTAB) iv[j] = *idum; } iy=iv[0]; } k=(*idum)/IQ; *idum=IA*(*idum-k*IQ)-IR*k; if (*idum < 0) *i ```