Seeding independent bits engine with boosts multiprecision

Question

I have a simple rejection sampling application which is wrapped in a class and used externally as the dummy example below shows. I was able to adapt this post to a boost::multiprecision use case. However I'm not sure how to appropriately seed the generator and can't find any random_device equivalent for boost.

The below code 'works', but if you run it multiple times in quick succession you will get the same random numbers which I don't want. Is there something more sensitive than time(NULL)?

#include <iostream>
#include <boost/multiprecision/cpp_dec_float.hpp>
#include <boost/multiprecision/cpp_int.hpp>
#include <boost/random.hpp>

using namespace boost::multiprecision; // used only for SO post
using namespace boost::random;

typedef independent_bits_engine<boost::mt19937, std::numeric_limits<cpp_dec_float_50>::digits, cpp_int> generator;


generator &gen()
{
  thread_local static generator genny(time(NULL));
  return genny;  
}

class dummy {
    public:
        dummy() = default;
        cpp_dec_float_50 rejectionSample() {
           uniform_real_distribution<cpp_dec_float_50> ur(0,1);
           cpp_dec_float_50 x = ur(gen());
           while (x > 0.1) 
               x = ur(gen());
           return x;
        }
};



int main()
{
    std::cout << std::setprecision(std::numeric_limits<cpp_dec_float_50>::digits10) << std::showpoint;

    dummy d;
    int nToGet = 5;
    for (int i = 0; i < nToGet; ++i) 
        std::cout << d.rejectionSample() << std::endl;
}

Have you thought about using chrono::steady_clock::now().time_since_epoch().count() as the seed? Eg. you could count milliseconds? Generally I believe chrono's clock is more precise than time(). — Marek Piotrowski, Aug 14 '20 at 12:12

score 1 · Accepted Answer · answered Aug 14 '20 at 13:22

This works:

#include <iostream>
#include <boost/multiprecision/cpp_dec_float.hpp>
#include <boost/multiprecision/cpp_int.hpp>
#include <random>
#include <boost/random.hpp>

using namespace boost::multiprecision;

typedef boost::random::independent_bits_engine<boost::mt19937, std::numeric_limits<cpp_dec_float_50>::digits, cpp_int> generator;

class dummy {
public:
    dummy()
    {
        ur_ = boost::random::uniform_real_distribution<cpp_dec_float_50>(0,1);
        std::random_device rd;
        gen_ = generator(rd());
        ur_(gen_);
    }

    cpp_dec_float_50 rejectionSample() {
        cpp_dec_float_50 x = ur_(gen_);
        while (x > 0.1)
        {
            x = ur_(gen_);
        }
        return x;
    }
private:
    boost::random::uniform_real_distribution<cpp_dec_float_50> ur_;
    generator gen_;
};



int main()
{
    std::cout << std::setprecision(std::numeric_limits<cpp_dec_float_50>::digits10) << std::showpoint;
    
    dummy d;
    int nToGet = 50;
    for (int i = 0; i < nToGet; ++i) {
        std::cout << d.rejectionSample() << std::endl;
    }
}

So your misunderstandings are:

You were calling a constructor to uniform_real_distribution on every call. Not wrong per se, but expensive.
You were re-seeding the generator on every call, with time(NULL). That's really bad; you're not getting really even getting a pseudorandom sequence that way; you're using the non-decreasing clock to seed a generator and using only the first value that comes out.
It seems like you were confused by the syntax of how the distributions work, since you had a call-like syntax to gen(). You are passing a callable to uniform_real_distribution.
You shouldn't use time(NULL) for seeding generators anymore; we have the rdseed assembly instruction which is much better. This is what is called by std::random_device.

Finally, doesn't uniform_real_distribution do rejection sampling internally? So why do you need to do it?

Thanks. This is basically what I wanted; a 'random' class that effectively wraps a given distribution type (`uniform_real` in this case) and can be used more flexibly. By rejection sample I just mean sampling some arbitrary 'space' via a rejection method. — algae, Aug 15 '20 at 12:33
@algae: I see; so the rejection sampling you were doing here was just a stand-in for some more complicated goal — user14717, Aug 15 '20 at 13:44
Yes. A minor question: Suppose `dummy` objects were instantiated and used in multiple places. Then, is it possible to 'globally' fix a specific seed (for debugging) so that each call to the program gives the same results? For instance you could overload the `dummy` constructor to accept a specific seed, but then every instance of `dummy` in your code would need to be changed. — algae, Aug 16 '20 at 02:32

Seeding independent bits engine with boosts multiprecision

1 Answers1