# QuantLib: Using Boost Random Number Generators

In the previous post on random number generators I introduced all the random number generators (RNGs) that are available in QuantLib. In most circumstances these generators will be sufficient. However, if you take a look at the boost/random package, you will find a wealth of random number generators with different properties. In some cases you might want to test if your results depend on the specific implementation of the random number generator and whether some other generator might give better, or more dependable results. In this post I will write two wrapper classes that will allow any Boost random number generator to be used with QuantLib.

### Anatomy of Boost Integer RNGs

To understand the differences between the Boost generators and the generators of QuantLib, I will first analyse the anatomy of the Boost random number generators. Boost knows two types of random number generators, integral and floating point generators. As the name suggests, integral generators return integral values and floating point generators return floating point values. Let’s look at the integral generators first. The integral type that a generator returns can have different bit length. For example the mt19937_64 is a 64 bit version of the Mersenne twister with periodicity of $$2^{19937}-1$$. Every generator class in the Boost library defines the result_type as the type that the generator returns. For example rand48 contains the line

typedef boost::uint32_t result_type;

All Boost random number generators override the function operator to access the next random number.

result_type operator()();

This will generate and return the next number in the sequence. Random number generators all have a default constructor and a constructor that takes a seed. For all integral engines the seed type is the same as the result type. In addition there is a seed method that re-initialises a random number generator with a seed. For example rand48 has the following constructors.

rand48();
rand48(result_type seed);

And the seed method is defined as follows.

seed(result_type seed);

There is also a seed function without arguments but we will only discuss the methods that are needed to create our generic wrapper class.

### A Generic QuantLib Wrapper for Boost’s Integral RNGs

The aim is to write a generic wrapper around all of Boost’s random number generators. The idea is to follow the pattern used in the Ranlux3UniformRng wrapper for Boost’s ranlux64_3_01 RNG. The difference will be that this wrapper will be templated with the type of the Boost RNG. I will start by creating the GenericBoostIntRNG class with the template argument RNG

template<class RNG>
class GenericBoostIntRNG {

Here RNG is the tye of the Boost random number generator. GenericBoostIntRNG needs to keep an instance of RNG as private member. Here we encounter the first little hurdle. Getting the next random number from a Boost RNG will modify the object. The reasoning is that the internal state of the RNG is changed and this manifests itself in the fact that the function operator is non-const. QuantLib, on the other hand, expects that getting the next random number does not modify the generator. The next method is a const method in all the QuantLib RNGs. The solution is to make the Boost RNG a mutable member of the wrapper.

  private:
mutable RNG generator;

Next we have to declare two local types. QuantLib demands the sample_type to define the return type of the next method. On the other hand, we also need to use result_type of RNG in many places and should therefore declare a shortcut.

  public:
typedef Sample<Real> sample_type;
typedef typename RNG::result_type result_type;

To make sure that we use GenericBoostIntRNG only with integral random number generators we should make sure that result_type is an integer type. This can be done with BOOST_STATIC_ASSERT and the standard numeric_limits type.

    BOOST_STATIC_ASSERT (std::numeric_limits<result_type>::is_integer);

At this point we should be defining the constructor. This is not straightforward and will require some extra code because we might need to create an automatic seed. We will postpone this and first focus on the easier methods that generate the random numbers.

    Real nextReal() const {
return (Real(generator()) + 0.5) / ceil_value;
}

This method creates a random integer number and convert it into a Real value. The ceil_value will be one greater than the maximum number returned by the random number generator. We define it as a  static private member of the class.

  private:
static Real ceil_value;

The value is set by using the max function of the Boost random number generator.

template<class RNG>
Real GenericBoostIntRNG<RNG>::ceil_value = Real(RNG::max()) + 1.0;

With this definition, the call to generator in the nextReal method will always return a value from 0 to ceil_value-1. This means that the value returned by nextReal will always fall within the closed interval.

$$\left[ \frac{0.5}{ceil}, \frac{ceil-0.5}{ceil}\right]$$
Mathematically this interval will never include 0 or 1 but we have to be careful when stating that the numerical value will lie in the open interval $$(0,1)$$. Rounding errors might cause the upper value to round to 1. This will be the case when the Real type has less bits that the integral type of the Boost RNG. If we stick to using 32 bit Boost RNGs and define Real as 64 bit double than all is fine and the numbers will always be less than 1. In any case the number will never be larger than 1 and an occasional value of 1 should not cause a problem. QuantLib random number generators are expected to define a next method that returns a sample_type. This is quickly done.

  public:
sample_type next() const {
return sample_type(nextReal(), 1.0);
}

Now we come back to creating the constructor for our generic wrapper class.

  public:
GenericBoostIntRNG(seed_type seed = 0)
: generator(getSeed(seed))
{}

This does not look complicated. We are initializing the Boost object generator with a seed that we obtain from a method called getSeed. The getSeed method will have to be implemented. If we specify a seed that is non-zero that we will simply pass it on. Otherwise we will use the QuantLib SeedGenerator class to create a suitable seed for the generator object.

  private:
result_type getSeed(result_type seed) {
unsigned long s = (seed != 0 ? seed : SeedGenerator::instance().get());
return (result_type) (s % RNG::max());
}

In the first line inside the method we assign the seed to an unsigned long. The SeedGenerator class in QuantLib is defined as a singleton. A call to SeedGenerator::instance().get() will return an unsigned long seed value based on the current time. Because unsigned long might be longer than the result_type of the Boost RNG we have to restrict the values by taking the modulus with the maximum value of the RNG.

With this the generic wrapper is complete. I have tested it by plugging various Boost random number generators into the EquityOption example of QuantLib. When using the Boost’s mt19937 generator, I get identical results to the QuantLib MersenneTwisterUniformRng.

### Wrapping Boost’s Floating Point Generators

Boost also specifies a number of floating point random number generators. For completeness I have also created a wrapper for using these within QuantLib. Using the experience gained from GenericBoostIntRNG the implementation of GenericBoostFloatRNG is relatively straightforward.

template<class RNG>
class GenericBoostFloatRNG {
private:
mutable RNG generator;
static Real ceil_value;
public:
typedef Sample<Real> sample_type;

typedef typename RNG::result_type result_type;
typedef boost::uint32_t seed_type;
BOOST_STATIC_ASSERT (std::numeric_limits<result_type>::is_float);

GenericBoostFloatRNG(seed_type seed = 0)
: generator(getSeed(seed))
{}

sample_type next() const {
return sample_type(nextReal(), 1.0);
}

result_type nextFloat() const {
return generator();
}

private:
seed_type getSeed(seed_type seed) {
unsigned long s = (seed != 0 ? seed : SeedGenerator::instance().get());
return (seed_type) (s % 0xffffffff);
}
};

There are three differences between this class and the GenericBoostIntRNG class. First, the nextFloat method can return the result of the generator directly. Second, we need to define the seed_type separately from the result_type. The seed_type will be integral while the result_type will be a floating point type. Third, the restriction of the seed to values accepted by the Boost RNG looks slightly different. All floating point random number generators accept uint32 as seed type. This means we can simply restrict the seed by applying a bit mask.

I hope someone will find this useful. The code for these wrapper classes is freely available under the MIT licence and available on Github. 