RAND

(Pseudo)random number generation

The library offers a set of random number generators. The random number generation is encapsulated in the class TRandomGenerator in RandomEng.pas that encapsulates various algorithms. There are versions to create integer based random numbers (for a specific range) as well as a floating point version that creates numbers within [0, 1) . This class features follwing algorithms:

• A thread save implementation of the standard Delphi Linear congrential generator algorithm. Thread savety has been achieved by putting the hidden seed variable into a class so there can be more than one instance. Interesting enough the Delphi integer random number generation used the high quadword of the multiplication as random number and keeps the lower quadword as next seed value. The floating point implementation only uses the low quadword to calculate the random floating point number.
• An implementation of the Mersenne Twister algorithm (with improved initialization). There is also a special function available to initialize that generator with an array of keys.
• Utilization of the Intel RDRAND instruction. CAUTION: the implementation is not yet tested and taken from Stackoverflow sources.
• Operating system provided random numbers: This swich utilizes the CryptGenRandom API from Microsoft to generate a set of random numbers. The class uses standard initialization values. There also exists an experimental (untested) OS dependent version for OSX that reads a file stream from /dev/random. Note that this function has not yet been tested so please handle with care.
• A random number generator based on the ChaCha (Salsa20 cipher) encryption algorithm. This algorithim is superfast since it utilises an AVX enhanced calculation.

The library also offers some functions to create different distributed random numbers:

• Gaussian distribution
• Poission
• Exponential distribution
• Erlang distribution

Units

The base random generator can be found in RandomEng.pas. The OS dependent generators in .\win\winRandomGen.pas and the OSX one in .\mac\MacOsRandomGen.pas.

Remarks

• If no hardware implementation is available the class falls back to the standard Delphi Implementation.
• If a seeding value of 0 is given the seed is initialized with the outcome of the current performance counter value
• Freepascal does not use the standard Delphi PRNG so I decided to disable some of the unit tests.
• The matrix library now also includes a call to the random generator class.
• The unit CPUFeatures.pas includes a check for the hardware support of RDRAND
• The standard seed for Mersenne Twister is 5489 . So call the init function if you need real random numbers.

Examples

// calculate 1000 random numbers
procedure TestMersenneTWister;
var gen : TRandomGenerator;
    start, stop : Int64;
    cnt : integer;
begin
     gen := TRandomGenerator.Create;
     gen.RandMethod := raSystem;

     start := MtxGetTime;
     for cnt := 0 to 1000- 1 do
         gen.Random;
     stop := MtxGetTime;

     gen.Free;
     Writeln(Format('Mersenne Twister: %.2f', [(stop - start)/mtxFreq*1000]));
end;

// calculate 1000 numbers in a specific range
procedure TestMersenneTWister;
var gen : TRandomGenerator;
    start, stop : Int64;
    cnt : integer;
begin
     gen := TRandomGenerator.Create;
     gen.RandMethod := raMersenneTwister;

     start := MtxGetTime;
     for cnt := 0 to 1000- 1 do
         gen.RandInt(100);  // values from 0 to 99
     stop := MtxGetTime;

     gen.Free;
     Writeln(Format('Mersenne Twister: %.2f', [(stop - start)/mtxFreq*1000]));
end;

// testing RDRND instruction
uses CPUFeatures;

procedure TestHardwareRandom;
var gen : TRandomGenerator;
    start, stop : Int64;
    cnt : integer;
begin
     if IsHardwareRNDSupport then
     begin
          gen := TRandomGenerator.Create;
          gen.RandMethod := raIntelRAND;

          start := MtxGetTime;
          for cnt := 0 to 1000000 - 1 do
              gen.Random;
          stop := MtxGetTime;

          gen.Free;
          Writeln(Format('Hardware random 1 milliont times: %.2f', [(stop - start)/mtxFreq*1000]));
     end
     else
     begin
          Writeln('No hardware random support detected' );
     end;
 end;

// initialize a matrix:

function InitRandMtx : IMatrix;
begin
     Result := TDoubleMatrix.CreateRand(100, 100, raOs, 0);
end;