Difference between revisions of "A simple implementation of the Mersenne twister"

From Free Pascal wiki
Jump to navigationJump to search
Line 7: Line 7:
 
It is basically taken from the C version in wikipedia and gererates the same output.<br>  
 
It is basically taken from the C version in wikipedia and gererates the same output.<br>  
 
Of course I added this one to wikipedia too. You don't have to create this class. it is fully static.<br>
 
Of course I added this one to wikipedia too. You don't have to create this class. it is fully static.<br>
Simply call Mer32.ExtractU32 to obtain a random value. Cast to integer for signed output.
+
Simply call Mer32.Random to obtain a random value. Cast to integer for signed output.
  
  
Line 15: Line 15:
 
interface
 
interface
 
type
 
type
   TMersenne = class sealed
+
   TMersenne32 = class sealed
 
   strict private
 
   strict private
 
   const  
 
   const  
Line 36: Line 36:
 
implementation  
 
implementation  
  
class constructor TMersenne.Create;
+
class constructor TMersenne32.Create;
 
begin
 
begin
 
   initialize(12345678);
 
   initialize(12345678);
 
end;
 
end;
  
class procedure TMersenne.Initialize(const seed:dword);inline;static;
+
class procedure TMersenne32.Initialize(const seed:dword);inline;static;
 
var  
 
var  
 
   i:dword;
 
   i:dword;
Line 51: Line 51:
 
end;
 
end;
 
    
 
    
class procedure TMersenne.Twist;inline;static;
+
class procedure TMersenne32.Twist;inline;static;
 
var  
 
var  
 
   i,x,xA:dword;
 
   i,x,xA:dword;
Line 66: Line 66:
 
end;
 
end;
  
class function TMersenne.ExtractU32:dword;inline;static;
+
class function TMersenne32.Random:dword;inline;static;
 
var
 
var
 
   y:dword;
 
   y:dword;

Revision as of 13:42, 12 March 2017

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 simple 32 bit implementation according to MT19937.

Enjoy,
Thaddy.

It is basically taken from the C version in wikipedia and gererates the same output.
Of course I added this one to wikipedia too. You don't have to create this class. it is fully static.
Simply call Mer32.Random to obtain a random value. Cast to integer for signed output.


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;
     class function Extractu32:dword;inline;static;
   end;
   
   mer32 = type TMersenne;

implementation 

class constructor TMersenne32.Create;
begin
  initialize(12345678);
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.Random: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;

end.

Reference

  1. Mersenne twister in WikiPedia

See also