aboutsummaryrefslogtreecommitdiffstats
path: root/util/random/mersenne32.cpp
blob: a7a4d17b2c49b798b74e1e5c6edec36fb56e7eb2 (plain) (blame)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
#include "mersenne32.h" 
 
#include <util/generic/array_size.h> 
#include <util/stream/input.h> 
 
using namespace NPrivate; 
 
#define M 397 
#define MATRIX_A 0x9908b0dfUL 
#define UPPER_MASK 0x80000000UL 
#define LOWER_MASK 0x7fffffffUL 
 
void TMersenne32::InitGenRand(ui32 s) noexcept {
    mt[0] = s; 
 
    for (mti = 1; mti < N; ++mti) { 
        mt[mti] = (1812433253UL * (mt[mti - 1] ^ (mt[mti - 1] >> 30)) + mti); 
    } 
} 
 
void TMersenne32::InitByArray(const ui32 init_key[], size_t key_length) noexcept {
    InitGenRand(19650218UL); 
 
    ui32 i = 1; 
    ui32 j = 0; 
    ui32 k = ui32(N > key_length ? N : key_length); 
 
    for (; k; k--) { 
        mt[i] = (mt[i] ^ ((mt[i - 1] ^ (mt[i - 1] >> 30)) * 1664525UL)) + init_key[j] + j; 
 
        ++i; 
        ++j; 
 
        if (i >= N) { 
            mt[0] = mt[N - 1]; 
            i = 1; 
        } 
 
        if (j >= key_length) { 
            j = 0; 
        } 
    } 
 
    for (k = N - 1; k; k--) { 
        mt[i] = (mt[i] ^ ((mt[i - 1] ^ (mt[i - 1] >> 30)) * 1566083941UL)) - i; 
 
        ++i; 
 
        if (i >= N) { 
            mt[0] = mt[N - 1]; 
            i = 1; 
        } 
    } 
 
    mt[0] = 0x80000000UL; 
} 
 
void TMersenne32::InitNext() noexcept {
    int kk; 
    ui32 y; 
    ui32 mag01[2] = { 
        0x0UL, 
        MATRIX_A, 
    }; 
 
    if (mti == N + 1) { 
        InitGenRand(5489UL); 
    } 
 
    for (kk = 0; kk < N - M; ++kk) { 
        y = (mt[kk] & UPPER_MASK) | (mt[kk + 1] & LOWER_MASK); 
        mt[kk] = mt[kk + M] ^ (y >> 1) ^ mag01[y & 0x1UL]; 
    } 
 
    for (; kk < N - 1; ++kk) { 
        y = (mt[kk] & UPPER_MASK) | (mt[kk + 1] & LOWER_MASK); 
        mt[kk] = mt[kk + (M - N)] ^ (y >> 1) ^ mag01[y & 0x1UL]; 
    } 
 
    y = (mt[N - 1] & UPPER_MASK) | (mt[0] & LOWER_MASK); 
    mt[N - 1] = mt[M - 1] ^ (y >> 1) ^ mag01[y & 0x1UL]; 
 
    mti = 0; 
} 
 
TMersenne32::TMersenne32(IInputStream& input)
    : mti(N + 1) 
{ 
    ui32 buf[128]; 
 
    input.LoadOrFail(buf, sizeof(buf)); 
    InitByArray(buf, Y_ARRAY_SIZE(buf));
} 
 
#undef M 
#undef MATRIX_A 
#undef UPPER_MASK 
#undef LOWER_MASK