Skip to main content

fcmaes_core/
rng.rs

1//! Random number generation for the optimizers.
2//!
3//! Wraps `rand_pcg::Pcg64` behind the small set of scalar and vector
4//! distributions required by the optimizers.
5//!
6//! # Examples
7//!
8//! Every optimizer takes an [`Rng`], so seeding it is what makes a run
9//! reproducible:
10//!
11//! ```
12//! use fcmaes_core::Rng;
13//!
14//! let mut first = Rng::new(42);
15//! let mut second = Rng::new(42);
16//! assert_eq!(first.uniform01(), second.uniform01());
17//!
18//! // A different seed gives an independent stream.
19//! let mut other = Rng::new(43);
20//! assert_ne!(Rng::new(42).uniform01(), other.uniform01());
21//! ```
22//!
23//! # Reference
24//!
25//! M. E. O'Neill, [“PCG: A Family of Simple Fast Space-Efficient Statistically
26//! Good Algorithms for Random Number
27//! Generation”](https://www.pcg-random.org/paper.html) (2014).
28
29use rand::{Rng as _, SeedableRng};
30use rand_distr::StandardNormal;
31use rand_pcg::Pcg64;
32
33/// Deterministic, seedable RNG shared by all fcmaes optimizers.
34#[derive(Clone, Debug)]
35pub struct Rng {
36    inner: Pcg64,
37}
38
39impl Rng {
40    /// Seed from a single 64-bit value (convenience for `seed`/`runid` pairs).
41    pub fn new(seed: u64) -> Self {
42        Self {
43            inner: Pcg64::seed_from_u64(seed),
44        }
45    }
46
47    /// Seed from a full 128-bit state + stream selector.
48    pub fn from_state_stream(state: u128, stream: u128) -> Self {
49        Self {
50            inner: Pcg64::new(state, stream),
51        }
52    }
53
54    /// Draw a uniform value in `[0, 1)`.
55    #[inline]
56    pub fn uniform01(&mut self) -> f64 {
57        self.inner.r#gen::<f64>()
58    }
59
60    /// Draw a full-width seed for an independent child operation.
61    #[inline]
62    pub fn next_u64(&mut self) -> u64 {
63        self.inner.r#gen::<u64>()
64    }
65
66    /// Standard normal sample `N(0, 1)`.
67    #[inline]
68    pub fn gaussian(&mut self) -> f64 {
69        self.inner.sample(StandardNormal)
70    }
71
72    /// Draw a normal value with mean `mu` and standard deviation `sdev`.
73    #[inline]
74    pub fn normreal(&mut self, mu: f64, sdev: f64) -> f64 {
75        self.gaussian() * sdev + mu
76    }
77
78    /// Draw an integer in `[0, max)` by scaling a uniform value.
79    ///
80    /// Returns zero when `max <= 0`.
81    #[inline]
82    pub fn int_below(&mut self, max: i64) -> i64 {
83        if max <= 0 {
84            return 0;
85        }
86        (max as f64 * self.uniform01()) as i64
87    }
88
89    /// Vector of `dim` uniform `[0, 1)` samples.
90    pub fn uniform_vec(&mut self, dim: usize) -> Vec<f64> {
91        (0..dim).map(|_| self.uniform01()).collect()
92    }
93
94    /// Vector of `dim` standard-normal samples.
95    pub fn normal_vec(&mut self, dim: usize) -> Vec<f64> {
96        (0..dim).map(|_| self.gaussian()).collect()
97    }
98}
99
100#[cfg(test)]
101mod tests {
102    use super::*;
103
104    #[test]
105    fn same_seed_reproduces_stream() {
106        let mut a = Rng::new(12345);
107        let mut b = Rng::new(12345);
108        for _ in 0..100 {
109            assert_eq!(a.uniform01(), b.uniform01());
110        }
111    }
112
113    #[test]
114    fn different_seed_diverges() {
115        let mut a = Rng::new(1);
116        let mut b = Rng::new(2);
117        // Extremely unlikely the first draws coincide.
118        assert_ne!(a.uniform01(), b.uniform01());
119    }
120
121    #[test]
122    fn full_width_seeds_reproduce_and_advance() {
123        let mut a = Rng::new(42);
124        let mut b = Rng::new(42);
125        let first = a.next_u64();
126        assert_eq!(first, b.next_u64());
127        assert_ne!(first, a.next_u64());
128    }
129
130    #[test]
131    fn uniform01_in_range() {
132        let mut rng = Rng::new(7);
133        for _ in 0..10_000 {
134            let u = rng.uniform01();
135            assert!((0.0..1.0).contains(&u));
136        }
137    }
138
139    #[test]
140    fn int_below_in_range() {
141        let mut rng = Rng::new(7);
142        for _ in 0..10_000 {
143            let v = rng.int_below(5);
144            assert!((0..5).contains(&v));
145        }
146        assert_eq!(rng.int_below(0), 0);
147        assert_eq!(rng.int_below(-3), 0);
148    }
149
150    #[test]
151    fn gaussian_statistics_are_sane() {
152        let mut rng = Rng::new(2024);
153        let n = 200_000;
154        let mut sum = 0.0;
155        let mut sumsq = 0.0;
156        for _ in 0..n {
157            let g = rng.gaussian();
158            sum += g;
159            sumsq += g * g;
160        }
161        let mean = sum / n as f64;
162        let var = sumsq / n as f64 - mean * mean;
163        assert!(mean.abs() < 0.02, "mean={mean}");
164        assert!((var - 1.0).abs() < 0.03, "var={var}");
165    }
166
167    #[test]
168    fn state_stream_normal_and_vector_helpers() {
169        let mut first = Rng::from_state_stream(123, 7);
170        let mut second = Rng::from_state_stream(123, 7);
171        assert_eq!(first.uniform_vec(8), second.uniform_vec(8));
172        let normal = first.normal_vec(16);
173        assert_eq!(normal.len(), 16);
174        assert!(normal.iter().all(|value| value.is_finite()));
175        let sample = first.normreal(10.0, 0.0);
176        assert_eq!(sample, 10.0);
177        assert!(first.uniform_vec(0).is_empty());
178        assert!(first.normal_vec(0).is_empty());
179    }
180}