From 16c465a30767404812d098c9c0d0ee4579c6bc71 Mon Sep 17 00:00:00 2001 From: Reiner Pope Date: Sat, 10 May 2025 21:54:04 -0700 Subject: [PATCH 1/8] Test f32/f64 generation. --- tests/smoke.rs | 33 +++++++++++++++++++++++++++++++++ 1 file changed, 33 insertions(+) diff --git a/tests/smoke.rs b/tests/smoke.rs index 7c92ee5..3dd1c1f 100644 --- a/tests/smoke.rs +++ b/tests/smoke.rs @@ -75,6 +75,39 @@ fn u128() { } } +#[test] +fn f32() { + let mut r = fastrand::Rng::with_seed(0); + let tiny = (-24.0f32).exp2(); + let mut count_tiny_nonzero = 0; + let mut count_top_half = 0; + for _ in 0..100_000_000 { + let x = r.f32(); + assert!(x >= 0.0 && x < 1.0); + if x > 0.0 && x < tiny { + count_tiny_nonzero += 1; + } else if x > 0.5 { + count_top_half += 1; + } + } + assert!(count_top_half >= 49_000_000); + assert!(count_tiny_nonzero > 0); +} + +#[test] +fn f64() { + let mut r = fastrand::Rng::with_seed(0); + let mut count_top_half = 0; + for _ in 0..100_000_000 { + let x = r.f64(); + assert!(x >= 0.0 && x < 1.0); + if x > 0.5 { + count_top_half += 1; + } + } + assert!(count_top_half >= 49_000_000); +} + #[test] fn fill() { let mut r = fastrand::Rng::new(); From 606d2c054a1cf3d2145be5d1c39584b8227677df Mon Sep 17 00:00:00 2001 From: Reiner Pope Date: Sat, 10 May 2025 21:54:04 -0700 Subject: [PATCH 2/8] Improve quality of f32/f64 generation. The previous int-to-float conversion had a bias of probability 2^-24 / 2^-53 for types f32 / f64 respectively. The new conversion has a bias of 2^-64, which is the same bias as the underlying WyRand generator. The new conversion is a slightly shorter instruction sequence on x86 and ARM, but executes as 1 more uop on x86. Seems unlikely to harm performance much, if at all. https://rust.godbolt.org/z/q3zMxEc3T --- src/lib.rs | 46 ++++++++++++++++++++++++++++++++++++++++------ tests/smoke.rs | 37 +++++++++++++++++++++++++++++++++++++ 2 files changed, 77 insertions(+), 6 deletions(-) diff --git a/src/lib.rs b/src/lib.rs index 51ec997..17d3eec 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -364,18 +364,52 @@ impl Rng { } } + /// Generates a random `f32` in range `0..=1`. + pub fn f32_inclusive(&mut self) -> f32 { + // Generate a number in 0..2^63 then convert to f32 and multiply by 2^(-63). + // + // Even though we're returning f32, we still generate u64 internally to make + // it possible to return nonzero numbers as small as 2^(-63). If we only + // generated u32 internally, the smallest nonzero number we could return + // would be 2^(-32). + // + // The integer we generate is in 0..2^63 rather than 0..2^64 to improve speed + // on x86-64, which has efficient i64->float conversion (cvtsi2ss) but for + // which u64->float conversion must be implemented in software. + // + // There is still some remaining bias in the int-to-float conversion, because + // nonzero numbers <=2^(-64) are never generated, even though they are + // expressible in f32. However, at this point the bias in int-to-float conversion + // is no larger than the bias in the underlying WyRand generator: since it only + // has a 64-bit state, it necessarily already have biases of at least 2^(-64) + // probability. + (self.u64(..) >> 1) as f32 * (-63.0f32).exp2() + } + /// Generates a random `f32` in range `0..1`. pub fn f32(&mut self) -> f32 { - let b = 32; - let f = core::f32::MANTISSA_DIGITS - 1; - f32::from_bits((1 << (b - 2)) - (1 << f) + (self.u32(..) >> (b - f))) - 1.0 + loop { + let x = self.f32_inclusive(); + if x < 1.0 { + return x; + } + } + } + + /// Generates a random `f64` in range `0..=1`. + pub fn f64_inclusive(&mut self) -> f64 { + // See the comment in f32_inclusive() for more details. + (self.u64(..) >> 1) as f64 * (-63.0f64).exp2() } /// Generates a random `f64` in range `0..1`. pub fn f64(&mut self) -> f64 { - let b = 64; - let f = core::f64::MANTISSA_DIGITS - 1; - f64::from_bits((1 << (b - 2)) - (1 << f) + (self.u64(..) >> (b - f))) - 1.0 + loop { + let x = self.f64_inclusive(); + if x < 1.0 { + return x; + } + } } /// Collects `amount` values at random from the iterable into a vector. diff --git a/tests/smoke.rs b/tests/smoke.rs index 3dd1c1f..43475fd 100644 --- a/tests/smoke.rs +++ b/tests/smoke.rs @@ -94,6 +94,29 @@ fn f32() { assert!(count_tiny_nonzero > 0); } +#[test] +fn f32_inclusive() { + let mut r = fastrand::Rng::with_seed(0); + let tiny = (-24.0f32).exp2(); + let mut count_top_half = 0; + let mut count_tiny_nonzero = 0; + let mut count_one = 0; + for _ in 0..100_000_000 { + let x = r.f32_inclusive(); + assert!(x >= 0.0 && x <= 1.0); + if x == 1.0 { + count_one += 1; + } else if x > 0.5 { + count_top_half += 1; + } else if x > 0.0 && x < tiny { + count_tiny_nonzero += 1; + } + } + assert!(count_top_half >= 49_000_000); + assert!(count_one > 0); + assert!(count_tiny_nonzero > 0); +} + #[test] fn f64() { let mut r = fastrand::Rng::with_seed(0); @@ -108,6 +131,20 @@ fn f64() { assert!(count_top_half >= 49_000_000); } +#[test] +fn f64_inclusive() { + let mut r = fastrand::Rng::with_seed(0); + let mut count_top_half = 0; + for _ in 0..100_000_000 { + let x = r.f64(); + assert!(x >= 0.0 && x <= 1.0); + if x > 0.5 { + count_top_half += 1; + } + } + assert!(count_top_half >= 49_000_000); +} + #[test] fn fill() { let mut r = fastrand::Rng::new(); From 25a58a19e06fc8cc890712ee21d38826b4b454fe Mon Sep 17 00:00:00 2001 From: Reiner Pope Date: Sat, 10 May 2025 22:36:46 -0700 Subject: [PATCH 3/8] Apply clippy suggestions --- tests/smoke.rs | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/tests/smoke.rs b/tests/smoke.rs index 43475fd..ec5e856 100644 --- a/tests/smoke.rs +++ b/tests/smoke.rs @@ -83,7 +83,7 @@ fn f32() { let mut count_top_half = 0; for _ in 0..100_000_000 { let x = r.f32(); - assert!(x >= 0.0 && x < 1.0); + assert!((0.0..1.0).contains(&x)); if x > 0.0 && x < tiny { count_tiny_nonzero += 1; } else if x > 0.5 { @@ -103,7 +103,7 @@ fn f32_inclusive() { let mut count_one = 0; for _ in 0..100_000_000 { let x = r.f32_inclusive(); - assert!(x >= 0.0 && x <= 1.0); + assert!((0.0..=1.0).contains(&x)); if x == 1.0 { count_one += 1; } else if x > 0.5 { @@ -123,7 +123,7 @@ fn f64() { let mut count_top_half = 0; for _ in 0..100_000_000 { let x = r.f64(); - assert!(x >= 0.0 && x < 1.0); + assert!((0.0..1.0).contains(&x)); if x > 0.5 { count_top_half += 1; } @@ -136,8 +136,8 @@ fn f64_inclusive() { let mut r = fastrand::Rng::with_seed(0); let mut count_top_half = 0; for _ in 0..100_000_000 { - let x = r.f64(); - assert!(x >= 0.0 && x <= 1.0); + let x = r.f64_inclusive(); + assert!((0.0..=1.0).contains(&x)); if x > 0.5 { count_top_half += 1; } From 0c328caa29cbeab3c96cc90cfc978af134ff8e6e Mon Sep 17 00:00:00 2001 From: Reiner Pope Date: Sat, 10 May 2025 22:38:07 -0700 Subject: [PATCH 4/8] Add some more commentary. --- src/global_rng.rs | 10 ++++++++++ src/lib.rs | 9 +++++++++ 2 files changed, 19 insertions(+) diff --git a/src/global_rng.rs b/src/global_rng.rs index 0db6c8d..11a4ff6 100644 --- a/src/global_rng.rs +++ b/src/global_rng.rs @@ -176,11 +176,21 @@ pub fn f32() -> f32 { with_rng(|r| r.f32()) } +/// Generates a random `f32` in range `0..=1`. +pub fn f32_inclusive() -> f32 { + with_rng(|r| r.f32_inclusive()) +} + /// Generates a random `f64` in range `0..1`. pub fn f64() -> f64 { with_rng(|r| r.f64()) } +/// Generates a random `f64` in range `0..=1`. +pub fn f64_inclusive() -> f64 { + with_rng(|r| r.f64_inclusive()) +} + /// Collects `amount` values at random from the iterable into a vector. pub fn choose_multiple(source: I, amount: usize) -> Vec { with_rng(|rng| rng.choose_multiple(source, amount)) diff --git a/src/lib.rs b/src/lib.rs index 17d3eec..ed4d5de 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -383,10 +383,16 @@ impl Rng { // is no larger than the bias in the underlying WyRand generator: since it only // has a 64-bit state, it necessarily already have biases of at least 2^(-64) // probability. + // + // See e.g. Section 3.1 of Thomas, David B., et al. "Gaussian random number generators, + // https://www.doc.ic.ac.uk/~wl/papers/07/csur07dt.pdf, for background. (self.u64(..) >> 1) as f32 * (-63.0f32).exp2() } /// Generates a random `f32` in range `0..1`. + /// + /// Function `f32_inclusive()` is a little simpler and faster, so default + /// to that if inclusive range is acceptable. pub fn f32(&mut self) -> f32 { loop { let x = self.f32_inclusive(); @@ -403,6 +409,9 @@ impl Rng { } /// Generates a random `f64` in range `0..1`. + /// + /// Function `f64_inclusive()` is a little simpler and faster, so default + /// to that if inclusive range is acceptable. pub fn f64(&mut self) -> f64 { loop { let x = self.f64_inclusive(); From 5b7946f89e96361437b6e21ed972ff0e93443f52 Mon Sep 17 00:00:00 2001 From: Reiner Pope Date: Sat, 10 May 2025 22:44:00 -0700 Subject: [PATCH 5/8] cargo fmt --- src/lib.rs | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/lib.rs b/src/lib.rs index ed4d5de..a224dbe 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -390,7 +390,7 @@ impl Rng { } /// Generates a random `f32` in range `0..1`. - /// + /// /// Function `f32_inclusive()` is a little simpler and faster, so default /// to that if inclusive range is acceptable. pub fn f32(&mut self) -> f32 { @@ -402,14 +402,14 @@ impl Rng { } } - /// Generates a random `f64` in range `0..=1`. + /// Generates a random `f64` in range `0..=1`. pub fn f64_inclusive(&mut self) -> f64 { // See the comment in f32_inclusive() for more details. (self.u64(..) >> 1) as f64 * (-63.0f64).exp2() } /// Generates a random `f64` in range `0..1`. - /// + /// /// Function `f64_inclusive()` is a little simpler and faster, so default /// to that if inclusive range is acceptable. pub fn f64(&mut self) -> f64 { From 9ae06081e71e93bfb9a8b17a7344b6c9c5a556ea Mon Sep 17 00:00:00 2001 From: Reiner Pope Date: Sat, 10 May 2025 22:47:29 -0700 Subject: [PATCH 6/8] Remove use of f32::exp2(). This is unavailable on thumbv7m-none-eabi. --- src/lib.rs | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/lib.rs b/src/lib.rs index a224dbe..b0be0a8 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -386,7 +386,8 @@ impl Rng { // // See e.g. Section 3.1 of Thomas, David B., et al. "Gaussian random number generators, // https://www.doc.ic.ac.uk/~wl/papers/07/csur07dt.pdf, for background. - (self.u64(..) >> 1) as f32 * (-63.0f32).exp2() + const MUL: f32 = 1.0 / (1u64 << 63) as f32; + (self.u64(..) >> 1) as f32 * MUL } /// Generates a random `f32` in range `0..1`. @@ -405,7 +406,8 @@ impl Rng { /// Generates a random `f64` in range `0..=1`. pub fn f64_inclusive(&mut self) -> f64 { // See the comment in f32_inclusive() for more details. - (self.u64(..) >> 1) as f64 * (-63.0f64).exp2() + const MUL: f64 = 1.0 / (1u64 << 63) as f64; + (self.u64(..) >> 1) as f64 * MUL } /// Generates a random `f64` in range `0..1`. From d571112899a21cfd6d13ecd0b5dd477c402f124e Mon Sep 17 00:00:00 2001 From: Reiner Pope Date: Sat, 10 May 2025 22:53:35 -0700 Subject: [PATCH 7/8] Add #[inline] --- src/lib.rs | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/lib.rs b/src/lib.rs index b0be0a8..a59316e 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -365,6 +365,7 @@ impl Rng { } /// Generates a random `f32` in range `0..=1`. + #[inline] pub fn f32_inclusive(&mut self) -> f32 { // Generate a number in 0..2^63 then convert to f32 and multiply by 2^(-63). // @@ -394,6 +395,7 @@ impl Rng { /// /// Function `f32_inclusive()` is a little simpler and faster, so default /// to that if inclusive range is acceptable. + #[inline] pub fn f32(&mut self) -> f32 { loop { let x = self.f32_inclusive(); @@ -404,6 +406,7 @@ impl Rng { } /// Generates a random `f64` in range `0..=1`. + #[inline] pub fn f64_inclusive(&mut self) -> f64 { // See the comment in f32_inclusive() for more details. const MUL: f64 = 1.0 / (1u64 << 63) as f64; @@ -414,6 +417,7 @@ impl Rng { /// /// Function `f64_inclusive()` is a little simpler and faster, so default /// to that if inclusive range is acceptable. + #[inline] pub fn f64(&mut self) -> f64 { loop { let x = self.f64_inclusive(); From 7c3b83b78140d6a74b1283dc4055d6b59307fc35 Mon Sep 17 00:00:00 2001 From: Reiner Pope Date: Sat, 10 May 2025 23:11:34 -0700 Subject: [PATCH 8/8] Switch from u64(..) to gen_u64(). They both end up compiling to the same after inlining and simplification, but the former requires more work from the optimizer to get there. --- src/lib.rs | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/lib.rs b/src/lib.rs index a59316e..ba9a61a 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -388,7 +388,7 @@ impl Rng { // See e.g. Section 3.1 of Thomas, David B., et al. "Gaussian random number generators, // https://www.doc.ic.ac.uk/~wl/papers/07/csur07dt.pdf, for background. const MUL: f32 = 1.0 / (1u64 << 63) as f32; - (self.u64(..) >> 1) as f32 * MUL + (self.gen_u64() >> 1) as f32 * MUL } /// Generates a random `f32` in range `0..1`. @@ -410,7 +410,7 @@ impl Rng { pub fn f64_inclusive(&mut self) -> f64 { // See the comment in f32_inclusive() for more details. const MUL: f64 = 1.0 / (1u64 << 63) as f64; - (self.u64(..) >> 1) as f64 * MUL + (self.gen_u64() >> 1) as f64 * MUL } /// Generates a random `f64` in range `0..1`.