Assignment 'Randomness of individual bits'

For our upcoming application, Monte Carlo integration, the uniformness of the generated random numbers is most important. There are applications, however, where individual bits of generated numbers are used or play an important role. In this assignment you will verify the serious non-randomness of the old UNIX's rand generator, and see a inexpensive trick how this problem was 'circumvented'.

The old UNIX's rand generator

The original BSD Unix's rand generator is a linear congruential generator with a = 1103515245, c = 12345 and m = 2³¹ = 2147483648.
The period of this generator is 2³¹. Still, this rand() generator is not very good, since the low bits of successive numbers are correlated.

Later implementations of rand() solve the correlation, we'll get back to that later. For now, you will use the GNU Scientific Library (GSL), since this library also provides the original, problematic version of rand.

Use GSL's rand implementation

Extend your EX_RNGs.cc and header file, such that there is an additional GSL_RAND generator available. In the GSL documentation on random number generation, find out how to initialise and use their generators.

Compile your new program and inspect some results in the generated randdata.txt. Do not forget to produce nonnormalised numbers.

Comparison with a 'fixed' rand

Just as you just did, add a generator SUN_RAND that employs the fixed version of rand, as shipped by Sun. Here is the code that should be added to LC_RNGs.cc (don't forget to adapt the header file):

/*
 * SUN_RAND class
 */
generators::SUN_RAND::SUN_RAND() : LC_RNG(1103515245, 12345, 1<<31) {}
inline unsigned int generators::SUN_RAND::getMaximum() { return 0x7fff; }

/* Sun's implementation of nextInt, to solve correlation of low bits. */
unsigned int generators::SUN_RAND::nextInt() {
    _seed = ( _a*_seed + _c ); // automatic mod 2^32
    return (((_seed) >> 16) & 0x7fff);
}

Important: make sure that your two new generators GSL_RAND and SUN_RAND are seeded with the same initial value.
Run your program and check whether the generated numbers differ from what you just saw.

A bit of a problem

A tiny bit of theory: the bitwise and (&) operator performs a logical and on each pair of bits of two numbers. What is the result of 0b1100 & 0b1010? Think of a way to 'extract' the bit value at position i from an integer variable using as little and as cheap operator(s) as possible.

Extend your main function, such that it prints each generated number on one line in your terminal, followed by the (last few) bits of its binary representation. Do this by writing your own code; do not use external bits of code such as the standard bitset class. Try to be as efficient as possible. An example:

0       : 16838  = ...000110
1       : 5758   = ...111110
2       : 10113  = ...000001

A bit pattern and a fix

Inspect the last few bits of the generated numbers, for both GSL_RAND and SUN_RAND. Which of them (or both?) has a problem, can you find a pattern in the generated samples?

Try to understand what the Sun implementation (see the SUN_RAND class above) does to avoid the problem you just found.

Hand in

Send your files rand.cc, LC_RNGs.cc, EX_RNGs.cc and relevant header (.h) files by e-mail. In the email, shortly note what you found to be the problem with the UNIX (and maybe also the SUN) RNG. Also answer these final questions in your e-mail: The 'fix' by Sun produces numbers that are bounded by some value. What is this value, and why? Can you give an upper bound for the period of this generator? Is the SUN RNG a linear congruential RNG?