A simple implementation of the Mersenne twister

TMersenne32

When you look at the source code that implements Free Pascal's random it may look overwelmingly complex. The algorithm itself, though, is very understandable.
Here is a simple 32 bit implementation according to MT19937. It is compatible with FPC's Random for 32 bit values. It is also about 45% faster than FPC's random.

Enjoy,
Thaddy.

It is verified against several mt19937 implementations in different languages, including Freepascal and C/C++.

Of course I added this one to wikipedia too. You don't have to create this class. it is fully static.
Simply call Rnd32.URand, URandf, Srand of Srandf to obtain a random value.

{ new version 2017-02-12 }
unit mersenne;
{$mode objfpc}
interface
type
   TMersenne32 = class sealed
   strict private
   const 
     // Define MT19937 constants (32-bit RNG)
     N = 624;M = 397;R = 31;A = $9908B0DF;F = 1812433253;
     U = 11;S = 7;B = $9D2C5680;T = 15;C = $EFC60000;L = 18;
     MASK_LOWER = (QWORD(1) << R) - 1;
     MASK_UPPER = QWORD(1) << R;
     class var mt:array[0..N-1] of dword;
     class var index:word;
     class procedure twist;inline;static;
   public
     class constructor create;
     class procedure initialize(const seed:dword);inline;static;
     { 32 bit unsigned integer }
     class function URand:dword;inline;static;
     { 32 bit float in range 0..1 }
     class function URandf:single;inline;static;
     { 32 bit signed integer }
     class function SRand:integer;inline;static;
     {32 bit float in the range -1..1 }
     class function SRandf:single;inline;static;
   end;
   
   Rnd32 = type TMersenne32;
   
implementation 

class constructor TMersenne32.Create;
begin
  initialize(5489); // default seed for mt199937
end;

class procedure TMersenne32.Initialize(const seed:dword);inline;static;
var 
  i:dword;
begin
  mt[0] := seed;
 for  i := 1 to pred(N) do 
   mt[i] := F * (mt[i - 1] xor (mt[i - 1] >> 30)) + i;
 index := N;
end;
   
class procedure TMersenne32.Twist;inline;static;
var 
  i,x,xA:dword;
begin
  for i := 0 to pred(N) do
  begin
    x := (mt[i] and MASK_UPPER) + (mt[(i + 1) mod N] and MASK_LOWER);
    xA := x >> 1;
    if ( x and $1 ) <> 0 then
	  xA := xA xor A;
    mt[i] := mt[(i + M) mod N] xor xA;
  end;
  index := 0;
end;

class function TMersenne32.URand:dword;inline;static;
var
  y:dword;
  i:integer;
begin
  i := index;
  if  index >= N then
  begin
    Twist;
    i := index;
  end;
  y := mt[i];index := i + 1;
  y := y xor (mt[i] >> U);
  y := y xor (y << S) and B;
  y := y xor (y << T) and C;
  y := y xor(y >> L);
  Result :=y;
end;

class function TMersenne32.URandf:single;inline;static;
begin
  Result := URand * 2.32830643653869628906e-10;  
end;

class function TMersenne32.SRandf:single;inline;static;
begin
  Result :=URand * 4.6566128730773926e-010;  
end;

class function TMersenne32.SRand:integer;inline;static;
begin
  {$ifopt D+}Result := 0;{$endif}
  Result := URand;  
end;
end.

C++ verification test:

Compile this with g++ mersplusplus.cpp -std=gnu++11 -o mersplusplus
Run ./mersplusplus and compare the output with TMersenne32. It should be (is) equal.

#include <iostream>
#include <random> 
using namespace std;
 int main()
{
    mt19937 mt_rand(5489u);
    for(int i=0;i<5;i++){
    cout << mt_rand() << endl;};
 
    return 0;
}

Left FPC, right C++

Verification against FPC's Random

program mersennetest2;
{$mode objfpc}
uses sysutils,dateutils,mersenne;
var 
  i:integer;
  start:TdateTime;
  m,f:int64;
begin
   Randseed :=5489; // default seed for mt199937
   start := Now;for i := 0 to 9999999 do mer32.urandf;write(mer32.urandf :2:10,' ');writeln('Mersenne32');m:=MilliSecondsBetween(Now,start);
   for i := 0 to 9999999 do Random ;write(single(Random) :2:10,' ');writeln('FPC Random');f:=MilliSecondsBetween(Now,start);  
   writeln('Mersenne32: ',m,' FPC Random ',f, ', speed difference (f *100)/m: ', 'Mersenne32 is ',f*100/m -100 :4:0 ,'% Faster than FPC's Random');
end.

Reference

1. Mersenne twister in WikiPedia

See also

• Dev random
• Functions for descriptive statistics
• Generating Random Numbers
• Marsaglia's pseudo random number generators