chj_unix_util/
xorshift.rs

1//! Avoid dependency on `rand` crate
2
3// struct Xorshift64(pub u64);
4
5// impl Xorshift64 {
6//     pub fn new(seed: u64) -> Self {
7//         if seed == 0 {
8//             panic!("seed must not be 0")
9//         }
10//         Self(seed)
11//     }
12
13//     pub fn get(&mut self) -> u64 {
14//         self.0 ^= self.0 << 13;
15//         self.0 ^= self.0 >> 17;
16//         self.0 ^= self.0 << 5;
17//         self.0
18//     }
19// }
20
21// https://handwiki.org/wiki/Xorshift
22
23pub struct Splitmix64(pub u64);
24
25impl Splitmix64 {
26    pub fn get(&mut self) -> u64 {
27        self.0 = self.0.wrapping_add(0x9E3779B97F4A7C15);
28        let mut result = self.0;
29        result = (result ^ (result >> 30)).wrapping_mul(0xBF58476D1CE4E5B9);
30        result = (result ^ (result >> 27)).wrapping_mul(0x94D049BB133111EB);
31        result ^ (result >> 31)
32    }
33}
34
35pub struct Xorshift128plus(pub [u64; 2]);
36
37impl Xorshift128plus {
38    pub fn new(seed: u64) -> Self {
39        let mut s = Splitmix64(seed);
40        Self([s.get(), s.get()])
41    }
42
43    pub fn get(&mut self) -> u64 {
44        let mut t = self.0[0];
45        let s = self.0[1];
46        self.0[0] = s;
47        t ^= t << 23;
48        t ^= t >> 18;
49        t ^= s ^ (s >> 5);
50        self.0[1] = t;
51        t.wrapping_add(s)
52    }
53}
54
55#[test]
56fn t_splitmix64() {
57    let mut s = Splitmix64(0);
58    assert_eq!(s.get(), 16294208416658607535);
59    assert_eq!(s.get(), 7960286522194355700);
60    assert_eq!(s.get(), 487617019471545679);
61    assert_eq!(s.get(), 17909611376780542444);
62    for _ in 0..1000000 - 4 {
63        s.get();
64    }
65    assert_eq!(s.get(), 14850574393604363050);
66    assert_eq!(s.get(), 1562119273537874705);
67    assert_eq!(s.get(), 1986060996022186059);
68    assert_eq!(s.get(), 17599710064536246394);
69    assert_eq!(s.get(), 9868005493142914005);
70}
71
72#[test]
73fn t_xorshift128plus() {
74    let mut s = Xorshift128plus::new(0);
75    assert_eq!(s.get(), 148304652509113927);
76    assert_eq!(s.get(), 6897519897668720478);
77    assert_eq!(s.get(), 8466708535677759538);
78    assert_eq!(s.get(), 4573841993332567017);
79    for _ in 0..1000000 - 4 {
80        s.get();
81    }
82    assert_eq!(s.get(), 7831888472505485542);
83    assert_eq!(s.get(), 9546132825291831751);
84    assert_eq!(s.get(), 10667650326493356499);
85    assert_eq!(s.get(), 13170528539071654594);
86    assert_eq!(s.get(), 16110659944319132175);
87}
88
89/* Yields the same results as:
90
91#include <cassert>
92
93// https://handwiki.org/wiki/Xorshift
94
95#include <stdint.h>
96
97struct splitmix64_state {
98    uint64_t s;
99};
100
101uint64_t splitmix64(struct splitmix64_state *state) {
102    uint64_t result = (state->s += 0x9E3779B97F4A7C15);
103    result = (result ^ (result >> 30)) * 0xBF58476D1CE4E5B9;
104    result = (result ^ (result >> 27)) * 0x94D049BB133111EB;
105    return result ^ (result >> 31);
106}
107
108
109struct xorshift128p_state {
110    uint64_t x[2];
111};
112
113/* The state must be seeded so that it is not all zero */
114uint64_t xorshift128p(struct xorshift128p_state *state)
115{
116    uint64_t t = state->x[0];
117    uint64_t const s = state->x[1];
118    state->x[0] = s;
119    t ^= t << 23;       // a
120    t ^= t >> 18;       // b -- Again, the shifts and the multipliers are tunable
121    t ^= s ^ (s >> 5);  // c
122    state->x[1] = t;
123    return t + s;
124}
125
126int main() {
127
128    struct splitmix64_state smstate = {0};
129
130    uint64_t tmp1 = splitmix64(&smstate);
131    uint64_t tmp2 = splitmix64(&smstate);
132
133    assert(tmp1 == 16294208416658607535ull);
134    assert(tmp2 == 7960286522194355700ull);
135
136    struct xorshift128p_state pstate = {{tmp1, tmp2}};
137
138    uint64_t a;
139    a = xorshift128p(&pstate); assert(a == 148304652509113927ull);
140    a = xorshift128p(&pstate); assert(a == 6897519897668720478ull);
141    a = xorshift128p(&pstate); assert(a == 8466708535677759538ull);
142    a = xorshift128p(&pstate); assert(a == 4573841993332567017ull);
143
144    for (int i = 0; i < 1000000 - 4; i++) {
145        xorshift128p(&pstate);
146    }
147
148    a = xorshift128p(&pstate); assert(a == 7831888472505485542ull);
149    a = xorshift128p(&pstate); assert(a == 9546132825291831751ull);
150    a = xorshift128p(&pstate); assert(a == 10667650326493356499ull);
151    a = xorshift128p(&pstate); assert(a == 13170528539071654594ull);
152    a = xorshift128p(&pstate); assert(a == 16110659944319132175ull);
153
154    return 0;
155}
156 */