Random Number Generators

Overview

The ability to generate random numbers quickly is a prerequisite for using randomization techniques. This article introduces common methods for generating random numbers.

Random numbers and pseudorandom numbers

Calling a single number a “random number” is meaningless, so below we use “random sequence” by default. Even when “random number” is mentioned, it refers to an element of a random sequence.

The computation performed by a modern computer is deterministic. Therefore, it is impossible to generate a truly unpredictable and non-repeatable random sequence using an algorithm alone.

However, in most situations we do not need such strong randomness. We only need the generated sequence to have the statistical characteristics of random sequences, such as a uniform distribution and mutual independence. Such a sequence is called a pseudorandom number sequence.

Examples of applications of random and pseudorandom numbers in real life and algorithms include:

• Pseudorandom numbers are often sufficient for sample surveys, because we care only about statistical characteristics.
• Network security often requires stronger randomness than the pseudorandom numbers just mentioned, because an attacker may exploit predictability.
• Randomized algorithms used in OI/ICPC generally need only pseudorandom numbers. This is because such algorithms often introduce randomness into the complexity analysis, thereby reducing the complexity. In essence, they still use only the statistical characteristics of random numbers.
• Some randomized algorithms, such as the Moser algorithm, use entropy-related properties of random numbers and therefore require true random numbers.

Implementation

rand

It is used to generate pseudorandom numbers. Its drawback is that it is relatively slow, and it requires #include<cstdlib>.

Calling the rand() function returns a random nonnegative integer in [0,RAND_MAX], where RAND_MAX is a macro in the standard library. On Linux systems, RAND_MAX equals $2^{31}-1$. The modulus operation can be used to limit the size of the generated number.

Using rand() requires a random seed. The srand(seed) function can be used to change the random seed to seed, although leaving it uninitialized is also possible.

If the same program is run twice with the same seed, on the same machine and with the same compiler, the generated results will be identical.

One option is to use the current system time as the random seed: srand(time(nullptr)).

Regarding the randomness of rand() and rand()%n:

• The C/C++ standard makes no requirements regarding any aspect of the quality of the random numbers generated by rand().
• The implementation of rand() used by GCC guarantees basic properties such as uniformity of the distribution, but has obvious defects such as a short period in the low bits. For example, on the author's machine, the sequence generated by rand()%2 has a period of about $2\cdot 10^6$.
• Even assuming that rand() is uniformly random, rand()%n cannot guarantee uniformity, because each number in [0,n) may occur a different number of times among 0%n,1%n,...,RAND_MAX%n.

Predefined random number generators

Several especially popular algorithms are predefined. Unless otherwise specified, they are defined in the <random> header.

mt19937

It is a random number generator class whose purpose is the same as rand(), and whose random-number range is the same as that of the unsigned int type.

Its advantages are high random-number quality (one indication is a longer period before repetition; in other respects it is at least no worse than rand()), and a much higher speed than rand(). It requires #include<random>.

mt19937 is based on the 32-bit Mersenne Twister, designed by Matsumoto and Nishimura in 1998². To use it, simply define a random-number generator with it: std::mt19937 myrand(seed). seed may be omitted; when seed is omitted, the default random seed is used.

mt19937 overloads operator (); call myrand() whenever a random number is needed, and it returns one.

Another similar generator is mt19937_64, based on the 64-bit Mersenne Twister and designed by Matsumoto and Nishimura in 2000. It is used in the same way as mt19937, but its random-number range is expanded to the range of the unsigned long long type.

#include <ctime>
#include <iostream>
#include <random>

using namespace std;

int main() {
  mt19937 myrand(time(nullptr));
  cout << myrand() << endl;
  return 0;
}

minstd_rand0

The linear congruential algorithm was discovered by Lewis, Goodman, and Miller in 1969, and adopted by Park and Miller in 1988 as the “minimum standard”.

The formula is as follows, where $A,C,M$ are predefined constants.

minstd_rand() is the newer “minimum standard”, recommended by Park, Miller, and Stockmeyer in 1993.

For minstd_rand0(), $s$ is a 32-bit unsigned integer, $A$ is 16807, $C$ is 0, and $M$ is 2147483647.

For minstd_rand(), $s$ is a 32-bit unsigned integer, $A$ is 48271, $C$ is 0, and $M$ is 2147483647.

random_shuffle

It is used to randomly shuffle a specified sequence. It requires #include<algorithm>.

To use it, pass the first and last pointers or iterators of the specified range (left-closed and right-open): std::random_shuffle(first, last) or std::random_shuffle(first, last, myrand).

The random-number generator used internally defaults to rand(). A custom random-number generator may also be passed in.

Regarding the randomness of random_shuffle:

• The C++ standard requires random_shuffle to choose among all possible permutations with equal probability, but the GCC³ compiler does not implement this strictly.
• One reason for the randomness defect in GCC's random_shuffle is that it uses an expression such as rand()%n. As mentioned earlier, this does not generate uniformly random integers.
• The second reason is that the value range of rand() is limited. If the length of the supplied range exceeds RAND_MAX, some permutations cannot be generated⁴.

shuffle

It has the same purpose as random_shuffle. It requires #include<algorithm>.

The difference is that it must use a custom random-number generator: std::shuffle(first, last, myrand).

The GCC³ implementation of shuffle meets the requirements of the C++ standard: it chooses among all possible permutations with equal probability.

The following is a data generator written with rand() and random_shuffle(). It generates random small data for “ZJOI2012” Disaster.

#include <algorithm>
#include <cstdlib>
#include <ctime>
#include <iostream>

int a[100];

int main() {
  srand(time(nullptr));
  int n = rand() % 99 + 1;
  for (int i = 1; i <= n; i++) a[i] = i;
  std::cout << n << '\n';
  for (int i = 1; i <= n; i++) {
    std::random_shuffle(a + 1, a + i);
    int cnt = rand() % i;
    for (int j = 1; j <= cnt; j++) std::cout << a[j] << ' ';
    std::cout << 0 << '\n';
  }
}

Here is the same data generator written with mt19937 and shuffle().

#include <algorithm>
#include <ctime>
#include <iostream>
#include <random>

int a[100];

int main() {
  std::mt19937 rng(time(nullptr));
  int n = rng() % 99 + 1;
  for (int i = 1; i <= n; i++) a[i] = i;
  std::cout << n << '\n';
  for (int i = 1; i <= n; i++) {
    std::shuffle(a + 1, a + i, rng);
    int cnt = rng() % i;
    for (int j = 1; j <= cnt; j++) std::cout << a[j] << ' ';
    std::cout << 0 << '\n';
  }
}

Here is an implementation that randomly permutes the first ten positive integers.

#include <algorithm>
#include <iostream>
#include <iterator>
#include <random>

int main() {
  std::vector<int> v = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10};

  std::random_device rd;
  std::mt19937 g(rd());

  std::shuffle(v.begin(), v.end(), g);

  std::copy(v.begin(), v.end(), std::ostream_iterator<int>(std::cout, " "));
  std::cout << "\n";
}

Nondeterministic uniformly distributed integer random-number generator

random_device is a hardware-based uniformly distributed random-number generator that can generate random numbers quickly before the entropy pool is exhausted. This class was defined in C++11 and requires the random header. Since performance drops sharply after the entropy pool is exhausted, it is recommended to use this method to generate seeds for pseudorandom numbers such as those from mt19937, rather than to generate numbers directly.

random_device is a nondeterministic uniform random-bit generator. If nondeterministic random-number generation is unsupported, an implementation is allowed to implement it using a pseudorandom engine. So far, the author has not received any report that NOIP judging machines do not support hardware-based uniformly distributed random-number generation. Out of caution, however, it is recommended to use this algorithm to generate random seeds.

The reference code is as follows.

#include <iostream>
#include <map>
#include <random>
#include <string>

int main() {
  std::random_device rd;
  std::map<int, int> hist;
  std::uniform_int_distribution<int> dist(0, 9);
  for (int n = 0; n < 20000; ++n) {
    ++hist[dist(rd)];  // Note: for demonstration only; once the entropy pool is exhausted,
                       // the performance of many random_device implementations drops sharply.
                       // In practice, random_device is usually used only to
                       // seed a pseudorandom number generator such as mt19937.
  }
  for (auto p : hist) {
    std::cout << p.first << " : " << std::string(p.second / 100, '*') << '\n';
  }
}

A possible output is shown below.

0 : ********************
1 : *******************
2 : ********************
3 : ********************
4 : ********************
5 : *******************
6 : ********************
7 : ********************
8 : *******************
9 : ********************

Random-number distributions

This section concerns distributions in which the generated random numbers appear with specified probabilities, such as equal probability, the Bernoulli distribution, the binomial distribution, the geometric distribution, and the standard normal (Gaussian) distribution.

For the exact class names, see the list under Pseudorandom number generation — Random number distributions.

Implementation

The following program simulates a six-sided die.

#include <iostream>
#include <random>

int main() {
  std::random_device rd;   // Used to obtain a seed for the random-number engine
  std::mt19937 gen(rd());  // Seed the standard mersenne_twister_engine
  std::uniform_int_distribution<> dis(1, 6);

  for (int n = 0; n < 10; ++n)
    // Transform the random unsigned int generated by gen into an int in [1, 6]
    std::cout << dis(gen) << ' ';
  std::cout << '\n';
}

Other implementation methods

Sometimes we need to implement our own random-number generator. Here are some commonly used methods for generating random numbers.

Linear congruential random-number generator

Use the following formula to generate the random-number sequence $\{R_i\}$:

where $A,B,P$ are all constants.

This method is easy to implement, but the generated random sequence has a short period (the maximum period is $P$, but in most cases it is shorter than $P$).

#include <iostream>
using namespace std;

struct myrand {
  int A, B, P, x;

  myrand(int A, int B, int P) {
    this->A = A;
    this->B = B;
    this->P = P;
  }

  // Generate the next random number in the random sequence
  int next() { return x = (A * x + B) % P; }
};

myrand rnd(3, 5, 97);  // Initialize a random-number generator

int main() {
  int x = rnd.next();
  cout << x << endl;
  return 0;
}

Lagged Fibonacci random-number generator

Use the following formula to generate the random-number sequence $\{R_i\}$ (where $0 < j < k$):

Here $P$ is usually a power of $2$ (commonly $2^{32}$ or $2^{64}$), and $\star$ denotes a binary operator. Addition, subtraction, multiplication, or XOR may be used.

Compared with the traditional linear congruential random-number generator, this method has a longer period, but its randomness is strongly affected by the initial conditions.

#include <iostream>
#include <vector>
using namespace std;

struct myrand {
  vector<unsigned> vec;
  int l, j, k, cur;

  myrand(int l, int j, int k) {
    this->l = l;
    this->j = j;
    this->k = k;
    cur = 0;
    for (int i = 0; i < l; i++) {
      vec.push_back(rand());  // First generate the first few elements of the random sequence using another method
    }
  }

  unsigned next() {
    vec[cur] = vec[(cur - j + l) % l] * vec[(cur - k + l) % l];
    // unsigned makes reduction modulo 2^32 automatic here
    return vec[cur++];
  }
};

myrand rnd(11, 4, 7);

int main() {
  unsigned x = rnd.next();
  cout << x << endl;
  return 0;
}

References and notes