C++11 类 <random> 的分布如何转换底层生成器？

Question

以下代码的行为似乎并不直观：

#include <random>
#include <iostream>
using namespace std;

int main()
{

  mt19937 MyGenerator(40);
  auto gauss = normal_distribution<double>(0,1);
  auto linear = uniform_real_distribution<double>(0,1);
  cout << gauss(MyGenerator) << endl; //line a
  cout << linear(MyGenerator) << endl; //line b
  cout << gauss(MyGenerator) << endl;
}

运行此代码给出输出

-0.816097
 0.705030
 0.303032.

如果现在交换 a 和 b 行的顺序，则输出变为

 0.644008
 0.338080
-0.639501.

前两个数字现在完全不同了，因为它们是由不同的分布产生的。然而，为什么第三个数字不同？ 在我的直觉中，分布应该抓取一个数字 c = MyGenerator() ，然后将其映射到特定范围内的随机数。在分发调用之后，随机数生成器将指向数字序列中的下一个数字。那么，在这两种情况下，第三次调用的结果不应该相同吗？

另一个观察： 实际上，向任一发行版添加第四次调用似乎确实重现了相同的数字。

Answer 1

libstdc++ 的normal_distribution 实现使用Marsaglia polar method。这种方法的有趣之处在于，每次传递都使用来自 URNG 的 两个 随机数来生成两个结果。

也就是说，对分布的第一次调用调用了 URNG 两次（可能更多次，因为它使用拒绝采样，但次数是偶数次）并返回一个结果；以下对分发的调用不会调用 URNG，但会返回保存的第二个结果。

这是extract from the source code，稍微重新格式化：

if (_M_saved_available)
{
    _M_saved_available = false;
    ret = _M_saved;
}
else
{
    result_type x, y, r2;
    do
    {
        x = result_type(2.0) * aurng() - 1.0;
        y = result_type(2.0) * aurng() - 1.0;
        r2 = x * x + y * y;
    }
    while (r2 > 1.0 || r2 == 0.0);

    const result_type mult = std::sqrt(-2 * std::log(r2) / r2);
    _M_saved = x * mult;
    _M_saved_available = true;
    ret = y * mult;
}

Answer 2

分发不需要为每个值调用底层生成器一次。一些分布最好通过组合多个随机值来计算。

例如，在 GNU 实现中，均匀分布的实现是

return (__aurng() * (__p.b() - __p.a())) + __p.a();

调用生成器__aurng一次；而正态分布的核心是：

do
{
    __x = result_type(2.0) * __aurng() - 1.0;
    __y = result_type(2.0) * __aurng() - 1.0;
    __r2 = __x * __x + __y * __y;
}
while (__r2 > 1.0 || __r2 == 0.0);

至少调用两次。