集群环境下的伪随机数生成器

Question

我如何在集群上生成独立的伪随机数，例如用于Monte Carlo模拟？我可以有许多计算节点（例如100个），并且我需要在每个节点上生成数百万个数字。我需要保证一个节点上的PRN序列不会与另一个节点上的PRN序列重叠。

• 我可以在根节点生成所有PRN，然后将它们发送到其他节点。但这会非常慢。
• 我可以在每个节点上跳转到序列中的已知距离。但是对于Mersenne-Twister或任何其他好的PRNG，是否存在这样的算法，可以在合理的时间和内存量内完成？
• 我可以在每个节点上使用不同的生成器。但是对于像Mersenne-Twister这样的好生成器是否可能？如何完成？
• 还有其他想法吗？

- Charles Brunet

---

相关：http://cstheory.stackexchange.com/questions/6903/parallel-pseudorandom-number-generators - Jukka Suomela

@Jukka Suomela：也许你应该在cstheory上发布自己的答案，提到我们已经在这里解决了这个问题。 - jopasserat

Answer 1

您不应该使用从相同原始流获得的潜在重叠的随机流。如果您没有测试结果交错的流的统计质量，您就无法知道其统计质量。

幸运的是，梅森旋转算法 (MT) 可以帮助您完成分布任务。使用其专用算法，称为 动态创建器 (DC)，您可以创建 独立的随机数生成器，将产生高度独立的随机流。

每个流将在使用它的节点上创建。基本上，可以将 DC 视为面向对象范例中创建 MT 的构造函数。每个不同的实例都设计为产生高度独立的随机序列。

您可以在此处找到 DC：http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/DC/dc.html
它很容易使用，您将能够修复不同的参数，例如您要获取的不同 MT 实例数量或这些 MT 的周期。根据其输入参数，DC 运行时会发生变化。

除了随 DC 一起提供的 README 外，还可以查看 DC 存档中的文件 example/new_example2.c。它显示了如何调用 给定不同的输入标识符 来获取独立序列，这基本上就是您必须标识集群作业的内容。

最后，如果您打算了解如何在并行或分布式环境中使用 PRNG，请阅读以下科学文章：

高性能计算随机流的实际分布，David RC Hill，在高性能计算和模拟国际会议 (HPCS)，2010年

- jopasserat

---

DCMT是我实际使用的工具。它使用节点编号来创建生成多项式，因此每个随机序列都是独立的。但是有没有证明呢？我记得在DCMT的原始文章中，他们并没有证明它能够正常工作，只是假设它应该可以工作。 - Charles Brunet

我希望有一个数学证明，但不幸的是并没有！我正在攻读关于高性能环境下随机模拟的博士学位，因此我广泛研究了这个问题。基本上，如果你没有内存限制（MT家族状态向量相当大），这种方法是确保序列之间高度独立的最佳方法。其他方法在我引用的文献中进行了调查，但作者支持DCMT的参数化方法。 - jopasserat

TinyMT 可能是一个不错的选择。 - jj1bdx

@max 像这样吗？http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt64.html - jopasserat

在DC中，get_mt_parameter函数的第一个参数允许您指定字长。 - jopasserat

Answer 2

好的，回答 #2 ;-)

我要说...保持简单。只需使用“短”种子来启动MT（假设这个种子是2³²位，为了限制不太好）。这假定短种子生成“足够分布”的MT起始状态（例如，在我的另一个答案中的代码中的init_genrand，希望如此）。这并不能保证等分布的起始状态，而是追求“足够好”，请参见下文。

每个节点将使用自己预先选择的种子序列（可以是传输的随机种子列表或类似于number_nodes * node_number * iteration的公式）。重要的是，初始的“短”种子永远不会在节点之间重新使用。

每个节点将使用一个以此种子初始化的MT PRNG，使用n次，其中n <<< MT_period / max_value_of_short_seed（TT800是2⁸⁰⁰-1和MT19937是2¹⁹⁹³⁷-1，因此n仍然可以是一个非常大的数字）。在n次之后，节点移动到所选列表中的下一个种子。

虽然我不能（也不可能）提供“保证”没有节点会同时（或根本）具有重复序列，但这里是AMD关于使用不同种子的说法：（显然，初始播种算法起着重要作用）。

在这里介绍的四种创建多个流的方法中，这是最不令人满意的...例如，如果初始值不足够远，则从不同起点生成的序列可能会重叠。如果所使用的生成器的周期很大，则减少了重叠序列的可能性。尽管无法保证序列的独立性，但由于其极长的周期，使用带有随机起始值的Mersenne Twister不太可能引起问题，特别是如果所需的序列数量很少...

祝编程愉快。

- user166390

---

AMD提出的四种解决方案仅适用于其库案例。但是，当您处理原始MT实现（例如@Charles Brunet所做的）时，通过动态创建器进行参数化是处理独立Mersenne Twister流的最合适方式。 - jopasserat

Answer 3

我可以在每个节点上跳到已知的距离。但是是否有这样一种算法适用于Mersenne-Twister或其他良好的伪随机数生成器，可以在合理的时间和内存范围内完成？

是的，请参见http://theo.phys.sci.hiroshima-u.ac.jp/~ishikawa/PRNG/mt_stream_en.html。这是获得独立随机数流的绝佳解决方案。通过进行比每个流中所需的随机数数量更大的跳跃来创建每个流的起始点，可以避免流重叠。

- M. S. B.

Answer 4

免责声明：我不确定在从任意“uint”（或x，其中x是一个较小的任意但唯一值）种子开始时，MT在循环重叠方面有什么保证，但这可能值得研究，因为如果有保证，则只需将每个节点以不同的“uint”种子启动，那么帖子的其余部分基本上变得无关紧要。（MT的周期长度/阶数是令人震惊的，将UINT_MAX除去后仍然留下一个难以理解的数字，除了在纸上。）

---

好的，下面是我的评论回答...

我喜欢第二种方法，使用预生成的状态集；然后在每个节点中初始化MT以给定的起始状态。

当然，只有初始状态必须被保留，一旦生成了这些状态，这些状态可以：

1. 如果满足要求，可以无限重复使用；或者
2. 在模拟运行时，可以在外部快速框架上向前生成下一个状态；或者
3. 如果可靠的消息传递和序列在节点之间以相同的速率使用，并且满足要求，则节点可以报告回结束状态。

考虑到MT的生成速度很快，我不建议采用上述第三种方法，因为它太复杂了，而且有许多附加条件。选项#1很简单，但可能不够动态。

选项#2似乎是一个非常好的可能性。服务器（一个“快速机器”，不一定是节点）只需要传输下一个“未使用的序列块”的起始状态（比如说，十亿个周期）--节点将在请求新块之前使用发生器进行十亿个周期。这将使其成为帖子中#1的混合体，具有非常少的消息传递。

在我的系统上，一台Core2 Duo电脑上，我可以使用下面提供的代码（在LINQPad中运行）在17秒内生成十亿个随机数。我不确定这是哪种MT变体。

void Main()
{
    var mt = new MersenneTwister();
    var start = DateTime.UtcNow;
    var ct = 1000000000;
    int n = 0;
    for (var i = 0; i < ct; i++) {      
        n = mt.genrand_int32();
    }
    var end = DateTime.UtcNow;
    (end - start).TotalSeconds.Dump();
}

// From ... and modified (stripped) to work in LINQPad.
// http://mathnet-numerics.googlecode.com/svn-history/r190/trunk/src/Numerics/Random/MersenneTwister.cs
// See link for license and copyright information.
public class MersenneTwister
{
    private const uint _lower_mask = 0x7fffffff;
    private const int _m = 397;
    private const uint _matrix_a = 0x9908b0df;
    private const int _n = 624;
    private const double _reciprocal = 1.0/4294967295.0;
    private const uint _upper_mask = 0x80000000;
    private static readonly uint[] _mag01 = {0x0U, _matrix_a};
    private readonly uint[] _mt = new uint[624];
    private int mti = _n + 1;

    public MersenneTwister() : this((int) DateTime.Now.Ticks)
    {
    }       
    public MersenneTwister(int seed)
    {
                init_genrand((uint)seed);
    }       

    private void init_genrand(uint s)
    {
        _mt[0] = s & 0xffffffff;
        for (mti = 1; mti < _n; mti++)
        {
            _mt[mti] = (1812433253*(_mt[mti - 1] ^ (_mt[mti - 1] >> 30)) + (uint) mti);
            _mt[mti] &= 0xffffffff;
        }
    }

    public uint genrand_int32()
    {
        uint y;
        if (mti >= _n)
        {
            int kk;

            if (mti == _n + 1) /* if init_genrand() has not been called, */
                init_genrand(5489); /* a default initial seed is used */

            for (kk = 0; kk < _n - _m; kk++)
            {
                y = (_mt[kk] & _upper_mask) | (_mt[kk + 1] & _lower_mask);
                _mt[kk] = _mt[kk + _m] ^ (y >> 1) ^ _mag01[y & 0x1];
            }
            for (; kk < _n - 1; kk++)
            {
                y = (_mt[kk] & _upper_mask) | (_mt[kk + 1] & _lower_mask);
                _mt[kk] = _mt[kk + (_m - _n)] ^ (y >> 1) ^ _mag01[y & 0x1];
            }
            y = (_mt[_n - 1] & _upper_mask) | (_mt[0] & _lower_mask);
            _mt[_n - 1] = _mt[_m - 1] ^ (y >> 1) ^ _mag01[y & 0x1];

            mti = 0;
        }

        y = _mt[mti++];

        /* Tempering */
        y ^= (y >> 11);
        y ^= (y << 7) & 0x9d2c5680;
        y ^= (y << 15) & 0xefc60000;
        y ^= (y >> 18);

        return y;
    }
}

编程愉快。

- user166390

---

当您使用uint种子生成MT时，它只是使用线性同余生成器生成其内部状态。因此，您是否知道使用两个不同的种子生成的两个状态之间的接近程度。在进行蒙特卡罗模拟时，我正在寻找非常罕见的事件。如果我需要1e-12的概率，则需要至少1e14个随机数。您的解决方案可能有效，但绝对不是最佳方法。 - Charles Brunet

Answer 5

TRNG是一个专门为并行集群环境设计的随机数生成器（特别是为德国TINA超级计算机构建的）。因此，它非常容易创建独立的随机数流，并生成非标准分布。这里有一个关于如何设置它的教程： http://www.lindonslog.com/programming/parallel-random-number-generation-trng/

- Lindon