From d0beebff8db18641bbb7b73e3ef7d014f3f8aaaa Mon Sep 17 00:00:00 2001 From: portlandHODL Date: Thu, 27 Aug 2026 07:45:51 -0700 Subject: [PATCH] positive_f64: saturate arithmetic at MIN_POSITIVE, define inf/inf as 1 PositiveF64 guarantees that it holds a positive, non-NaN f64, but its arithmetic could escape this invariant at the extremes: inf / inf produced NaN (silently breaking the Eq contract), a finite value divided by inf produced 0.0, and repeated division (e.g. through deeply-nested-or policies) could underflow a probability to 0.0. Harden the arithmetic as described in issue #1033: - Multiplication and division results that would underflow to a subnormal number or to 0.0 now saturate at f64::MIN_POSITIVE, exposed as PositiveF64::MIN_POSITIVE. Once this value is reached, further division (or multiplication by values less than one) is a no-op. - inf / inf is defined to be 1.0, making NaN unreachable and the Eq implementation sound. - A finite value divided by inf underflows to 0.0 and therefore saturates at MIN_POSITIVE like any other underflow. - NormalizedIterator uses the same hardened division, closing the compiler-side underflow path. Document the invariant and saturation semantics on the type (replacing the stale "Ordered f64 for comparison" doc line left over from OrdF64), document that NormalizedIterator yields mathematically incorrect results when its input contains multiple infinities, and add targeted unit tests for the extreme-value behavior, including the reversed-iteration path the compiler uses. Fixes #1033 --- src/primitives/positive_f64.rs | 134 +++++++++++++++++++++++++++++++-- 1 file changed, 129 insertions(+), 5 deletions(-) diff --git a/src/primitives/positive_f64.rs b/src/primitives/positive_f64.rs index fd7dc7b61..f7995b439 100644 --- a/src/primitives/positive_f64.rs +++ b/src/primitives/positive_f64.rs @@ -7,7 +7,22 @@ use core::{cmp, hash, ops}; use crate::Threshold; -/// Ordered f64 for comparison. +/// A positive floating-point number. +/// +/// This type guarantees that the contained value is positive: it will never +/// hold 0.0, a negative number, `-inf` or `NaN`. (Positive infinity *is* a +/// permissible value.) This guarantee makes it safe to implement [`Eq`], even +/// though the underlying [`PartialEq`] implementation passes through to `f64`. +/// +/// To uphold the guarantee, arithmetic on this type saturates below at +/// [`PositiveF64::MIN_POSITIVE`]: any operation whose result would underflow +/// to a subnormal number or to 0.0 instead yields [`PositiveF64::MIN_POSITIVE`]. +/// This means that once you obtain [`PositiveF64::MIN_POSITIVE`], dividing it +/// further (or multiplying it by values less than one) is a no-op. +/// +/// Division involving infinity is also defined so as to avoid `NaN` and 0.0: +/// `inf / inf` is defined to be 1.0, and ` / inf` underflows to 0.0 +/// and therefore yields [`PositiveF64::MIN_POSITIVE`]. #[derive(Copy, Clone, PartialEq, Debug)] pub struct PositiveF64(f64); @@ -15,6 +30,13 @@ impl PositiveF64 { /// The constant one. pub const ONE: Self = Self(1.0); + /// The smallest value that arithmetic on [`PositiveF64`] can produce. + /// + /// This is the smallest positive normal `f64`. Operations whose results + /// would underflow below this value saturate here instead; see the type + /// documentation for more detail. + pub const MIN_POSITIVE: Self = Self(f64::MIN_POSITIVE); + /// Constant used in unit tsets #[cfg(test)] pub const ONE_QUARTER: Self = Self(0.25); @@ -40,6 +62,14 @@ impl PositiveF64 { /// Internally clones the iterator and runs it twice, so best to only use /// this with reference-based iterators obtained with e.g. `slice.iter()` /// rather than "owning" iterators like you'd get from `vec.into_iter()`. + /// + /// Note that the result is mathematically incorrect if the iterator + /// contains multiple infinities: since `inf / inf` is defined to be 1.0, + /// every infinite item normalizes to 1.0 and the items total to the + /// number of infinities rather than to 1. (The correct behavior would be + /// to yield `1 / <# infinities>` for each one, but that would require an + /// extra counting pass, and this case never occurs in this crate's own + /// usage.) pub fn normalized_iter(iter: I) -> NormalizedIterator where I: Iterator + Clone, @@ -103,6 +133,36 @@ impl From for PositiveF64 { fn from(value: NonZeroU32) -> Self { Self(f64::from(u32::from(value))) } } +/// Clamps the result of an arithmetic operation on positive floats to the +/// range of permissible [`PositiveF64`] values. +/// +/// Positive floats cannot produce `NaN` or negative values when added or +/// multiplied together, and can only produce `NaN` when divided in the +/// `inf / inf` case (which is special-cased in [`positive_div`]). They can, +/// however, underflow to a subnormal number or to 0.0; in this case we +/// saturate at `f64::MIN_POSITIVE`. +fn clamp_positive(f: f64) -> f64 { + // As described above, `f` being NaN is impossible here. (If we did get a + // NaN, `f64::max` would turn it into `f64::MIN_POSITIVE` anyway.) + f.max(f64::MIN_POSITIVE) +} + +/// Multiplies two positive floats, saturating below at `f64::MIN_POSITIVE`. +fn positive_mul(a: f64, b: f64) -> f64 { clamp_positive(a * b) } + +/// Divides two positive floats, saturating below at `f64::MIN_POSITIVE`. +/// +/// Defines `inf / inf` to be 1.0 to avoid producing `NaN`. A finite value +/// divided by `inf` underflows to 0.0 and therefore saturates at +/// `f64::MIN_POSITIVE`, like any other underflowing division. +fn positive_div(a: f64, b: f64) -> f64 { + if a.is_infinite() && b.is_infinite() { + 1.0 + } else { + clamp_positive(a / b) + } +} + macro_rules! impl_op { ($trait:ident, $op:ident, $expr:expr) => { impl ops::$trait for PositiveF64 { @@ -128,9 +188,14 @@ macro_rules! impl_op { } impl_op!(Add, add, f64::add); -impl_op!(Mul, mul, f64::mul); -impl_op!(Div, div, f64::div); +impl_op!(Mul, mul, positive_mul); +impl_op!(Div, div, positive_div); +/// Iterator over [`PositiveF64`]s normalized to total to 1. +/// +/// Constructed by [`PositiveF64::normalized_iter`]; see its documentation for +/// details, including a note on mathematically incorrect behavior when the +/// input contains multiple infinities. pub struct NormalizedIterator { iter: I, /// Sum must be nonnegative, and may only be zero if `iter` is empty. @@ -143,7 +208,9 @@ where { type Item = I::Item; fn next(&mut self) -> Option { - self.iter.next().map(|x| PositiveF64(x.0 / self.sum)) + self.iter + .next() + .map(|x| PositiveF64(positive_div(x.0, self.sum))) } fn size_hint(&self) -> (usize, Option) { self.iter.size_hint() } @@ -154,7 +221,9 @@ where I: Iterator + DoubleEndedIterator, { fn next_back(&mut self) -> Option { - self.iter.next_back().map(|x| PositiveF64(x.0 / self.sum)) + self.iter + .next_back() + .map(|x| PositiveF64(positive_div(x.0, self.sum))) } } @@ -167,3 +236,58 @@ impl FusedIterator for NormalizedIterator where I: Iterator + FusedIterator { } + +#[cfg(test)] +mod tests { + use super::PositiveF64; + + #[test] + fn div_inf_by_inf_is_one() { + let inf = PositiveF64::new(f64::INFINITY).unwrap(); + assert_eq!(inf / inf, PositiveF64::ONE); + } + + #[test] + fn arithmetic_saturates_at_min_positive() { + let inf = PositiveF64::new(f64::INFINITY).unwrap(); + let two = PositiveF64::new(2.0).unwrap(); + let min = PositiveF64::MIN_POSITIVE; + // A finite value divided by inf underflows to 0.0 and saturates. + assert_eq!(PositiveF64::ONE / inf, min); + // Dividing MIN_POSITIVE further would underflow to a subnormal + // number; the division is instead a no-op. + assert_eq!(min / two, min); + // Multiplication saturates the same way. + assert_eq!(min * min, min); + } + + #[test] + fn ordinary_arithmetic_unchanged() { + let two = PositiveF64::new(2.0).unwrap(); + let three = PositiveF64::new(3.0).unwrap(); + assert_eq!(two * three, PositiveF64::new(6.0).unwrap()); + assert_eq!(three / two, PositiveF64::new(1.5).unwrap()); + } + + #[test] + fn normalized_iter_extreme_values() { + use super::PositiveF64 as P; + + // inf and inf: the sum is inf, and each item normalizes to 1.0. + // This is mathematically incorrect (the items total to 2, not 1); + // see the documented limitation on `normalized_iter`. + let infs = [ + P::new(f64::INFINITY).unwrap(), + P::new(f64::INFINITY).unwrap(), + ]; + let normalized: Vec

= P::normalized_iter(infs.iter().copied()).collect(); + assert_eq!(normalized, vec![P::ONE, P::ONE]); + + // A tiny value next to an infinite one: saturates at MIN_POSITIVE. + // Collected in reverse to also exercise `next_back`, which the + // compiler uses via `rev()`. + let mixed = [P::MIN_POSITIVE, P::new(f64::INFINITY).unwrap()]; + let normalized: Vec

= P::normalized_iter(mixed.iter().copied()).rev().collect(); + assert_eq!(normalized, vec![P::ONE, P::MIN_POSITIVE]); + } +}