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
|
#include "mersenne_twister_64.h"
#include "time.h"
enum {
MM = 156,
MATRIX_A = 0xb5026f5aa96619e9ull,
UM = 0xffffffff80000000ull,
LM = 0x7fffffffull
};
void kit_mt64_init_array(kit_mt64_state_t *const state,
ptrdiff_t const size,
uint64_t const *const seed) {
for (ptrdiff_t i = 0; i < size && i < KIT_MT64_N; i++)
state->mt[i] = seed[i];
for (state->index = size; state->index < KIT_MT64_N; state->index++)
state->mt[state->index] = (6364136223846793005ull *
(state->mt[state->index - 1] ^
(state->mt[state->index - 1] >>
62u)) +
state->index);
}
void kit_mt64_init(kit_mt64_state_t *const state,
uint64_t const seed) {
kit_mt64_init_array(state, 1, &seed);
}
uint64_t kit_mt64_generate(kit_mt64_state_t *const state) {
static uint64_t const mag01[2] = { 0ull, MATRIX_A };
int i;
uint64_t x;
if (state->index >= KIT_MT64_N) {
for (i = 0; i < KIT_MT64_N - MM; i++) {
x = (state->mt[i] & UM) | (state->mt[i + 1] & LM);
state->mt[i] = state->mt[i + MM] ^ (x >> 1u) ^
mag01[(int) (x & 1ull)];
}
for (; i < KIT_MT64_N - 1; i++) {
x = (state->mt[i] & UM) | (state->mt[i + 1] & LM);
state->mt[i] = state->mt[i + (MM - KIT_MT64_N)] ^ (x >> 1u) ^
mag01[(int) (x & 1ull)];
}
x = (state->mt[KIT_MT64_N - 1] & UM) | (state->mt[0] & LM);
state->mt[KIT_MT64_N - 1] = state->mt[MM - 1] ^ (x >> 1u) ^
mag01[(int) (x & 1ull)];
state->index = 0;
}
x = state->mt[state->index++];
x ^= (x >> 29u) & 0x5555555555555555ull;
x ^= (x << 17u) & 0x71d67fffeda60000ull;
x ^= (x << 37u) & 0xfff7eee000000000ull;
x ^= (x >> 43u);
return x;
}
uint64_t kit_mt64_seed() {
static uint64_t n = 0;
struct timespec t;
timespec_get(&t, TIME_UTC);
uint64_t seed[2] = { (n++) ^ (uint64_t) t.tv_sec,
(uint64_t) t.tv_nsec };
kit_mt64_state_t s;
kit_mt64_init_array(&s, sizeof seed / sizeof *seed, seed);
return kit_mt64_generate(&s);
}
|
