boost 随机数生成器的丢弃行为

Question

我正在使用适配器 class 包装增强随机数生成器以实现 Monte Carlo 例程。在class的成员函数上写单元测试时，我假设.discard(unsigned int N)的行为是抽取N个随机数而不存储它们，从而推进rng的状态。提升代码是：

void discard(boost::uintmax_t z)
{
    if(z > BOOST_RANDOM_MERSENNE_TWISTER_DISCARD_THRESHOLD) {
        discard_many(z);
    } else {
        for(boost::uintmax_t j = 0; j < z; ++j) {
            (*this)();
        }
    }
}

这支持了我的假设。但是，我发现 .discard(1) 产生的序列与没有丢弃的相同序列没有一个数字不同。代码：

#include <iostream>
#include <iomanip>
#include <random>
#include <boost/random.hpp>

int main()
{
    boost::mt19937 uGenOne(1);
    boost::variate_generator<boost::mt19937&, boost::normal_distribution<> > distOne(uGenOne, boost::normal_distribution<>());

    boost::mt19937 uGenTwo(1);
    boost::variate_generator<boost::mt19937&, boost::normal_distribution<> > distTwo(uGenTwo, boost::normal_distribution<>());
    distTwo.engine().discard(1);

    unsigned int M = 10;
    std::vector<double> variatesOne(M);
    std::vector<double> variatesTwo(M);

    for (unsigned int m = 0; m < M; ++m) {
        variatesOne[m] = distOne();
        variatesTwo[m] = distTwo();
    }

    for (unsigned int m = 0; m < M; ++m)
        std::cout << std::left << std::setw(15) << variatesOne[m] << variatesTwo[m] << std::endl;

    return 0;
}

产出

2.28493        0.538758  
-0.668627      -0.0017866
0.00680682     0.619191  
0.26211        0.26211   
-0.806832      -0.806832 
0.751338       0.751338  
1.50612        1.50612   
-0.0631903     -0.0631903
0.785654       0.785654  
-0.923125      -0.923125

我对 .discard 如何运作的解释不正确吗？为什么两个序列在前三个输出中不同，然后相同？

（此代码是在 msvc 19.00.23918 和 g++ 4.9.2 上在 cygwin 上编译的，结果相同）。

Answer 1

这里的问题似乎是引擎没有被正确修改或者发行版正在添加一些额外的工作。如果我们像

那样直接使用引擎

int main()
{
    boost::mt19937 uGenOne(1);

    boost::mt19937 uGenTwo(1);
    uGenTwo.discard(1);

    unsigned int M = 10;
    std::vector<double> variatesOne(M);
    std::vector<double> variatesTwo(M);

    for (unsigned int m = 0; m < M; ++m) {
        variatesOne[m] = uGenOne();
        variatesTwo[m] = uGenTwo();
    }

    for (unsigned int m = 0; m < M; ++m)
        std::cout << std::left << std::setw(15) << variatesOne[m] << variatesTwo[m] << std::endl;

    return 0;
}

它产生

1.7911e+09     4.28288e+09
4.28288e+09    3.09377e+09
3.09377e+09    4.0053e+09
4.0053e+09     491263
491263         5.5029e+08
5.5029e+08     1.29851e+09
1.29851e+09    4.29085e+09
4.29085e+09    6.30312e+08
6.30312e+08    1.01399e+09
1.01399e+09    3.96591e+08

因为我们丢弃了第一个输出，所以这是一个 1 移位序列。

所以你对 discard 的工作原理是正确的。我不确定为什么在通过 boost::variate_generator 执行此操作时会出现差异。我不明白为什么前三个数字不同，但所有其余输出都匹配。

Answer 2

只是添加上一个答案评论中的一个重要细节。正如@NathanOliver 提到的，.discard 增加生成器，生成器将制服发送到正态分布，将制服转换为正态分布。 boost::normal_distribution 使用 Ziggurat algorithm 这是一种 "acceptance/rejection" 算法。它随机绘制一个制服，对其进行操作，然后检查它是否处于所需的分布中。如果不是，则拒绝并随机生成新的制服。

for(;;) {
        std::pair<RealType, int> vals = generate_int_float_pair<RealType, 8>(eng);
        int i = vals.second;
        int sign = (i & 1) * 2 - 1;
        i = i >> 1;
        RealType x = vals.first * RealType(table_x[i]);
        if(x < table_x[i + 1]) return x * sign;
        if(i == 0) return generate_tail(eng) * sign;
        RealType y = RealType(table_y[i]) + uniform_01<RealType>()(eng) * RealType(table_y[i + 1] - table_y[i]);
        if (y < f(x)) return x * sign;
    }

关键的一点是，如果最后一个if失败，那么for循环会再次启动，并且会再次触发对generate_int_float_pair的调用。这意味着底层生成器递增的次数是未知的。

因此，正常序列将具有不同的编号，直到子序列中拒绝总和相同的点，此时剩余的均匀序列相同。这发生在问题中发布的示例的第三个位置。（它实际上有点微妙，因为底层生成器可以在 Ziggurat 算法中被调用一次或两次，但本质是相同的——一旦序列同步，它们就永远不会产生不同的变量）。

boost 随机数生成器的丢弃行为

Behavior of discard for boost random number generators

c++

boost-random