diff options
Diffstat (limited to 'tests/rand.c')
-rw-r--r-- | tests/rand.c | 83 |
1 files changed, 83 insertions, 0 deletions
diff --git a/tests/rand.c b/tests/rand.c new file mode 100644 index 0000000..925d1ba --- /dev/null +++ b/tests/rand.c @@ -0,0 +1,83 @@ +// Taken from +// http://www.math.sci.hiroshima-u.ac.jp/m-mat/MT/MT2002/CODES/mt19937ar.c + +#define _XOPEN_SOURCE 700 + +#include <assert.h> +#include <inttypes.h> +#include <stdbool.h> +#include <stdio.h> + + +#define N 624U +static const uint16_t M = 397; +static const uint32_t MATRIX_A = 0x9908b0dfUL; +static const uint32_t UPPER_MASK = 0x80000000UL; +static const uint32_t LOWER_MASK = 0x7fffffffUL; + +static uint32_t mt[N]; +static uint16_t mti; +static bool initialized = false; + +static void +init(uint32_t s) { + initialized = true; + mt[0] = s & 0xffffffffUL; + for (mti = 1; mti < N; mti++) { + mt[mti] = (1812433253 * (mt[mti - 1] ^ (mt[mti - 1] >> 30))) + mti; + } +} + +/// Generates a random number on the [ 0, 0xffffffff ]-interval. +static uint32_t +generate(void) { + assert(initialized == true); + static uint32_t mag01[2] = { 0, MATRIX_A }; + + if (mti >= N) { + uint16_t kk; + + for (kk = 0; kk < N - M; kk++) { + const uint32_t y0 = + (mt[kk] & UPPER_MASK) | + (mt[kk + 1] & LOWER_MASK); + mt[kk] = mt[kk + M] ^ (y0 >> 1) ^ mag01[y0 & 1U]; + } + + for (; kk < N - 1; kk++) { + const uint32_t y1 = + (mt[kk] & UPPER_MASK) | + (mt[kk + 1] & LOWER_MASK); + mt[kk] = mt[kk + (M - N)] ^ (y1 >> 1) ^ mag01[y1 & 1U]; + } + + const uint32_t y2 = + (mt[N - 1] & UPPER_MASK) | + (mt[0] & LOWER_MASK); + mt[N - 1] = mt[M - 1] ^ (y2 >> 1) ^ mag01[y2 & 1U]; + mti = 0; + } + + uint32_t y = mt[mti++]; + + y ^= (y >> 11); + y ^= (y << 7) & 0x9d2c5680UL; + y ^= (y << 15) & 0xefc60000UL; + y ^= (y >> 18); + + return y; +} + + +int +main(void) { + init(2222); + + for (uint16_t i = 0; i < 10000; i++) { + assert( + printf("%" PRIu32 "\n", generate()) > 0 + ); + } + + return 0; +} |