blob: bec518af80596a43c710185f17126d48b9d7a555 (
plain)
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
99
100
101
102
103
104
105
106
|
#include "mt19937.h"
#include "mt19937-jump.h"
void mt19937_seed(mt19937_state *state, uint32_t seed) {
int pos;
seed &= 0xffffffffUL;
/* Knuth's PRNG as used in the Mersenne Twister reference implementation */
for (pos = 0; pos < RK_STATE_LEN; pos++) {
state->key[pos] = seed;
seed = (1812433253UL * (seed ^ (seed >> 30)) + pos + 1) & 0xffffffffUL;
}
state->pos = RK_STATE_LEN;
}
/* initializes mt[RK_STATE_LEN] with a seed */
static void init_genrand(mt19937_state *state, uint32_t s) {
int mti;
uint32_t *mt = state->key;
mt[0] = s & 0xffffffffUL;
for (mti = 1; mti < RK_STATE_LEN; mti++) {
/*
* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier.
* In the previous versions, MSBs of the seed affect
* only MSBs of the array mt[].
* 2002/01/09 modified by Makoto Matsumoto
*/
mt[mti] = (1812433253UL * (mt[mti - 1] ^ (mt[mti - 1] >> 30)) + mti);
/* for > 32 bit machines */
mt[mti] &= 0xffffffffUL;
}
state->pos = mti;
return;
}
/*
* initialize by an array with array-length
* init_key is the array for initializing keys
* key_length is its length
*/
void mt19937_init_by_array(mt19937_state *state, uint32_t *init_key,
int key_length) {
/* was signed in the original code. RDH 12/16/2002 */
int i = 1;
int j = 0;
uint32_t *mt = state->key;
int k;
init_genrand(state, 19650218UL);
k = (RK_STATE_LEN > key_length ? RK_STATE_LEN : key_length);
for (; k; k--) {
/* non linear */
mt[i] = (mt[i] ^ ((mt[i - 1] ^ (mt[i - 1] >> 30)) * 1664525UL)) +
init_key[j] + j;
/* for > 32 bit machines */
mt[i] &= 0xffffffffUL;
i++;
j++;
if (i >= RK_STATE_LEN) {
mt[0] = mt[RK_STATE_LEN - 1];
i = 1;
}
if (j >= key_length) {
j = 0;
}
}
for (k = RK_STATE_LEN - 1; k; k--) {
mt[i] = (mt[i] ^ ((mt[i - 1] ^ (mt[i - 1] >> 30)) * 1566083941UL)) -
i; /* non linear */
mt[i] &= 0xffffffffUL; /* for WORDSIZE > 32 machines */
i++;
if (i >= RK_STATE_LEN) {
mt[0] = mt[RK_STATE_LEN - 1];
i = 1;
}
}
mt[0] = 0x80000000UL; /* MSB is 1; assuring non-zero initial array */
}
void mt19937_gen(mt19937_state *state) {
uint32_t y;
int i;
for (i = 0; i < N - M; i++) {
y = (state->key[i] & UPPER_MASK) | (state->key[i + 1] & LOWER_MASK);
state->key[i] = state->key[i + M] ^ (y >> 1) ^ (-(y & 1) & MATRIX_A);
}
for (; i < N - 1; i++) {
y = (state->key[i] & UPPER_MASK) | (state->key[i + 1] & LOWER_MASK);
state->key[i] = state->key[i + (M - N)] ^ (y >> 1) ^ (-(y & 1) & MATRIX_A);
}
y = (state->key[N - 1] & UPPER_MASK) | (state->key[0] & LOWER_MASK);
state->key[N - 1] = state->key[M - 1] ^ (y >> 1) ^ (-(y & 1) & MATRIX_A);
state->pos = 0;
}
extern inline uint64_t mt19937_next64(mt19937_state *state);
extern inline uint32_t mt19937_next32(mt19937_state *state);
extern inline double mt19937_next_double(mt19937_state *state);
void mt19937_jump(mt19937_state *state) { mt19937_jump_state(state); }
|