Implement MT19937 for challenge 21.

src/mt19937.rs

@@ -0,0 +1,80 @@
+use std::num::Wrapping;
+
+#[derive(Debug)]
+pub struct MT19937 {
+    mt: Vec<u32>,
+    index: usize,
+}
+
+const N: usize = 624;
+
+impl MT19937 {
+    pub fn new() -> MT19937 {
+        let mut mt = MT19937 {
+            mt: vec![0; N],
+            index: N + 1,
+        };
+        mt.seed(5489);
+        mt
+    }
+
+    pub fn seed(&mut self, seed: u32) {
+        self.index = N;
+        const F: u32 = 1812433253;
+        self.mt[0] = seed;
+        for i in 1..N {
+            self.mt[i] = (Wrapping(F) * Wrapping(self.mt[i - 1] ^ (self.mt[i - 1] >> 30))
+                + Wrapping(i as u32))
+            .0;
+        }
+    }
+
+    pub fn extract_number(&mut self) -> u32 {
+        if self.index >= N {
+            if self.index > N {
+                panic!("Generator was never seeded");
+            }
+            self.twist();
+        }
+
+        const S: u32 = 7;
+        const T: u32 = 15;
+        const U: u32 = 11;
+        const B: u32 = 0x9D2C5680;
+        const C: u32 = 0xEFC60000;
+        const D: u32 = 0xFFFFFFFF;
+        const L: u32 = 18;
+
+        let mut y = self.mt[self.index];
+        y = y ^ ((y >> U) & D);
+        y = y ^ ((y << S) & B);
+        y = y ^ ((y << T) & C);
+        y = y ^ (y >> L);
+        self.index += 1;
+        y
+    }
+
+    fn twist(&mut self) {
+        const M: usize = 397;
+        const R: u32 = 31;
+        const A: u32 = 0x9908B0DF;
+        const LOWER_MASK: u32 = (1 << R) - 1; // 0x7fffffff
+        const UPPER_MASK: u32 = !LOWER_MASK; // 0x80000000
+        const FIRST_HALF: usize = N - M;
+
+        for i in 0..FIRST_HALF {
+            let bits = self.mt[i] & UPPER_MASK | self.mt[i + 1] & LOWER_MASK;
+            self.mt[i] = self.mt[i + M] ^ (bits >> 1) ^ ((bits & 1) * A);
+        }
+
+        for i in FIRST_HALF..(N - 1) {
+            let bits = self.mt[i] & UPPER_MASK | self.mt[i + 1] & LOWER_MASK;
+            self.mt[i] = self.mt[i - FIRST_HALF] ^ (bits >> 1) ^ ((bits & 1) * A);
+        }
+
+        let i = N - 1;
+        let bits = self.mt[i] & UPPER_MASK | self.mt[0] & LOWER_MASK;
+        self.mt[i] = self.mt[M - 1] ^ (bits >> 1) ^ ((bits & 1) * A);
+        self.index = 0;
+    }
+}