Skip to main content

base/
rng.rs

1// rng -- deterministic pseudo-random numbers for simulations, tests and games (NOT for keys:
2// those come from plat::fill_random)
3//
4//   SplitMix64::new(seed)      next_u64; derives seeds (Vigna's splitmix64)
5//   Pcg32::new(seed, stream)   PCG-XSH-RR 64/32 (O'Neill), the reference pcg32_srandom_r: the
6//                              same seed and stream give the same numbers on every target;
7//                              different streams are independent sequences
8//   Pcg32::derive(seed, id)    a generator for one of many objects (seed and id mixed through
9//                              SplitMix64), so adding an object does not shift the others
10//   next_u32 / next_u64 / below(n) (unbiased) / range(lo, hi) / f64 (0..1, 53 bits) /
11//   chance(num, den) / chance_ppm(ppm) / fill(bytes) / shuffle(slice) / advance(delta)
12
13/// SplitMix64: a 64-bit generator mostly used to expand one seed into many
14#[derive(Clone, Debug, PartialEq, Eq)]
15pub struct SplitMix64 {
16    state: u64,
17}
18
19impl SplitMix64 {
20    pub const fn new(seed: u64) -> SplitMix64 {
21        SplitMix64 { state: seed }
22    }
23    pub fn next_u64(&mut self) -> u64 {
24        self.state = self.state.wrapping_add(0x9E37_79B9_7F4A_7C15);
25        let mut z = self.state;
26        z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9);
27        z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB);
28        z ^ (z >> 31)
29    }
30}
31
32const PCG_MUL: u64 = 6_364_136_223_846_793_005;
33
34/// PCG32 (XSH-RR, 64-bit state, 32-bit output)
35#[derive(Clone, Debug, PartialEq, Eq)]
36pub struct Pcg32 {
37    state: u64,
38    inc: u64,
39}
40
41impl Pcg32 {
42    /// the generator of `seed` on `stream` (pcg32_srandom_r(seed, stream))
43    pub fn new(seed: u64, stream: u64) -> Pcg32 {
44        let mut g = Pcg32 { state: 0, inc: (stream << 1) | 1 };
45        g.step();
46        g.state = g.state.wrapping_add(seed);
47        g.step();
48        g
49    }
50    /// a generator for object `id` under a common seed
51    pub fn derive(seed: u64, id: u64) -> Pcg32 {
52        let mut m = SplitMix64::new(seed ^ id.wrapping_mul(0xD6E8_FEB8_6659_FD93));
53        let s = m.next_u64();
54        let st = m.next_u64();
55        Pcg32::new(s, st)
56    }
57    fn step(&mut self) {
58        self.state = self.state.wrapping_mul(PCG_MUL).wrapping_add(self.inc);
59    }
60    pub fn next_u32(&mut self) -> u32 {
61        let old = self.state;
62        self.step();
63        let xs = (((old >> 18) ^ old) >> 27) as u32;
64        let rot = (old >> 59) as u32;
65        xs.rotate_right(rot)
66    }
67    pub fn next_u64(&mut self) -> u64 {
68        let hi = self.next_u32() as u64;
69        (hi << 32) | self.next_u32() as u64
70    }
71    /// uniform in 0..n without modulo bias (n = 0 gives 0)
72    pub fn below(&mut self, n: u32) -> u32 {
73        if n == 0 {
74            return 0;
75        }
76        let threshold = n.wrapping_neg() % n;
77        loop {
78            let r = self.next_u32();
79            if r >= threshold {
80                return r % n;
81            }
82        }
83    }
84    /// uniform in 0..n for a 64-bit n without modulo bias (n = 0 gives 0)
85    fn below64(&mut self, n: u64) -> u64 {
86        if n == 0 {
87            return 0;
88        }
89        let threshold = n.wrapping_neg() % n;
90        loop {
91            let r = self.next_u64();
92            if r >= threshold {
93                return r % n;
94            }
95        }
96    }
97    /// uniform in lo..=hi
98    pub fn range(&mut self, lo: i64, hi: i64) -> i64 {
99        if hi <= lo {
100            return lo;
101        }
102        let span = hi as i128 - lo as i128 + 1; // 2 ..= 2^64
103        if span <= u32::MAX as i128 {
104            return lo + self.below(span as u32) as i64;
105        }
106        if span == 1i128 << 64 {
107            return self.next_u64() as i64;
108        }
109        // wide ranges: rejection on 64 bits, as below()
110        let span = span as u64;
111        let threshold = span.wrapping_neg() % span;
112        loop {
113            let r = self.next_u64();
114            if r >= threshold {
115                return (lo as i128 + (r % span) as i128) as i64;
116            }
117        }
118    }
119    /// uniform in [0, 1) with 53 random bits
120    pub fn f64(&mut self) -> f64 {
121        (self.next_u64() >> 11) as f64 * (1.0 / 9_007_199_254_740_992.0)
122    }
123    /// true with probability num / den
124    pub fn chance(&mut self, num: u32, den: u32) -> bool {
125        num >= den || (den > 0 && self.below(den) < num)
126    }
127    /// true with probability ppm / 1 000 000 (loss and error rates)
128    pub fn chance_ppm(&mut self, ppm: u32) -> bool {
129        self.chance(ppm, 1_000_000)
130    }
131    pub fn fill(&mut self, out: &mut [u8]) {
132        for c in out.chunks_mut(4) {
133            let r = self.next_u32().to_le_bytes();
134            c.copy_from_slice(&r[..c.len()]);
135        }
136    }
137    /// Fisher-Yates
138    pub fn shuffle<T>(&mut self, v: &mut [T]) {
139        let mut i = v.len();
140        while i > 1 {
141            let j = if i <= u32::MAX as usize { self.below(i as u32) as usize } else { self.below64(i as u64) as usize };
142            i -= 1;
143            v.swap(i, j);
144        }
145    }
146    /// skips `delta` outputs in O(log delta) (Brown's jump ahead for LCGs)
147    pub fn advance(&mut self, mut delta: u64) {
148        let (mut acc_mul, mut acc_add) = (1u64, 0u64);
149        let (mut cur_mul, mut cur_add) = (PCG_MUL, self.inc);
150        while delta > 0 {
151            if delta & 1 != 0 {
152                acc_mul = acc_mul.wrapping_mul(cur_mul);
153                acc_add = acc_add.wrapping_mul(cur_mul).wrapping_add(cur_add);
154            }
155            cur_add = cur_mul.wrapping_add(1).wrapping_mul(cur_add);
156            cur_mul = cur_mul.wrapping_mul(cur_mul);
157            delta >>= 1;
158        }
159        self.state = acc_mul.wrapping_mul(self.state).wrapping_add(acc_add);
160    }
161}