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

From Lazarus wiki
Jump to navigationJump to search
m
Line 107: Line 107:
 
Left FPC, right C++<br>
 
Left FPC, right C++<br>
 
[[File:Selection_019.png]]
 
[[File:Selection_019.png]]
 +
== verification against FPC's Random ==
 +
<syntaxhighlight>
 +
program mersennetest2;
 +
{$mode objfpc}
 +
uses mersenne;
 +
var
 +
  i:integer;
 +
begin
 +
  Randseed :=5489;
 +
  for i := 0 to 4 do write(mer32.random * 2.32830643653869628906e-10 :2:10,' ');writeln('Mersenne32');
 +
  for i := 0 to 4 do write(single(Random) :2:10,' ');writeln('FPC Random'); 
 +
end.
 +
</syntaxhighlight>
  
 
== Reference ==
 
== Reference ==

Revision as of 23:15, 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 a simple 32 bit implementation according to MT19937. It is compatible with FPC's Random for 32 bit values. (integer,cardinal/dword, single)

Enjoy,
Thaddy.

It is basically taken from the C version in wikipedia and gererates the same output. It is also verified against the mt19937 class in GNU C++11.

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

implementation 

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

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++
Selection 019.png

verification against FPC's Random

program mersennetest2;
{$mode objfpc}
uses mersenne;
var 
  i:integer;
begin
   Randseed :=5489;
   for i := 0 to 4 do write(mer32.random * 2.32830643653869628906e-10 :2:10,' ');writeln('Mersenne32');
   for i := 0 to 4 do write(single(Random) :2:10,' ');writeln('FPC Random');  
end.

Reference

  1. Mersenne twister in WikiPedia

See also