Many problems in engineering, finance and statistics can not be solved by direct methods, but a great number of them can be solved approximately using randomized algorithms. All those algorithms need flexible and efficient pseudorandom number generators. An effective implementation of PRNG in the F# language is somewhat tricky.
The F# language is based on .Net framework, which already provides a fast pseudorandom generator class
System.Random
, but it meets only certain statistical requirements for randomness. However, there are
two classes of generators that have good performance:
 Combined multiple recursive generators (e.g., L’Ecuyer’s
MRG32k3a
)  Twisted general feedback shift register generators (e.g., Matsumoto and Nishimura’s
MT19937ar
)
Despite the fact that MT19937ar
has longer period than MRG32k3a
, the state space in
MT19937ar
is very large and it is somewhat slow. Thus, MRG32k3a
becomes a better
alternative for many use cases.
Algorithm
MRG32k3a
is just one of several good parameters sets that L’Ecuyer obtained for CMRGs ^{1}.
It has two components of order 3. The state at step \(n\) is a pair of vectors
$S_{1,n}=(s_{1,n},s_{1,n+1}, s_{1,n+2})$ and $S_{2,n}=(s_{2,n},s_{2,n+1}, s_{2,n+2})$, which
change according to the linear recurrences:
The uniformly distributed output \(u_n\) is obtained by
\[z_n = (s_{1,n}s_{2,n}) \mathop{mod} 4294967087\] \[u_n=\begin{cases} z_n / 4294967088, & \text{if } z > 0 \\ 4294967087 / 4294967088, & \text{if } z = 0 \end{cases}\]Implementation
Two key requirement for the random number generator are the computational performance and
minimal memory consumption. MRG32k3a
algorithm requires only a pair of 3 element
double valued arrays for storing the current generator state, but F# arrays are heap allocated.
Fortunately, the F# language supports structure types, which can be allocated on the stack:
[<Struct>]
type State32k3a =
val mutable S10 : double
val mutable S11 : double
val mutable S12 : double
val mutable S20 : double
val mutable S21 : double
val mutable S22 : double
new (s10, s11, s12, s20, s21, s22) =
{ S10 = s10; S11 = s11; S12 = s12;
S20 = s20; S21 = s21; S22 = s22 }
The implementation of MRG32k3a
gets a reference to State32k3a
structure and updates its
fields. In this way no intermediate object allocation and copy happens, and the computational
performance is close to the C implementation.
let inline private (%.) p m =
let k = int32(p / m)
let r = p  double(k) * m
if r < 0.0 then r + m else r
let MRG32k3a (state : State32k3a byref) =
let p1 = (1403580.0 * state.S11  810728.0 * state.S10) %. 4294967087.0
state.S10 < state.S11; state.S11 < state.S12; state.S12 < p1
let p2 = (527612.0 * state.S21  1370589.0 * state.S20) %. 4294944443.0
state.S20 < state.S21; state.S21 < state.S22; state.S22 < p2
if p1 <= p2 then p1  p2 + 4294967087.0 else p1  p2
Nonetheless, it is possible to create a nice functional wrapper around the low level implementation:
let makeMRG32k2a (seed : int32) =
let seed = float(seed)
let state = ref (State32k3a(seed, seed, seed, seed, seed, seed))
let gen = (fun () > let mutable s = !state in let r = (MRG32k3a &s) in state := s; r)
gen() > ignore; gen() > ignore; // warm up the generator
gen
This implementaion returns a closure that can be used for generating a sequence of random numbers:
let gen = MRG.defaultMRG32k2a
for i in 1..10 do
printfn "%A" (gen())
Though, the price for this nice wrapper is high enough. It involves a heap allocated box for the structure and additional copy operations for reading/writing the structure from the box.
Conslusions
Despite the fact that pure functions and immutable structures are preferable in the most cases because they reduce the number of mistakes, it is important to have clearly designed low level functions for advanced optimizations. Unlike R, Python and other scripting languages, F# has no language border between high level and low level implementations, and, thus, facilitates reuse of low level implementations.

L’Ecuyer Good parameters and implementations for combined multiple recursive random number generators // Operations Research.  1999.  47(1)  pp. 159164. ↩