C++ gamma_distribution 返回无穷大

Question

我最近写了一些代码，需要从 Gamma 分布中抽取大量数据。我使用标准的 gamma_distribution 方法实现了这个，但我发现偶尔(难得一见)这会返回值 1。#INF 没有明显的原因。

下面是我发现的一个最简单的例子，它展示了这个问题。在我的测试中，我发现问题通常发生在第 10 亿次迭代左右。

#include "stdafx.h"
#include <iostream>
#include <random>

std::random_device rd;
std::default_random_engine generator(rd());

int _tmain(int argc, _TCHAR* argv[])
{
    // create gamma distribution and random variable x
    std::gamma_distribution<double> rgamma(2.0,1.0);
    double x;

    // loop through a large number of iterations
    for (unsigned long int i=0; i<int(4e9); i++) {

        // print update to console every million iterations
        if ((i+1)%int(1e6)==0)
            std::cout << "iteration: " << (i+1)/1e6 << " million\n";

        // draw new value of x from gamma distribution
        x = rgamma(generator);

        // if x==infinity then break
        if ((1.0/x)==0) {
            std::cout << "Error at iteration " << (i+1) << ": x=" << x << "\n";
            std::cin.get();
            exit(1);
        }
    }
    // print message if reach end of loop
    std::cout << "end\n";
    std::cin.get();
    return 0;
}

我很想知道其他人是否可以复制这个问题。我不知道这是否相关，但上面的程序是在 Visual Studio 2010 中作为 Win32 应用程序编写的，并在具有 8 核 Intel 处理器的 Windows 机器上运行。

目前我已经通过捕获无限值并将它们变成大数字来解决这个问题。但是，如果有人对发生这种情况的原因/方式有任何见解，我们将不胜感激!

Answer 1

我通过以下方式使随机种子具有确定性:

//std::random_device rd;
std::default_random_engine generator(-246744094);

并且可以在 8.32 亿次迭代中始终如一地重现此错误。我用不同的编译器(VS2013 32/64 位和 Intel C++ 2016 64 位)和不同的机器(我的桌面 i7 和集群节点 Xeon E5)得到了相同的结果。我的 rgamma.max() 也同意 Bob 的。

然后我用 MinGW (g++ main.cpp -o main.exe -std=c++0x -O3) 编译它 - 它运行时没有产生任何无限。

所以 - 我假设:这是通过 Visual Studio 链接的库中的某种错误，无论您使用的是 MS 还是 Intel 编译器。不确定是否可以报告/修复 - 或者是否可以更快地找到一些可靠的 gammas 源代码。