分布和内部状态

Distributions and internal state

在 Whosebug 上,有很多关于从先验未知范围生成均匀分布整数的问题。例如

典型的解决方案是这样的:

inline std::mt19937 &engine()
{
  thread_local std::mt19937 eng;
  return eng;
}

int get_int_from_range(int from, int to)
{
  std::uniform_int_distribution<int> dist(from, to);
  return dist(engine());
}

鉴于分布应该是一个轻量级对象并且不存在多次重新创建它的性能问题,似乎即使是简单的分布也可能很好并且通常会有 some internal state.

所以我想知道是否通过不断重置它来干扰分布的工作方式(即在每次调用 get_int_from_range 时重新创建分布)我得到正确分布的结果。

Pete Becker 和 Steve Jessop 之间有一个 long discussion,但没有最后一句话。 在另一个问题中(Should I keep the random distribution object instance or can I always recreate it?)内部状态的"problem"似乎不​​是很重要。

C++标准对这个问题有什么保证吗?

下面的实现(来自 N4316 - std::rand replacement)是否更可靠一些?

int get_int_from_range(int from, int to)
{
  using distribution_type = std::uniform_int_distribution<int>;
  using param_type = typename distribution_type::param_type;

  thread_local std::uniform_int_distribution<int> dist;
  return dist(engine(), param_type(from, to));    
}

编辑

这重用了一个可能的分布内部状态,但它很复杂,我不确定它是否值得这么麻烦:

int get_int_from_range(int from, int to)
{
  using range_t = std::pair<int, int>;
  using map_t = std::map<range_t, std::uniform_int_distribution<int>>;

  thread_local map_t range_map;

  auto i = range_map.find(range_t(from, to));
  if (i == std::end(range_map))
    i = range_map.emplace(
          std::make_pair(from, to),
          std::uniform_int_distribution<int>{from, to}).first;

  return i->second(engine());
}

(来自)

有趣的问题。

So I was wondering if interfering with how the distribution works by constantly resetting it (i.e. recreating the distribution at every call of get_int_from_range) I get properly distributed results.

我已经用 uniform_int_distributionpoisson_distribution 编写了测试代码。如果您愿意,可以很容易地扩展它来测试另一个发行版。答案好像是.

样板代码:

#include <random>
#include <memory>
#include <chrono>
#include <utility>

typedef std::mt19937_64 engine_type;

inline size_t get_seed()
    { return std::chrono::system_clock::now().time_since_epoch().count(); }

engine_type& engine_singleton()
{  
    static std::unique_ptr<engine_type> ptr;

    if ( !ptr ) 
        ptr.reset( new engine_type(get_seed()) );
    return *ptr;
}

// ------------------------------------------------------------------------

#include <cmath>
#include <cstdio>
#include <vector>
#include <string>
#include <algorithm>

void plot_distribution( const std::vector<double>& D, size_t mass = 200 )
{
    const size_t n = D.size();
    for ( size_t i = 0; i < n; ++i ) 
    {
        printf("%02ld: %s\n", i, 
            std::string(static_cast<size_t>(D[i]*mass),'*').c_str() );
    }
}

double maximum_difference( const std::vector<double>& x, const std::vector<double>& y )
{
    const size_t n = x.size(); 

    double m = 0.0;
    for ( size_t i = 0; i < n; ++i )
        m = std::max( m, std::abs(x[i]-y[i]) );

    return m;
}

实际测试代码:

#include <iostream>
#include <vector>
#include <cstdio>
#include <random>
#include <string>
#include <cmath>

void compare_uniform_distributions( int lo, int hi )
{
    const size_t sample_size = 1e5;

    // Initialize histograms
    std::vector<double> H1( hi-lo+1, 0.0 ), H2( hi-lo+1, 0.0 );

    // Initialize distribution
    auto U = std::uniform_int_distribution<int>(lo,hi);

    // Count!
    for ( size_t i = 0; i < sample_size; ++i )
    {
        engine_type E(get_seed());

        H1[ U(engine_singleton())-lo ] += 1.0;
        H2[ U(E)-lo ] += 1.0;
    }

    // Normalize histograms to obtain "densities"
    for ( size_t i = 0; i < H1.size(); ++i )
    {
        H1[i] /= sample_size; 
        H2[i] /= sample_size; 
    }

    printf("Engine singleton:\n"); plot_distribution(H1);
    printf("Engine creation :\n"); plot_distribution(H2);
    printf("Maximum difference: %.3f\n", maximum_difference(H1,H2) );
    std::cout<< std::string(50,'-') << std::endl << std::endl;
}

void compare_poisson_distributions( double mean )
{
    const size_t sample_size = 1e5;
    const size_t nbins = static_cast<size_t>(std::ceil(2*mean));

    // Initialize histograms
    std::vector<double> H1( nbins, 0.0 ), H2( nbins, 0.0 );

    // Initialize distribution
    auto U = std::poisson_distribution<int>(mean);

    // Count!
    for ( size_t i = 0; i < sample_size; ++i )
    {
        engine_type E(get_seed());
        int u1 = U(engine_singleton());
        int u2 = U(E);

        if (u1 < nbins) H1[u1] += 1.0;
        if (u2 < nbins) H2[u2] += 1.0;
    }

    // Normalize histograms to obtain "densities"
    for ( size_t i = 0; i < H1.size(); ++i )
    {
        H1[i] /= sample_size; 
        H2[i] /= sample_size; 
    }

    printf("Engine singleton:\n"); plot_distribution(H1);
    printf("Engine creation :\n"); plot_distribution(H2);
    printf("Maximum difference: %.3f\n", maximum_difference(H1,H2) );
    std::cout<< std::string(50,'-') << std::endl << std::endl;

}

// ------------------------------------------------------------------------

int main()
{
    compare_uniform_distributions( 0, 25 );
    compare_poisson_distributions( 12 );
}

运行它here.


Does the C++ standard make any guarantee regarding this topic?

据我所知没有。但是,我要说的是,该标准隐含地建议不要每次都重新创建引擎;对于任何分布 DistribDistrib::operator() 的原型采用引用 URNG& 而不是 const 引用。这是可以理解的要求,因为引擎可能需要更新其内部状态,但这也意味着代码看起来像这样

auto U = std::uniform_int_distribution(0,10);
for ( <something here> ) U(engine_type());

不编译,这对我来说是一个明确的动机,不要写这样的代码。


我相信有很多关于如何正确使用随机库的建议。如果您必须处理使用 random_devices 的可能性并允许用于测试目的的确定性播种,它确实会变得复杂,但我认为也可以提出我自己的建议:

#include <random>
#include <chrono>
#include <utility>
#include <functional>

inline size_t get_seed()
    { return std::chrono::system_clock::now().time_since_epoch().count(); }

template <class Distrib>
using generator_type = std::function< typename Distrib::result_type () >;

template <class Distrib, class Engine = std::mt19937_64, class... Args>
inline generator_type<Distrib> get_generator( Args&&... args )
{ 
    return std::bind( Distrib( std::forward<Args>(args)... ), Engine(get_seed()) ); 
}

// ------------------------------------------------------------------------

#include <iostream>

int main()
{
    auto U = get_generator<std::uniform_int_distribution<int>>(0,10);
    std::cout<< U() << std::endl;
}

运行它here。希望这对您有所帮助!

编辑 我的第一个推荐是错误的,对此我深表歉意;我们不能像上面的测试那样使用单例引擎,因为这意味着两个统一的 int 分布会产生相同的随机序列。相反,我依赖于 std::bind 使用自己的种子在 std::function 中本地复制新创建的引擎这一事实,这会产生预期的行为;具有相同分布的不同生成器会产生不同的随机序列。