commit 4dda41b3c51551872d4c650f6fcf55859c1abe9c
parent bf02706fd5fb6adea3f1c66d7fc798888d9d8cbb
Author: Jack Sangdahl <jack@pngu.org>
Date: Thu, 23 Oct 2025 11:04:38 -0600
challenge 21: implement the mt19937 mersenne twister rng
Diffstat:
3 files changed, 115 insertions(+), 0 deletions(-)
diff --git a/lib/mt19937.c b/lib/mt19937.c
@@ -0,0 +1,69 @@
+#include "mt19937.h"
+
+struct mt19937_state mt19937_init(uint32_t seed) {
+ struct mt19937_state obj;
+ memset(&obj, 0, sizeof(struct mt19937_state));
+ obj.mt[0] = seed;
+ for (size_t i = 1; i < MT19937_N; i++)
+ obj.mt[i] = 0x6C078965 * (obj.mt[i - 1] ^ (obj.mt[i - 1] >> 30)) + i;
+ return obj;
+}
+
+void mt19937_regenerate(struct mt19937_state *obj) {
+ for (size_t i = 0; i < MT19937_N; i++) {
+ uint32_t y = obj->mt[i] & 0x80000000;
+ y += obj->mt[(i + 1) % MT19937_N] & 0x7FFFFFFF;
+ uint32_t z = obj->mt[(i + 397) % MT19937_N];
+ obj->mt[i] = z ^ (y >> 1);
+ if (y % 2)
+ obj->mt[i] ^= 0x9908B0DF;
+ }
+}
+
+uint32_t mt19937_temper(uint32_t y) {
+ y ^= y >> 11;
+ y ^= (y << 7) & MT19937_B;
+ y ^= (y << 15) & MT19937_C;
+ y ^= y >> 18;
+ return y;
+}
+
+uint32_t mt19937_extract(struct mt19937_state *obj) {
+ if (!obj->idx)
+ mt19937_regenerate(obj);
+ uint32_t y = obj->mt[obj->idx];
+ obj->idx++;
+ obj->idx %= MT19937_N;
+ return mt19937_temper(y);
+}
+
+uint32_t unshift_right(uint32_t val, unsigned shift) {
+ uint32_t z = 0;
+ for (size_t i = 0; i * shift < 32; i++) {
+ uint32_t mask = (0xFFFFFFFF << (32 - shift)) >> (shift * i);
+ uint32_t part = val & mask;
+ val ^= part >> shift;
+ z |= part;
+ }
+ return z;
+}
+
+uint32_t unshift_left_mask(uint32_t val, unsigned shift, uint32_t mask) {
+ uint32_t z = 0;
+ for (size_t i = 0; i * shift < 32; i++) {
+ uint32_t part_mask = (0xFFFFFFFF >> (32 - shift)) << (shift * i);
+ uint32_t part = val & part_mask;
+ val ^= (part << shift) & mask;
+ z |= part;
+ }
+ return z;
+}
+
+uint32_t mt19937_untemper(uint32_t y) {
+ uint32_t res = y;
+ res = unshift_right(res, 18);
+ res = unshift_left_mask(res, 15, MT19937_C);
+ res = unshift_left_mask(res, 7, MT19937_B);
+ res = unshift_right(res, 11);
+ return res;
+}
diff --git a/lib/mt19937.h b/lib/mt19937.h
@@ -0,0 +1,21 @@
+#pragma once
+
+#include <stdint.h>
+#include <string.h>
+
+#define MT19937_N 624
+#define MT19937_B 0x9D2C5680
+#define MT19937_C 0xEFC60000
+
+struct mt19937_state {
+ uint32_t mt[MT19937_N];
+ uint32_t idx;
+};
+
+struct mt19937_state mt19937_init(uint32_t seed);
+void mt19937_regenerate(struct mt19937_state *obj);
+uint32_t mt19937_temper(uint32_t y);
+uint32_t mt19937_extract(struct mt19937_state *obj);
+uint32_t unshift_right(uint32_t val, unsigned shift);
+uint32_t unshift_left_mask(uint32_t val, unsigned shift, uint32_t mask);
+uint32_t mt19937_untemper(uint32_t y);
diff --git a/src/21_implement_mt19937_rng.c b/src/21_implement_mt19937_rng.c
@@ -0,0 +1,25 @@
+// Implement the MT19937 Mersenne Twister RNG
+#include "mt19937.h"
+#include <assert.h>
+#include <inttypes.h>
+#include <stdio.h>
+
+/*
+ * you can get the pseudocode for this from wikipedia. if you're writing in
+ * python, ruby, or (gah) php, your language is probably already giving you
+ * mt19937 as 'rand()'. don't use it. write the rng yourself.
+ */
+
+int main(void) {
+ static uint32_t seed = 5489;
+ static uint32_t known_vals[] = {
+ -795755684, 581869302, -404620562, -708632711, 545404204,
+ -133711905, -372047867, 949333985, -1579004998, 1323567403};
+ struct mt19937_state mt = mt19937_init(seed);
+ for (int i = 0; i < 10; i++) {
+ uint32_t extracted = mt19937_extract(&mt);
+ printf("%" PRIu32 "\n", extracted);
+ assert(extracted == known_vals[i]);
+ }
+ return 0;
+}