-- ---
-- title: Random numbers
-- ---
--
-- Generating random numbers in a purely functional language requires
-- us to manually track the state of the random number generator
-- (RNG). The following code shows the idea. Note that the random
-- numbers we will generate are of mediocre quality at best, but
-- they're good enough for small visualisations and simulations and
-- such. Random number generation is a very complex subject, but most
-- of the complexities are the same in any language. Here we focus
-- only on the ones that are particular to Futhark.
--
-- First we define a [module type](abstract-data-types.html)
-- describing the interface of a random number generator.
module type rand = {
-- It consists of a random number generator *state*:
type rng
-- These states can be initialised from a *seed*:
val init : i32 -> rng
-- To generate a uniformly distributed random integer, we pass in an
-- RNG state, which will give us back an integer and a new state:
val rand : rng -> (rng, i32)
-- It is the user's responsibility not to use the same state more than
-- once (unless that is what is desired).
--
-- We will sometimes need to use one seed to generate multiple random
-- streams. The `split` function splits one RNG state
-- in twain:
val split : rng -> (rng, rng)
-- And `split_n` splits it arbitrary many ways, which is useful
-- when we wish to compute many parallel streams in parallel:
val split_n : (n: i64) -> rng -> [n]rng
}
-- We can implement the `rand` module type using various different
-- algorithms. For this example, we'll go with a simple [linear
-- congruential
-- generator](https://en.wikipedia.org/wiki/Linear_congruential_generator)
-- (LCG):
module lcg : rand = {
-- The RNG state is just 32 bits.
type rng = u32
-- An LCG is defined by the three constants *a*, *c*, and *m*. We use
-- the following values, which are the same that are used by glibc:
let a : u32 = 1103515245
let c : u32 = 12345
let m : u32 = 1<<31
-- Note that with these constants, we will never generate a negative
-- random integer.
let rand rng =
let rng' = (a * rng + c) % m
in (rng', i32.u32 rng')
-- Initialisation consists of a few rounds of a [hash function found
-- on
-- StackOverflow](https://stackoverflow.com/questions/664014/what-integer-hash-function-are-good-that-accepts-an-integer-hash-key/12996028#12996028),
-- as is tradition.
let init (x: i32): u32 =
let x = u32.i32 x
let x = ((x >> 16) ^ x) * 0x45d9f3b
let x = ((x >> 16) ^ x) * 0x45d9f3b
let x = ((x >> 16) ^ x)
in x
-- The main weakness of this initialisation is that a zero seed will
-- result in a zero state, which will be "stuck", producing only
-- zeroes. So don't do that.
--
-- The splitting functions creates RNG states by slightly twiddling
-- the states.
let split (rng: rng) =
(init (i32.u32 (rand rng).0),
init (i32.u32 rng))
let split_n n rng =
tabulate n (\i -> init (i32.u32 rng ^ i32.i64 i))
}
-- This is how we might use it:
let test seed =
let rng = lcg.init seed
let (rng, x) = lcg.rand rng
let (_, y) = lcg.rand rng
in (x, y)
-- Usually we don't need integers in the full `i32` range, so let's
-- define a handy function for generating integers in the range `[0,
-- bound-1]`:
let rand_i32 (rng: lcg.rng) (bound: i32) =
let (rng, x) = lcg.rand rng
in (rng, x % bound)
-- Note that the resulting integers are not fully uniformly
-- distributed unless the span of the RNG is divisible by `bound`. In
-- practice, as long as `bound` isn't too large, it should not matter
-- much.
--
-- We can also define a function that gives us a random `f64` in the
-- interval `[0,1]`:
let rand_f64 (rng: lcg.rng) =
let (rng, x) = lcg.rand rng
in (rng, f64.i32 x / f64.i32 i32.highest)
-- For a real random number package, we'd then go on defining various
-- other function for generating normally distributed numbers. We
-- would also write these functions as part of a parametric module,
-- rather than hard-coding the use of `lcg`. All this is done in the
-- [cpprandom](https://github.com/diku-dk/cpprandom) package, which
-- you should use for real Futhark programs.
--
-- # See also
--
-- The [cpprandom](https://github.com/diku-dk/cpprandom) package.