diff options
-rw-r--r-- | source/kit/CMakeLists.txt | 4 | ||||
-rw-r--r-- | source/kit/mersenne_twister_64.c | 62 | ||||
-rw-r--r-- | source/kit/mersenne_twister_64.h | 37 | ||||
-rw-r--r-- | source/test/unittests/CMakeLists.txt | 3 | ||||
-rw-r--r-- | source/test/unittests/mersenne_twister_64.test.c | 33 |
5 files changed, 137 insertions, 2 deletions
diff --git a/source/kit/CMakeLists.txt b/source/kit/CMakeLists.txt index 9cc1f7e..5891851 100644 --- a/source/kit/CMakeLists.txt +++ b/source/kit/CMakeLists.txt @@ -3,6 +3,7 @@ target_sources( PRIVATE atomic.c input_buffer.c input_stream.c lower_bound.c string_ref.c async_function.c allocator.c array_ref.c dynamic_array.c + mersenne_twister_64.c PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}/atomic.h> $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}/allocator.h> @@ -12,4 +13,5 @@ target_sources( $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}/input_stream.h> $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}/input_buffer.h> $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}/lower_bound.h> - $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}/array_ref.h>) + $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}/array_ref.h> + $<BUILD_INTERFACE:${CMAKE_CURRENT_LIST_DIR}/mersenne_twister_64.h>) diff --git a/source/kit/mersenne_twister_64.c b/source/kit/mersenne_twister_64.c new file mode 100644 index 0000000..4be4313 --- /dev/null +++ b/source/kit/mersenne_twister_64.c @@ -0,0 +1,62 @@ +#include "mersenne_twister_64.h" + +#include <time.h> + +void kit_mt64_init(kit_mt64_state_t *state, uint64_t seed) { + state->mt[0] = seed; + for (state->index = 1; 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); +} + +uint64_t kit_mt64_generate(kit_mt64_state_t *state) { + static uint64_t const mag01[2] = { 0ull, 0xB5026F5AA96619E9ull }; + + int i; + uint64_t x; + + if (state->index >= KIT_MT64_N) { + for (i = 0; i < KIT_MT64_N - 156; i++) { + x = (state->mt[i] & 0xFFFFFFFF80000000ull) | + (state->mt[i + 1] & 0x7FFFFFFFull); + state->mt[i] = state->mt[i + 156] ^ (x >> 1u) ^ + mag01[(int) (x & 1ull)]; + } + + for (; i < KIT_MT64_N - 1; i++) { + x = (state->mt[i] & 0xFFFFFFFF80000000ull) | + (state->mt[i + 1] & 0x7FFFFFFFull); + state->mt[i] = state->mt[i + (156 - KIT_MT64_N)] ^ (x >> 1u) ^ + mag01[(int) (x & 1ull)]; + } + + x = (state->mt[KIT_MT64_N - 1] & 0xFFFFFFFF80000000ull) | + (state->mt[0] & 0x7FFFFFFFull); + state->mt[KIT_MT64_N - 1] = state->mt[156 - 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() { + struct timespec t; + timespec_get(&t, TIME_UTC); + + kit_mt64_state_t s; + kit_mt64_init(&s, (uint64_t) t.tv_nsec); + + return kit_mt64_generate(&s); +} diff --git a/source/kit/mersenne_twister_64.h b/source/kit/mersenne_twister_64.h new file mode 100644 index 0000000..605a25f --- /dev/null +++ b/source/kit/mersenne_twister_64.h @@ -0,0 +1,37 @@ +#ifndef KIT_MERSENNE_TWISTER_64_H +#define KIT_MERSENNE_TWISTER_64_H + +#include <stddef.h> +#include <stdint.h> + +#ifdef __cplusplus +extern "C" { +#endif + +enum { + KIT_MT64_N = 312, +}; + +typedef struct { + uint64_t mt[KIT_MT64_N]; + uint64_t index; +} kit_mt64_state_t; + +void kit_mt64_init(kit_mt64_state_t *state, uint64_t seed); + +uint64_t kit_mt64_generate(kit_mt64_state_t *state); + +uint64_t kit_mt64_seed(); + +#ifndef LAPLACE_DISABLE_SHORT_NAMES +# define mt64_state_t kit_mt64_state_t +# define mt64_init kit_mt64_init +# define mt64_generate kit_mt64_generate +# define mt64_seed kit_mt64_seed +#endif + +#ifdef __cplusplus +} +#endif + +#endif diff --git a/source/test/unittests/CMakeLists.txt b/source/test/unittests/CMakeLists.txt index 3dda690..5ec950d 100644 --- a/source/test/unittests/CMakeLists.txt +++ b/source/test/unittests/CMakeLists.txt @@ -3,4 +3,5 @@ target_sources( PRIVATE async_function.test.c test_duration.test.c main.test.c string_ref.test.c atomic.test.c array_ref.test.c input_stream.test.c - lower_bound.test.c input_buffer.test.c dynamic_array.test.c) + lower_bound.test.c mersenne_twister_64.test.c input_buffer.test.c + dynamic_array.test.c) diff --git a/source/test/unittests/mersenne_twister_64.test.c b/source/test/unittests/mersenne_twister_64.test.c new file mode 100644 index 0000000..37c82b4 --- /dev/null +++ b/source/test/unittests/mersenne_twister_64.test.c @@ -0,0 +1,33 @@ +#include "../../kit/mersenne_twister_64.h" + +#define TEST_FILE mersenne_twister_64 +#include "../../kit_test/test.h" + +enum { SIZE = 1000 }; + +TEST("mt64 same seeds") { + uint64_t seed = mt64_seed(); + + mt64_state_t foo, bar; + mt64_init(&foo, seed); + mt64_init(&bar, seed); + + int ok = 1; + for (ptrdiff_t i = 0; i < SIZE; i++) + ok = ok && mt64_generate(&foo) == mt64_generate(&bar); + + REQUIRE(ok); +} + +TEST("mt64 different seeds") { + mt64_state_t foo, bar; + mt64_init(&foo, 42); + mt64_init(&bar, 4242424242); + + ptrdiff_t difference_count = 0; + for (ptrdiff_t i = 0; i < SIZE; i++) + if (mt64_generate(&foo) != mt64_generate(&bar)) + difference_count++; + + REQUIRE(difference_count > SIZE / 2); +} |