[M3devel] Mersenne Twister in Modula-3, Error?

[Rodney M. Bates]

Can you post the C code you are comparing against?

On 10/18/2011 07:16 PM, Ken Durocher wrote:
> I implemented the Mersenne Twister in Modula-3 some time ago, however when I test it out (the same test data that the authors of the Mersenne Twister used), the output is different.  I tested this code on both 32 and 64 bit machines, and they both give me the same results, but the results are different from that of the original Mersenne Twister code.
>
> Here's the code:
>
> INTERFACE MT19937ar;
>
> IMPORT Word;
>
> PROCEDURE initGenrand(s: Word.T);
> PROCEDURE initByArray(key: ARRAY OF Word.T; length: INTEGER);
>
> PROCEDURE genrandInt32(): Word.T;
> PROCEDURE genrandInt31(): INTEGER;
> PROCEDURE genrandReal1(): REAL;
> PROCEDURE genrandReal2(): REAL;
> PROCEDURE genrandReal3(): REAL;
>
> END MT19937ar.
>
> (*********)
>
> MODULE MT19937ar;
>
> IMPORT Word;
> FROM Word IMPORT And, Xor, Or, RightShift, LeftShift;
>
> CONST
>    n = 624;
>    m = 397;
>    matrix_a = 16_9908b0df;
>    upper_mask = 16_80000000;
>    lower_mask = 16_7fffffff;
>
> VAR
>    mt: ARRAY [0..n] OF Word.T;
>    mti := n + 1;
>
> PROCEDURE initGenrand(s: Word.T) =
>    BEGIN
>      mt[0] := And(s, 16_ffffffff);
>      FOR mti := 1 TO n - 1 DO
>        mt[mti] := 1812433253 * (Xor(mt[mti-1], RightShift(mt[mti-1], 30))) + mti;
>        mt[mti] := And(mt[mti], 16_ffffffff);
>      END;
>    END initGenrand;
>
> PROCEDURE initByArray(key: ARRAY OF Word.T; length: INTEGER) =
>    VAR i, j, k: INTEGER;
>    BEGIN
>      i := 1; j := 0;
>      initGenrand(19650218);
>      IF n > length THEN k := n ELSE k := length; END;
>      WHILE k # 0 DO
>        mt[i] := Xor(mt[i], Xor(mt[i-1], RightShift(mt[i-1], 30)) * 1664525) +
>                     key[j] + j;
>        mt[i] := And(mt[i], 16_ffffffff);
>        INC(i); INC(j);
>        IF i >= n THEN
>          mt[0] := mt[n-1];
>          i := 1;
>        END;
>
>        IF j >= length THEN
>          j := 0;
>        END;
>        DEC(k);
>      END;
>
>      k := n - 1;
>      WHILE k # 0 DO
>        mt[i] := Xor(mt[i], Xor(mt[i-1], RightShift(mt[i-1], 30)) * 1566083941) - i;
>        mt[i] := And(mt[i], 16_ffffffff);
>        INC(i);
>        IF i >= n THEN
>          mt[0] := mt[n-1];
>          i := 1;
>        END;
>        DEC(k);
>      END;
>
>      mt[0] := 16_80000000;
>    END initByArray;
>
> PROCEDURE genrandInt32(): Word.T =
>    VAR y: Word.T;
>        k: INTEGER := 0;
>        mag01 := ARRAY [0..1] OF Word.T {0, matrix_a};
>    BEGIN
>      IF mti >= n THEN
>        IF mti = (n + 1) THEN
>          initGenrand(5489);
>        END;
>
>        WHILE k < (n - m) DO
>          y := Or(And(mt[k], upper_mask), And(mt[k+1], lower_mask));
>          mt[k] := Xor(Xor(mt[k+m], RightShift(y, 1)), mag01[And(y, 1)]);
>          INC(k);
>        END;
>        WHILE k < (n - 1) DO
>          y := Or(And(mt[k], upper_mask), And(mt[k+1], lower_mask));
>          mt[k] := Xor(Xor(mt[k+(m-n)], RightShift(y, 1)), mag01[And(y, 1)]);
>          INC(k);
>        END;
>
>        y := Or(And(mt[n-1], upper_mask), And(mt[0], lower_mask));
>        mt[n-1] := Xor(Xor(mt[m-1], RightShift(y, 1)), mag01[And(y, 1)]);
>        mti := 0;
>      END;
>
>      y := mt[mti];
>      INC(mti);
>
>      y := Xor(y, RightShift(y, 11));
>      y := Xor(y, And(LeftShift(y, 7), 16_9d2c5680));
>      y := Xor(y, And(LeftShift(y, 15), 16_efc60000));
>      y := Xor(y, RightShift(y, 18));
>
>      RETURN y;
>    END genrandInt32;
>
> PROCEDURE genrandInt31(): INTEGER =
>    BEGIN
>      RETURN RightShift(genrandInt32(), 1);
>    END genrandInt31;
>
> PROCEDURE genrandReal1(): REAL =
>    BEGIN
>      RETURN FLOAT(genrandInt32(), REAL) * (1.0 / 4294967295.0);
>    END genrandReal1;
>
> PROCEDURE genrandReal2(): REAL =
>    BEGIN
>      RETURN FLOAT(genrandInt32(), REAL) * (1.0 / 4294967296.0);
>    END genrandReal2;
>
> PROCEDURE genrandReal3(): REAL =
>    BEGIN
>      RETURN (FLOAT(genrandInt32(), REAL) + 0.5) * (1.0 / 4294967296.0);
>    END genrandReal3;
>
> BEGIN
> END MT19937ar.
>
> (********)
>
> MODULE Main;
>
> IMPORT IO, Fmt, MT19937ar;
>
> VAR init := ARRAY [0..3] OF INTEGER {16_123, 16_234, 16_345, 16_456};
>      len := 4;
>
> BEGIN
>    MT19937ar.initByArray(init, len);
>    IO.Put("Random integers: ");
>    FOR i := 1 TO 999 DO
>      IO.Put(Fmt.Unsigned(MT19937ar.genrandInt32(), base := 10));
>      IO.Put(" ");
>      IF (i MOD 5) = 4 THEN IO.Put("\n"); END;
>    END;
>    IO.Put("\n");
> END Main.
>
> And here is the first few numbers it outputs:
>
> 3499211612 581869302 3890346734 3586334585 ...
>
> Compare that to the first few numbers of output from the original Mersenne Twister code, in C:
>
> 1067595299 955945823 477289528 4107218783 ...
>
> Anyone know why they are different? Is the code wrong, or is it an implementation difference, or what?