为 64 位 LCG 找到更多独立的种子值(MMIX(由 Knuth 提供)) [英] find more indipendent seed value for a 64 bit LCG (MMIX (by Knuth))

问题描述

我使用的是 64 位 LCG(MMIX(由 Knuth)).它在我的代码中生成特定的随机数块，使用它们来执行一些操作.我的代码在单核中工作，我想并行化工作以减少执行时间.
在开始考虑这种意义上的更高级方法之前，我想简单地并行执行更多相同的代码(实际上，代码在一定数量的独立模拟中重复相同的任务，因此我可以简单地将模拟数量分成更多相同的代码并并行运行).
我现在唯一的问题是为每个代码找到一个种子；特别是，为了避免在不同代码中生成的数据之间出现不需要的非平凡相关性的可能性，我必须确保在各种代码中生成的随机数不会重叠.为此，从第一个代码中的某个种子开始，我必须找到一种方法来找到一个非常远的值(下一个种子)，不是绝对值而是伪随机序列(所以，这样，从第一个到第二个种子，我需要大量的 LCG 步骤).
我的第一次尝试是这样的:
从序列中生成的2个连续数之间的LCG关系开始

所以，原则上，我可以用 n = 2^40 和 I_0 等于第一个种子的值来计算上述关系，并从第一个随机 CLG 序列中获得一个距离 2^40 步的新种子.
问题是，这样做，我必须进行溢出计算 a^n.事实上，对于 MMIX(由 Knuth 提供)a~2^62，我使用 unsigned long long int(64 位).请注意，这里唯一的问题是上述关系中的分数.如果只有和和乘法，我可以忽略由于以下模块化属性引起的溢出问题(实际上我使用 2^64 作为 c(64 位生成器)):

I'm using a 64 bit LCG (MMIX (by Knuth)). It generate a certain block of random numbers inside my code, which use them to perform some operations. My code works in single core and I would like to parallelize the work to reduce the execution time.
Before start thinking to more advanced methods in this sense I'd like to simply execute more identical codes in parallel (in fact the code repeats the same task over a certain numbers of indipendent simulation, so I can simply split the number of simulation between more identical codes and run them in parallel).
My only problem now is to find a seed for each code; in particular, to avoid the possibility of unwanted non trivial correlation between data generated in different codes, I have to be sure that the random number generated in the various codes don't overlap. To do so, starting from a certain seed in the first code I have to find a way to find a value (the next seed) very distant not in absolute value but in the pseudo-random sequence (so, such that, to go from the first to the second seed, I need a huge number of steps of LCG).
My first attempt was this:
starting from the LCG relation between 2 consecutive numbers generated in the sequence

So, in principle, I could calculate the above relation with, say, n = 2^40 and I_0 equal to the value of the first seed, and obtain a new seed distant 2^40 steps in the random CLG sequence from the first one.
The problem is that, doing so, I necessary go in overflow calculating a^n. In fact for MMIX (by Knuth) a~2^62 and i use unsigned long long int (64 bit). Note that the only problem here is the fraction in the above relation. If there only were sum and multiplication I could ignore the overflow problem due to the following modular properties (in fact I'm using 2^64 as c (64 bit generator)):

那么，从某个值(第一个种子)开始，如何在 LC 伪随机序列中找到距离大量步长的第二个值?

So, starting from a certain value (first seed), how can I find a second one distant a huge number of step in the LC pseudo-random sequence?

[编辑]
r3mainer 解决方案非常适合 Python 代码.我现在正在尝试使用 unsigned __int128 变量在 c 中实现它.我只有一个问题:原则上我应该计算:

[EDIT]
r3mainer solution is perfectly suited for python codes. I'm trying now to implement it in c using unsigned __int128 variables. I have only one problem: in principle I should compute:

比如说，为了简单起见，我想计算:

Say, for simplicity, I want to compute:

n = 2^40 和 c(a-1)~2^126.我继续一个循环.从 temp = a 开始，在每次迭代中我计算 temp = temp*temp，然后我计算 temp%c(a-1).问题出在第二步 (temp = temp*temp).temp 实际上可以是，原则上任何数字 temp 是一个大数字，请说 >2^64，我会溢出，在下一个模块操作之前达到2^128 - 1.那么有没有办法避免呢?目前，我看到的唯一解决方案是使用位循环执行每个乘法，如下所示:c 代码:防止在具有巨大模块(接近溢出阈值的模块)的模块化操作中溢出在乘法过程中还有其他方法可以进行模运算吗?
(请注意，作为 c = 2^64，使用 mod(c) 操作我没有同样的问题，因为溢出点(对于 ull int 变量)与模块重合)

with n = 2^40 and c(a-1)~2^126. I proceed with a cycle.Starting with temp = a, in each iteration I compute temp = temp*temp, then I compute temp%c(a-1). The problem is in the second step (temp = temp*temp). temp in fact could be, in principle any number < c(a-1)~2^126. If temp is a big number, say > 2^64, I'll go in overflow, reaching 2^128 - 1, before the next module operation. So is there a way to avoid it? For now the only solution I see is to perform each multiplication with a loop over bit, as suggested here: c code: prevent overflow in modular operation with huge modules (modules near the overflow treshold) Is there another way to perform module operation during the multiplication?
(note that being c = 2^64, with mod(c) operation I don't have the same problem because the overflow point (for ull int variables) coincides with the module)

推荐答案

x[n+1] = (x[n] * a + c) % m 形式的任何 LCG 都可以很快跳到任意位置.

Any LCG of the form x[n+1] = (x[n] * a + c) % m can be skipped to an arbitrary position very quickly.

从零种子值开始，LCG 的前几次迭代将为您提供以下序列:

Starting with a seed value of zero, the first few iterations of the LCG will give you this sequence:

x₀ = 0
x₁ = c % m
x₂ = (c(a + 1)) % m
x₃ = (c(a² + a + 1)) % m
x₄ = (c(a³ + a² + a + 1)) % m

很容易看出每个项实际上是几何级数的总和，可以使用 简单公式:

It's pretty easy to see that each term is actually the sum of a geometric series, which can be calculated with a simple formula:

x_n = (c(a^{n-1} + a^{n-2} + ... + a + 1)) % m
    = (c * (a^n - 1) / (a - 1)) % m

(a^n - 1) 项可以通过 快速计算模幂，但除以(a-1) 有点棘手，因为(a-1) 和m 都是甚至(即，不是互质)，所以我们无法计算 模乘逆的 模乘逆code>(a-1) mod m 直接.

The (a^n - 1) term can be calculated quickly by modular exponentiation, but dividing by (a-1) is a bit tricky because (a-1) and m are both even (i.e., not coprime), so we can't calculate the modular multiplicative inverse of (a-1) mod m directly.

相反，计算(a^n-1) mod m*(a-1)，然后将结果直接(非模)除以a-1.在 Python 中，计算过程如下:

Instead, calculate (a^n-1) mod m*(a-1), then perform a straightforward (non-modular) division of the result by a-1. In Python, the calculation would go something like this:

def lcg_skip(m, a, c, n):
    # Calculate nth term of LCG sequence with parameters m (modulus),
    # a (multiplier) and c (increment), assuming an initial seed of zero
    a1 = a - 1
    t = pow(a, n, m * a1) - 1
    t = (t * c // a1) % m
    return t

def test(nsteps):
    m = 2**64
    a = 6364136223846793005
    c = 1442695040888963407
    #
    print("Calculating by brute force:")
    seed = 0
    for i in range(nsteps):
        seed = (seed * a + c) % m
    print(seed)
    #
    print("Calculating by fast method:")
    # Calculate nth term by modular exponentiation
    print(lcg_skip(m, a, c, nsteps))

test(1000000)

因此要创建具有非重叠输出序列的 LCG，您需要做的就是使用由 lcg_skip() 生成的初始种子值和 n 的值距离足够远.

So to create LCGs with non-overlapping output sequences, all you would need to do is use initial seed values generated by lcg_skip() with values of n that are far enough apart.

这篇关于为 64 位 LCG 找到更多独立的种子值(MMIX(由 Knuth 提供))的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！