Navigation

Fast skipping in a linear congruential generator

Introduction

The LCG is perhaps the simplest pseudorandom number generator (PRNG) algorithm. The generator state consists of a single integer, and each step involves a multiplication, addition, and modulo operation. The behavior of an LCG is defined by the following recurrence relation:

The integer constants \(a\), \(b\), \(m\) are chosen carefully to give good properties. For example: having a full period of \(m\), having \(m\) be a power of 2 (which is efficient on computers), and/or having \(b\) be 0 (to save an operation). The seed value can be chosen (almost) arbitrarily for each random number generator instance. The sequence of values \(x_0, x_1, x_2, \ldots, x_k, \ldots\) represents the pseudorandom numbers produced by the generator. These values can be used as-is, or they can be filtered and have their distribution altered.

Fast skipping

Suppose we know \(x_0\) and we want to determine \(x_n\) for some \(n \in \mathbb{N}\). Clearly it is possible to iterate the recurrence formula \(n\) times to get our desired result, which takes \(\Theta(n)\) time. But with some clever arithmetic, we can get the result much faster, in only \(\Theta(\log n)\) time.

First, consider the special simple case where \(b = 0\). Then \(x_n \equiv a^n x_0 \text{ mod } m\). We can compute \(a^n x_0 \text{ mod } m\) quickly using the well-known modular exponentiation algorithm, which is exponentiation by squaring with a reduction of each intermediate result modulo \(m\).

Now for the general case, consider the LCG recurrence relation without the modulo and with \(x \in \mathbb{R}\):

Backward iteration

As long as \(a\) and \(m\) are coprime, the multiplicative inverse of \(a\) exists and can be computed by the extended Euclidean algorithm.

The formula above calculates one step backward. We can do fast backward skipping using the forward forward formula, but replacing its \(a\) variables to \((a^{-1} \text{ mod } m)\) and its \(b\) variables to \((-(a^{-1} b) \text{ mod } m)\).

Source code

The fast skipping and backward iteration concepts discussed in the text are implemented in the following runnable programs:

Notes

[0]: Finding this closed-form solution requires some ingenuity, a table of algebraic expansions, or a computer algebra system (CAS). But now that the solution has been stated, you can prove its correctness by induction quite easily.

[1]: Proof of \(\displaystyle\frac{a^n - 1}{a - 1} = \sum_{i=0}^{n-1} a^i \) is left as an exercise for the reader.

[2]: This operation is correct if we revisit the definition of the modulus operator: \(r = x \text{ mod } m\) finds the unique pair of integers \(q\) and \(r\) such that \(0 \le r < m\) and \(x = qm + r\).