rand_pcg/
pcg128.rs

1// Copyright 2018 Developers of the Rand project.
2// Copyright 2017 Paul Dicker.
3// Copyright 2014-2017 Melissa O'Neill and PCG Project contributors
4//
5// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
6// https://www.apache.org/licenses/LICENSE-2.0> or the MIT license
7// <LICENSE-MIT or https://opensource.org/licenses/MIT>, at your
8// option. This file may not be copied, modified, or distributed
9// except according to those terms.
10
11//! PCG random number generators
12
13// This is the default multiplier used by PCG for 128-bit state.
14const MULTIPLIER: u128 = 0x2360_ED05_1FC6_5DA4_4385_DF64_9FCC_F645;
15
16use core::fmt;
17use rand_core::{impls, le, RngCore, SeedableRng};
18#[cfg(feature = "serde")]
19use serde::{Deserialize, Serialize};
20
21/// A PCG random number generator (XSL RR 128/64 (LCG) variant).
22///
23/// Permuted Congruential Generator with 128-bit state, internal Linear
24/// Congruential Generator, and 64-bit output via "xorshift low (bits),
25/// random rotation" output function.
26///
27/// This is a 128-bit LCG with explicitly chosen stream with the PCG-XSL-RR
28/// output function. This combination is the standard `pcg64`.
29///
30/// Despite the name, this implementation uses 32 bytes (256 bit) space
31/// comprising 128 bits of state and 128 bits stream selector. These are both
32/// set by `SeedableRng`, using a 256-bit seed.
33///
34/// Note that two generators with different stream parameters may be closely
35/// correlated.
36#[derive(Clone, PartialEq, Eq)]
37#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
38pub struct Lcg128Xsl64 {
39    state: u128,
40    increment: u128,
41}
42
43/// [`Lcg128Xsl64`] is also officially known as `pcg64`.
44pub type Pcg64 = Lcg128Xsl64;
45
46impl Lcg128Xsl64 {
47    /// Multi-step advance functions (jump-ahead, jump-back)
48    ///
49    /// The method used here is based on Brown, "Random Number Generation
50    /// with Arbitrary Stride,", Transactions of the American Nuclear
51    /// Society (Nov. 1994).  The algorithm is very similar to fast
52    /// exponentiation.
53    ///
54    /// Even though delta is an unsigned integer, we can pass a
55    /// signed integer to go backwards, it just goes "the long way round".
56    ///
57    /// Using this function is equivalent to calling `next_64()` `delta`
58    /// number of times.
59    #[inline]
60    pub fn advance(&mut self, delta: u128) {
61        let mut acc_mult: u128 = 1;
62        let mut acc_plus: u128 = 0;
63        let mut cur_mult = MULTIPLIER;
64        let mut cur_plus = self.increment;
65        let mut mdelta = delta;
66
67        while mdelta > 0 {
68            if (mdelta & 1) != 0 {
69                acc_mult = acc_mult.wrapping_mul(cur_mult);
70                acc_plus = acc_plus.wrapping_mul(cur_mult).wrapping_add(cur_plus);
71            }
72            cur_plus = cur_mult.wrapping_add(1).wrapping_mul(cur_plus);
73            cur_mult = cur_mult.wrapping_mul(cur_mult);
74            mdelta /= 2;
75        }
76        self.state = acc_mult.wrapping_mul(self.state).wrapping_add(acc_plus);
77    }
78
79    /// Construct an instance compatible with PCG seed and stream.
80    ///
81    /// Note that the highest bit of the `stream` parameter is discarded
82    /// to simplify upholding internal invariants.
83    ///
84    /// Note that two generators with different stream parameters may be closely
85    /// correlated.
86    ///
87    /// PCG specifies the following default values for both parameters:
88    ///
89    /// - `state = 0xcafef00dd15ea5e5`
90    /// - `stream = 0xa02bdbf7bb3c0a7ac28fa16a64abf96`
91    pub fn new(state: u128, stream: u128) -> Self {
92        // The increment must be odd, hence we discard one bit:
93        let increment = (stream << 1) | 1;
94        Lcg128Xsl64::from_state_incr(state, increment)
95    }
96
97    #[inline]
98    fn from_state_incr(state: u128, increment: u128) -> Self {
99        let mut pcg = Lcg128Xsl64 { state, increment };
100        // Move away from initial value:
101        pcg.state = pcg.state.wrapping_add(pcg.increment);
102        pcg.step();
103        pcg
104    }
105
106    #[inline]
107    fn step(&mut self) {
108        // prepare the LCG for the next round
109        self.state = self
110            .state
111            .wrapping_mul(MULTIPLIER)
112            .wrapping_add(self.increment);
113    }
114}
115
116// Custom Debug implementation that does not expose the internal state
117impl fmt::Debug for Lcg128Xsl64 {
118    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
119        write!(f, "Lcg128Xsl64 {{}}")
120    }
121}
122
123impl SeedableRng for Lcg128Xsl64 {
124    type Seed = [u8; 32];
125
126    /// We use a single 255-bit seed to initialise the state and select a stream.
127    /// One `seed` bit (lowest bit of `seed[8]`) is ignored.
128    fn from_seed(seed: Self::Seed) -> Self {
129        let mut seed_u64 = [0u64; 4];
130        le::read_u64_into(&seed, &mut seed_u64);
131        let state = u128::from(seed_u64[0]) | (u128::from(seed_u64[1]) << 64);
132        let incr = u128::from(seed_u64[2]) | (u128::from(seed_u64[3]) << 64);
133
134        // The increment must be odd, hence we discard one bit:
135        Lcg128Xsl64::from_state_incr(state, incr | 1)
136    }
137}
138
139impl RngCore for Lcg128Xsl64 {
140    #[inline]
141    fn next_u32(&mut self) -> u32 {
142        self.next_u64() as u32
143    }
144
145    #[inline]
146    fn next_u64(&mut self) -> u64 {
147        self.step();
148        output_xsl_rr(self.state)
149    }
150
151    #[inline]
152    fn fill_bytes(&mut self, dest: &mut [u8]) {
153        impls::fill_bytes_via_next(self, dest)
154    }
155}
156
157/// A PCG random number generator (XSL 128/64 (MCG) variant).
158///
159/// Permuted Congruential Generator with 128-bit state, internal Multiplicative
160/// Congruential Generator, and 64-bit output via "xorshift low (bits),
161/// random rotation" output function.
162///
163/// This is a 128-bit MCG with the PCG-XSL-RR output function, also known as
164/// `pcg64_fast`.
165/// Note that compared to the standard `pcg64` (128-bit LCG with PCG-XSL-RR
166/// output function), this RNG is faster, also has a long cycle, and still has
167/// good performance on statistical tests.
168#[derive(Clone, PartialEq, Eq)]
169#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
170pub struct Mcg128Xsl64 {
171    state: u128,
172}
173
174/// A friendly name for [`Mcg128Xsl64`] (also known as `pcg64_fast`).
175pub type Pcg64Mcg = Mcg128Xsl64;
176
177impl Mcg128Xsl64 {
178    /// Multi-step advance functions (jump-ahead, jump-back)
179    ///
180    /// The method used here is based on Brown, "Random Number Generation
181    /// with Arbitrary Stride,", Transactions of the American Nuclear
182    /// Society (Nov. 1994).  The algorithm is very similar to fast
183    /// exponentiation.
184    ///
185    /// Even though delta is an unsigned integer, we can pass a
186    /// signed integer to go backwards, it just goes "the long way round".
187    ///
188    /// Using this function is equivalent to calling `next_64()` `delta`
189    /// number of times.
190    #[inline]
191    pub fn advance(&mut self, delta: u128) {
192        let mut acc_mult: u128 = 1;
193        let mut acc_plus: u128 = 0;
194        let mut cur_mult = MULTIPLIER;
195        let mut cur_plus: u128 = 0;
196        let mut mdelta = delta;
197
198        while mdelta > 0 {
199            if (mdelta & 1) != 0 {
200                acc_mult = acc_mult.wrapping_mul(cur_mult);
201                acc_plus = acc_plus.wrapping_mul(cur_mult).wrapping_add(cur_plus);
202            }
203            cur_plus = cur_mult.wrapping_add(1).wrapping_mul(cur_plus);
204            cur_mult = cur_mult.wrapping_mul(cur_mult);
205            mdelta /= 2;
206        }
207        self.state = acc_mult.wrapping_mul(self.state).wrapping_add(acc_plus);
208    }
209
210    /// Construct an instance compatible with PCG seed.
211    ///
212    /// Note that PCG specifies a default value for the parameter:
213    ///
214    /// - `state = 0xcafef00dd15ea5e5`
215    pub fn new(state: u128) -> Self {
216        // Force low bit to 1, as in C version (C++ uses `state | 3` instead).
217        Mcg128Xsl64 { state: state | 1 }
218    }
219}
220
221// Custom Debug implementation that does not expose the internal state
222impl fmt::Debug for Mcg128Xsl64 {
223    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
224        write!(f, "Mcg128Xsl64 {{}}")
225    }
226}
227
228/// We use a single 126-bit seed to initialise the state and select a stream.
229/// Two `seed` bits (lowest order of last byte) are ignored.
230impl SeedableRng for Mcg128Xsl64 {
231    type Seed = [u8; 16];
232
233    fn from_seed(seed: Self::Seed) -> Self {
234        // Read as if a little-endian u128 value:
235        let mut seed_u64 = [0u64; 2];
236        le::read_u64_into(&seed, &mut seed_u64);
237        let state = u128::from(seed_u64[0]) | (u128::from(seed_u64[1]) << 64);
238        Mcg128Xsl64::new(state)
239    }
240}
241
242impl RngCore for Mcg128Xsl64 {
243    #[inline]
244    fn next_u32(&mut self) -> u32 {
245        self.next_u64() as u32
246    }
247
248    #[inline]
249    fn next_u64(&mut self) -> u64 {
250        self.state = self.state.wrapping_mul(MULTIPLIER);
251        output_xsl_rr(self.state)
252    }
253
254    #[inline]
255    fn fill_bytes(&mut self, dest: &mut [u8]) {
256        impls::fill_bytes_via_next(self, dest)
257    }
258}
259
260#[inline(always)]
261fn output_xsl_rr(state: u128) -> u64 {
262    // Output function XSL RR ("xorshift low (bits), random rotation")
263    // Constants are for 128-bit state, 64-bit output
264    const XSHIFT: u32 = 64; // (128 - 64 + 64) / 2
265    const ROTATE: u32 = 122; // 128 - 6
266
267    let rot = (state >> ROTATE) as u32;
268    let xsl = ((state >> XSHIFT) as u64) ^ (state as u64);
269    xsl.rotate_right(rot)
270}
rand_pcg/pcg128.rs

rand_pcg/
pcg128.rs
