Skip to main content

rand_float/
standard.rs

1//! Shift-and-scale conversion.
2//!
3//! The technique used by most language runtimes and libraries: take as many
4//! bits as the target type has of precision (53 for `f64`, 24 for `f32`) and
5//! scale them by the corresponding power of two. Every arithmetic step is
6//! exact, so the result is equispaced: the 2⁵³ multiples of 2⁻⁵³ in [0 . . 1)
7//! for `f64`, each with probability 2⁻⁵³.
8//!
9//! This is by far the cheapest conversion, but most representable floats never
10//! occur: no nonzero value below 2⁻⁵³ can appear.
11//!
12//! Note that using a larger number of bits (e.g., Go used 63) leads to a
13//! nonuniform distribution because of the round-to-even of IEEE 754.
14
15/// 2⁻⁵³.
16const TWO_M53: f64 = 1.0 / (1u64 << 53) as f64;
17/// 2⁻²⁴.
18const TWO_M24: f32 = 1.0 / (1u32 << 24) as f32;
19
20/// Returns a random `f64` uniform on the 2⁵³ multiples of 2⁻⁵³ in [0 . . 1):
21/// the top 53 bits of one word, scaled by 2⁻⁵³.
22#[inline]
23pub fn f64_53bits(mut bits: impl FnMut() -> u64) -> f64 {
24    (bits() >> 11) as f64 * TWO_M53
25}
26
27/// Returns a random `f32` uniform on the 2²⁴ multiples of 2⁻²⁴ in [0 . . 1):
28/// the top 24 bits of one word, scaled by 2⁻²⁴.
29#[inline]
30pub fn f32_24bits(mut bits: impl FnMut() -> u64) -> f32 {
31    (bits() >> 40) as f32 * TWO_M24
32}
33
34#[cfg(test)]
35mod tests {
36    use super::*;
37    use crate::sources::Weyl;
38
39    #[test]
40    fn test_range_and_distribution() {
41        let mut rng = Weyl(42);
42        for _ in 0..100_000 {
43            let x = f64_53bits(|| rng.next_u64());
44            assert!((0.0..1.0).contains(&x), "f64_53bits: {x}");
45            // Every value is a multiple of 2^-53.
46            assert_eq!(x, (x / TWO_M53).round() * TWO_M53);
47            let y = f32_24bits(|| rng.next_u64());
48            assert!((0.0..1.0).contains(&y), "f32_24bits: {y}");
49            assert_eq!(y, (y / TWO_M24).round() * TWO_M24);
50        }
51    }
52
53    #[test]
54    fn test_extremes() {
55        assert_eq!(f64_53bits(|| 0), 0.0);
56        assert_eq!(f64_53bits(|| !0u64), 1.0 - TWO_M53);
57        assert_eq!(f32_24bits(|| 0), 0.0);
58        assert_eq!(f32_24bits(|| !0u64), 1.0 - TWO_M24);
59    }
60}