tango.math.random.Random

Random number generators This is an attempt at having a good flexible and easy to use random number generator. ease of use:
  • shared generator for quick usage available through the "rand" object 
    1

    int i=rand.uniformR(10); // a random number from [0;10)
  • simple Random (non threadsafe) and RandomSync (threadsafe) types to create new generators (for heavy use a good idea is one Random object per thread)
  • several distributions can be requested like this 
    1

    rand.distributionD!(type)(paramForDistribution)
    the type can often be avoided if the parameters make it clear. From it single numbers can be generated with .getRandom(), and variables initialized either with call style (var) or with .randomize(var). Utility functions to generate numbers directly are also available. The choice to put all the distribution in a single object that caches them has made (for example) the gamma distribution very easy to implement.
  • sample usage: 
    1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58

    auto r=new Random();
int i; float f; real rv; real[100] ar0; real[] ar=ar0[];
// initialize with uniform distribution
i=r.uniform!(int);
f=r.uniform!(float);
rv=r.uniform!(real);
foreach (ref el;ar)
  el=r.uniform!(real);
// another way to do all the previous in one go:
r(i)(f)(rv)(ar);
// unfortunetely one cannot use directly ar0...
// uniform distribution 0..10
i=r.uniformR(10);
f=r.uniformR(10.0f);
rv=r.uniformR(10.0L);
foreach (ref el;ar)
  el=r.uniformR(10.0L);
// another way to do all the previous in one go:
r.uniformRD(10)(i)(f)(r)(ar);
// uniform numbers in [5;10)
i=r.uniformR2(5,10);
// uniform numbers in (5;10)
f=r.uniformR2(5.0f,10.0f);
rv=r.uniformR2(5.0L,10.0L);
foreach (ref el;ar)
  el=r.uniformR2(5.0L,10.0L);
// another way to do all the previous in one go:
r.uniformR2D(5.0L,10.0L)(i)(f)(r)(ar);
// uniform distribution -10..10
i=r.uniformRSymm(10);
// well you get it...
r.uniformRSymmD(10)(i)(f)(r)(ar);
// any distribution can be stored
auto r2=r.uniformRSymmD(10);
// and used later
r2(ar);
// complex distributions (normal,exp,gamma) are produced for the requested type
r.normalSource!(float)()(f);
// with sigma=2
r.normalD(2.0f)(f);
// and can be used also to initialize other types
r.normalSource!(float)()(r)(ar);
r.normalD(2.0f)(r)(ar);
// but this is different from
r.normalSource!(real)()(i)(r)(ar);
r.normalD(2.0L)(i)(r)(ar);
// as the source generates numbers of its type that then are simply cast to
// the type needed.
// Uniform distribution (as its creation for different types has no overhead)
// is never cast, so that (for example) bounds exclusion for floats is really
// guaranteed.
// For the other distribution using a distribution of different type than
// the variable should be done with care, as underflow/overflow might ensue.
//
// Some utility functions are also available
int i2=r.uniform!(int)();
int i2=r.randomize(i); // both i and i2 are initialized to the same value
float f2=r.normalSigma(3.0f);

flexibility:

  • easily swappable basic source 
    1
2

    // a random generator that uses the system provided random generator:
auto r=RandomG!(Urandom)();
    One could also build an engine that can be changed at runtime (that calls a delegate for example), but this adds a little overhead, and changing engine is not something done often, so this is not part of the library.
  • ziggurat generator can be easily adapted to any decreasing derivable distribution, the hard parametrization (to find xLast) can be done automatically
  • several distributions available "out of the box"

Quality:

  • the default Source combines two surces that pass all statistical tests (KISS+CMWC) (P. L'Ecuyer and R. Simard, ACM Transactions on Mathematical Software (2007), 33, 4, Article 22, for KISS, see CMWC engine for the other)
  • floating point uniform generator always initializes the full mantissa, the only flaw is a (*very* small) predilection of 0 as least important bit (IEEE rounds to 0 in case of tie). Using a method that initializes the full mantissa was shown to improve the quality of subsequntly derived normal distribued numbers (Thomas et al. Gaussian random number generators. Acm Comput Surv (2007) vol. 39 (4) pp. 11)
  • Ziggurat method, a very fast and accurate method was used for both Normal and exp distributed numbers.
  • gamma distribued numbers uses a method recently proposed by Marsaglia and Tsang. The method is very fast, and should be good. My (Fawzi's) feeling is that the transformation h(x)=(1+d*x)^3 might lose a couple of bits of precision in some cases, but it is unclear if this might become visible in (*very* extensive) tests or not.
    • the basic source can be easily be changed with something else

    Efficiency:

  • very fast methods have been used, and some effort has been put into optimizing some of them, but not all, but the interface has been choosen so that close to optimal implementation can be provided through the same interface.
  • Normal and Exp sources allocated only upon request: no memory waste, but a (*very* small) speed hit, that can be avoided by storing the source in a variable and using it (not going through the RandomG)
    • )

    Annoyances:

    • I have added two "next" methods to RandomG for backward compatibility reasons, and the .instance from Random has been replaced by the "rand" object. The idea behind this is that RandomG is a template and rand it should be shared across all templates. If the name rand is considered bad one could change it. I kept .instance static method that returns rand, so this remain a dropin replacement of the old random.
    • You cannot initialize a static array directly, this because randomize is declared like this: 
      1

      U randomize(U)(ref U a) { }
      and a static array cannot be passed by reference. Removing the ref would make arrays initialized, and scalar not, which is much worse.

    License:

    BSD style: see license.txt

    Version:

    Initial release: July 2008

    Author:

    Fawzi Mohamed
    template isFloat(T)
    if T is a float
    alias KissCmwc_32_1 DefaultEngine
    The default engine, a reasonably collision free, with good statistical properties not easy to invert, and with a relatively small key (but not too small)
    class RandomG(SourceT = DefaultEngine) [final]
    Class that represents a random number generator. Normally you should get random numbers either with call-like interface: auto r=new Random(); r(i)(j)(k); or with randomize r.randomize(i); r.randomize(j); r.randomize(k); if you use this you should be able to easily switch distribution later, as all distributions support this interface, and can be built on the top of RandomG auto r2=r.NormalSource!(float)(); r2(i)(j)(k); there are utility methods within random for the cases in which you do not want to build a special distribution for just a few numbers
    this(bool randomInit = true)
    Creates and seeds a new generator
    RandomG seed()
    if source.canSeed seeds the generator using the shared rand generator (use urandom directly if available?)
    RandomG seed(scope uint delegate() seedSource)
    if source.canSeed seeds the generator using the given source of uints
    uint next()
    uint next(uint to)
    uint next(uint from, uint to)
    RandomG!(Sync!(DefaultEngine)) instance() [static]
    compatibility with old Random, deprecate??
    T uniform(T, bool boundCheck = true)() [@property]
    uniform distribution on the whole range of integer types, and on the (0;1) range for floating point types. Floating point guarantees the initialization of the full mantissa, but due to rounding effects it might have *very* small dependence due to rounding effects on the least significant bit (in case of tie 0 is favored). if boundCheck is false in the floating point case bounds might be included (but with a lower propability than other numbers)
    T uniformR(T, bool boundCheck = true)(T to)
    uniform distribution on the range [0;to) for integer types, and on the (0;to) range for floating point types. Same caveat as uniform(T) apply
    T uniformRSymm(T, bool boundCheck = true, bool excludeZero = isFloat!(T))(T to, int iter = 2000)
    uniform distribution on the range (-to;to) for integer types, and on the (-to;0)(0;to) range for floating point types if boundCheck is true. If boundCheck=false the range changes to [-to;0)u(0;to] with a slightly lower propability at the bounds for floating point numbers. excludeZero controls if 0 is excluded or not (by default float exclude it, ints no). Please note that the probability of 0 in floats is very small due Cannot be used on unsigned types.
    In here there is probably one of the few cases where c handling of modulo of negative numbers is handy
    T uniformR2(T, bool boundCheck = true)(T from, T to)
    uniform distribution [from;to) for integers, and (from;to) for floating point numbers. if boundCheck is false the bounds are included in the floating point number distribution. the range for int and long is limited to only half the possible range (it could be worked around using long aritmethic for int, and doing a carry by hand for long, but I think it is seldomly needed, for int you are better off using long when needed)
    T uniformEl(T)(const(T[]) arr)
    returns a random element of the given array (which must be non empty)
    U randomizeUniform(U, bool boundCheck)(ref U a)
    randomizes the given array and returns it (for some types this is potentially more efficient, both from the use of random numbers and speedwise)
    U randomizeUniformR(U, V, bool boundCheck = true)(ref U a, V to)
    randomizes the given array and returns it (for some types this is potentially more efficient, both from the use of random numbers and speedwise)
    U randomizeUniformR2(U, V, W, bool boundCheck = true)(ref U a, V from, W to)
    randomizes the given variable and returns it (for some types this is potentially more efficient, both from the use of random numbers and speedwise)
    U randomizeUniformRSymm(U, V, bool boundCheck = true, bool excludeZero = isFloat!(BaseTypeOfArrays!(U)))(ref U a, V to)
    randomizes the given variable like uniformRSymm and returns it (for some types this is potentially more efficient, both from the use of random numbers and speedwise)
    RandG spawn(RandG = RandomG)()
    returns another (mostly indipendent, depending on seed size) random generator
    struct UniformDistribution(T, bool boundCheck)
    uniform distribution on the whole range for integers, and on (0;1) for floats with boundCheck=true this is equivalent to r itself, here just for completness
    UniformDistribution opCall(U, S...)(ref U a, S args)
    chainable call style initialization of variables (thorugh a call to randomize)
    T getRandom()
    returns a random number
    U randomize(U)(ref U a)
    initialize el
    struct UniformRDistribution(T, bool boundCheck)
    uniform distribution on the subrange [0;to) for integers, (0;to) for floats
    UniformRDistribution create(RandomG r, T to) [static]
    initializes the probability distribution
    UniformRDistribution opCall(U)(ref U a)
    chainable call style initialization of variables (thorugh a call to randomize)
    T getRandom()
    returns a random number
    U randomize(U)(ref U a)
    initialize el
    struct UniformRSymmDistribution(T, bool boundCheck = true, bool excludeZero = isFloat!(T))
    uniform distribution on the subrange (-to;to) for integers, (-to;0)u(0;to) for floats excludeZero controls if the zero should be excluded, boundCheck if the boundary should be excluded for floats
    UniformRSymmDistribution create(RandomG r, T to) [static]
    initializes the probability distribution
    UniformRSymmDistribution opCall(U)(ref U a)
    chainable call style initialization of variables (thorugh a call to randomize)
    T getRandom()
    returns a random number
    U randomize(U)(ref U a)
    initialize el
    struct UniformR2Distribution(T, bool boundCheck)
    uniform distribution on the subrange (-to;to) for integers, (0;to) for floats
    UniformR2Distribution create(RandomG r, T from, T to) [static]
    initializes the probability distribution
    UniformR2Distribution opCall(U, S...)(ref U a, S args)
    chainable call style initialization of variables (thorugh a call to randomize)
    T getRandom()
    returns a random number
    U randomize(U)(ref U a)
    initialize a
    struct GammaDistribution(T)
    gamma distribution f=x^(alpha-1)*exp(-x/theta)/(gamma(alpha)*theta^alpha) alpha has to be bigger than 1, for alpha<1 use gammaD(alpha)=gammaD(alpha+1)*pow(r.uniform!(T),1/alpha) from Marsaglia and Tsang, ACM Transaction on Mathematical Software, Vol. 26, N. 3 2000, p 363-372
    GammaDistribution opCall(U, S...)(ref U a, S args)
    chainable call style initialization of variables (thorugh a call to randomize)
    T getRandom(T a = alpha, T t = theta)
    returns a single random number
    U randomize(U)(ref U b, T a = alpha, T t = theta)
    initializes b with gamma distribued random numbers
    U randomizeOp(U, S)(S delegate(T) op, ref U b, T a = alpha, T t = theta)
    maps op on random numbers (of type T) and initializes b with it
    NormalSource!(RandomG,T) normalSource(T)()
    generators of normal numbers (sigma=1,mu=0) of the given type f=exp(-x*x/(2*sigma^2))/(sqrt(2 pi)*sigma)
    ExpSource!(RandomG,T) expSource(T)()
    generators of exp distribued numbers (beta=1) of the given type f=1/beta*exp(-x/beta)
    NormalSource!(RandomG,T).NormalDistribution normalD(T)(T sigma = cast(T)1, T mu = cast(T)0)
    generators of normal numbers with a different default sigma/mu f=exp(-x*x/(2*sigma^2))/(sqrt(2 pi)*sigma)
    ExpSource!(RandomG,T).ExpDistribution expD(T)(T beta)
    exponential distribued numbers with a different default beta f=1/beta*exp(-x/beta)
    GammaDistribution!(T) gammaD(T)(T alpha = cast(T)1, T theta = cast(T)1)
    gamma distribued numbers with the given default alpha
    UniformDistribution!(T,true) uniformD(T)()
    uniform distribution on the whole integer range, and on (0;1) for floats should return simply this??
    UniformDistribution!(T,false) uniformBoundsD(T)()
    uniform distribution on the whole integer range, and on [0;1] for floats
    UniformRDistribution!(T,true) uniformRD(T)(T to)
    uniform distribution [0;to) for ints, (0:to) for reals
    UniformRDistribution!(T,false) uniformRBoundsD(T)(T to)
    uniform distribution [0;to) for ints, [0:to] for reals
    UniformRSymmDistribution!(T,true,isFloat!(T)) uniformRSymmD(T)(T to)
    uniform distribution (-to;to) for ints and (-to;0)u(0;to) for reals
    UniformRSymmDistribution!(T,false,isFloat!(T)) uniformRSymmBoundsD(T)(T to)
    uniform distribution (-to;to) for ints and [-to;0)u(0;to] for reals
    UniformR2Distribution!(T,true) uniformR2D(T)(T from, T to)
    uniform distribution [from;to) for ints and (from;to) for reals
    UniformR2Distribution!(T,false) uniformR2BoundsD(T)(T from, T to)
    uniform distribution [from;to) for ints and [from;to] for reals
    T normal(T)()
    returns a normal distribued number
    T normalSigma(T)(T sigma)
    returns a normal distribued number with the given sigma
    T normalSigmaMu(T)(T sigma, T mu)
    returns a normal distribued number with the given sigma and mu
    T exp(T)()
    returns an exp distribued number
    T expBeta(T)(T beta)
    returns an exp distribued number with the given scale beta
    T gamma(T)(T alpha = cast(T)1, T sigma = cast(T)1)
    returns a gamma distribued number from Marsaglia and Tsang, ACM Transaction on Mathematical Software, Vol. 26, N. 3 2000, p 363-372
    string toString() [override]
    writes the current status in a string
    size_t fromString(const(char[]) s)
    reads the current status from a string (that should have been trimmed) returns the number of chars read
    RandomG opCall(U)(ref U a)
    chainable call style initialization of variables (thorugh a call to randomize)
    T getRandom(T)()
    returns a random number
    U randomize(U)(ref U a)
    initialize el
    alias RandomG!() Random
    make the default random number generator type (a non threadsafe random number generator) easily available you can safely expect a new instance of this to be indipendent from all the others
    alias RandomG!(Sync!(DefaultEngine)) RandomSync
    default threadsafe random number generator type
    RandomSync rand [gshared, static]
    shared locked (threadsafe) random number generator initialized with urandom if available, with time otherwise