diff --git a/float/CHANGELOG.md b/float/CHANGELOG.md index 4d450cc2..99c8d2bb 100644 --- a/float/CHANGELOG.md +++ b/float/CHANGELOG.md @@ -8,6 +8,14 @@ it, and the complex literal macro uses it. ### Fix +- `powi`'s overflow guard no longer reports a spurious overflow when the base is very close to 1 + (a large significand with a large negative exponent). It estimated `log2(base)` with `log2_est`, + which for such a base is the difference of two ~1e3-magnitude terms and catastrophically cancels + to ~1e-4 of `f32` noise; scaled by a large exponent that noise crossed the overflow threshold + — which on 32-bit targets is only `isize::MAX·log2(B) ≈ 7e9` — and returned a spurious ±inf. This + made high-precision base conversion (`FBig::with_base`) panic ("arithmetic operations with the + infinity are not allowed!") on 32-bit targets (wasm32, i686). The guard now uses the + bit-length-based `log2_bounds`, which does not cancel (#95). - `to_f64`/`to_f32` now round the source once, directly to the target's precision at its own magnitude (fewer than 53/24 bits for subnormals), instead of through a fixed 53/24-bit intermediate that re-rounds into the subnormal grid. This removes a 1-ULP double-rounding error diff --git a/float/src/convert.rs b/float/src/convert.rs index a73d1e1e..99f75e3d 100644 --- a/float/src/convert.rs +++ b/float/src/convert.rs @@ -1043,4 +1043,46 @@ mod tests { let int_val = FBig::::new(Repr::new(IBig::from(42), 0), Context::new(8)); assert_eq!(IBig::try_from(int_val), Ok(IBig::from(42))); } + + #[test] + fn with_base_high_precision_no_overflow() { + // Regression for issue #95: converting a high-precision base-2 float to base + // 10 panicked on 32-bit targets ("arithmetic operations with the infinity are + // not allowed!"). The base conversion evaluates exp(r) as `sum^(B^n)` through + // `powi` with a huge exponent (B^n) on a base (sum) very close to 1; `powi`'s + // overflow guard estimated log2(base) with the catastrophically-canceling + // `log2_est`, and the ~1e-4 of f32 noise scaled by the exponent crossed the + // (much smaller on 32-bit) isize threshold, yielding a spurious ±inf that then + // panicked when shifted. See `powi` in exp.rs for the fix. + use crate::round::mode::Zero; + use core::str::FromStr; + + // The reporter's input: -1.1111…0011 in binary (578 significant bits), written + // in the hex form dashu accepts for base-2 floats (`0x1.…`). The value is + // identical to the raw binary literal. + let num = FBig::::from_str( + "-0x1.fffdc8d645194a5a95df4be063472d4406dd096339dd7dc2a8527d208b3da7b9e5c36b4f49a7982cb2ad20a4e7e4c016f858fe8cddea011a6d01fe3823189c4ed4f57a7babc331498", + ) + .unwrap(); + + // at the original 578-bit precision the conversion succeeds … + let a = num + .clone() + .with_precision(578) + .value() + .with_base::<10>() + .value(); + assert!(a.repr().is_finite()); + // … and so does a slightly higher precision (586), which panicked on 32-bit + // (wasm32 / i686). The result matches the value computed on 64-bit. Compared + // by value (FBig equality ignores context) rather than via string formatting, + // so this works under no_std too. + let b = num.with_precision(586).value().with_base::<10>().value(); + assert!(b.repr().is_finite()); + let expected = FBig::::from_str( + "-1.9999661944503703041843468850635057967553124154072485151176192294480158424234268438137612977886891381228704640656094986435381057574477216648567249609280392009533217665484389886", + ) + .unwrap(); + assert_eq!(b, expected); + } } diff --git a/float/src/exp.rs b/float/src/exp.rs index 9d8152f3..dd8ca617 100644 --- a/float/src/exp.rs +++ b/float/src/exp.rs @@ -133,32 +133,50 @@ impl Context { return Ok(repr.map(|v| FBig::new(v, *self))); } - // Guard against exponent overflow for astronomically large results: the result - // magnitude has log2 ≈ exp·log2(base); if that exceeds the isize exponent range, - // return ±inf (|base| > 1) or 0 (|base| < 1) instead of overflowing mid-computation. - let base_log2 = base.log2_est() as f64; + // Guard against exponent overflow/underflow for astronomically large results. + // The result magnitude has log2 ≈ exp·log2(base); if that leaves the isize + // exponent range, return ±inf (|base|>1) or signed 0 (|base|<1) instead of + // overflowing the exponent mid-computation. + // + // Use the *bounds* of log2(base), never the point estimate `log2_est`: when + // base is very close to 1 (a large significand with a large negative exponent), + // log2(base) is the difference of two ~1e3-magnitude terms and suffers + // catastrophic cancellation — `log2_est` returns ~1e-4 of f32 noise rather than + // ~0. Scaled by a large exponent that noise crosses the overflow threshold, + // which on 32-bit is only isize::MAX·log2(B) ≈ 7e9 (vs ≈3e19 on 64-bit), so the + // guard fires spuriously and returns ±inf — see issue #95. The bounds are + // derived from the exact significand bit length, so they don't cancel. Declare + // overflow only on the lower bound (no false positives) and underflow only on + // the upper bound (no false underflows); anything in between is computed. + let (base_log2_lb, base_log2_ub) = base.log2_bounds(); + let base_log2_lb = base_log2_lb as f64; + let base_log2_ub = base_log2_ub as f64; let threshold = (isize::MAX as f64) * (B.log2_est() as f64); let exp_f64 = i64::try_from(&exp).ok().map(|e| e as f64); let overflows = match exp_f64 { - Some(e) => e * base_log2 > threshold, - None => base_log2 != 0.0, // exp doesn't fit i64: overflows unless |base| == 1 + Some(e) => e * base_log2_lb > threshold, + None => base_log2_lb > 0.0, // exp doesn't fit i64: overflows iff base > 1 }; if overflows { - return if base_log2 > 0.0 { - Err(FpError::Overflow(if base.sign() == Sign::Negative { - Sign::Negative - } else { - Sign::Positive - })) + return Err(FpError::Overflow(if base.sign() == Sign::Negative { + Sign::Negative } else { - // |base| < 1 and exponent huge → underflow to signed zero - let underflow_sign = if base.sign() == Sign::Negative && exp.bit(0) { - Sign::Negative - } else { - Sign::Positive - }; - Err(FpError::Underflow(underflow_sign)) + Sign::Positive + })); + } + let underflows = !base.significand.is_zero() + && match exp_f64 { + Some(e) => e * base_log2_ub < -threshold, + None => base_log2_ub < 0.0, // exp doesn't fit i64: underflows iff base < 1 + }; + if underflows { + // |base| < 1 and exponent huge → underflow to signed zero + let underflow_sign = if base.sign() == Sign::Negative && exp.bit(0) { + Sign::Negative + } else { + Sign::Positive }; + return Err(FpError::Underflow(underflow_sign)); } let work_context = if self.is_limited() {