algorithm - 最小标准随机生成器基础知识

Question

最小标准随机数生成器。

我正在阅读有关最低标准数字生成器的信息，如下所示

Given a random integer xn, the next random integer in a random sequence is given by computing xn+1 = a xn (mod m), where a = 7 ^ 5 = 16807 and m = 2^31 − 1 = 2147483647; as a check on your implementation, if x0 = 1, then x10000 = 1043618065.

Park and Miller chose m as the largest Mersenne prime less than 2^32; the smallest primitive root of m is 7, and since 5 is also prime, 7^5 is also a primitive root, hence their choice of a. Because a is a primitive root of m, all values in the range 1 to m − 1 inclusive will be generated before any repeat, so the random number generator has full period. The multiplier a = 16807 has been shown to have good randomness properties.

Subsequent to their original paper, Park and Miller recommended 48271 as an improvement, and some people use 69621, but we’ll continue to use 16807. The easiest way to implement that is obvious: just multiply a by the current value of x and compute the modulus.

But that may cause overflow in the intermediate multiplication, rendering the results incorrect.

A trick of Linus Schrage allows that multiplication to be done without overflow: Compute q = ⌊m / a⌋ and r = m (mod a) so that m = a q + r. Then

a new x can be computed by hi = ⌊x / q⌋, lo = x (mod q), x = a · lo − r · hi, then adding m to x if x ≤ 0.

我的问题是作者如何根据 hi = floor(x/q) 和 lo = x(modq) 计算新的 x。我在这里寻找步骤。请解释。

Answer 1

让我们简化符号。设置 H = hi 和 L = lo。我们有 m = a * q + r。一个简单的计算表明 q = 127773 和 r = 2836。我们观察到 a < q。

现在给定 x_{n} 并计算 H = x_{n}/q 和 L = x_{n} % q。所以，x_{n} = q * H + L 其中 L < q。

根据定义 x_{n+1} = a * x_{n} mod m。计算右侧(在减少模 m 之前)我们得到 a * x_{n} = a * (q * H + L) = a * q * H + a * L = (m - r) * H + a * L = m * H - r * H + a * L。

现在让我们考虑 r * H。显然 0 <= r * H < a * (x_{n}/q)。由于 x_{n} < m 并且如上观察到的 a < q，a * (x_{n}/q) < m。特别是，它不会溢出。

同样 0 < a * L < a * q < m。所以，再次没有溢出。

我们得出结论，x_{n+1} = m * H - r * H + a * L。将后者对模 m 进行约化，我们得到 x_{n+1} = -r * H + a * L，两者均不存在右边的两个表达式溢出 m。如果总和为负，我们加上 m 就完成了。