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}