From 4be3e260e51ef0cfc0a84869dfe49e8143f0d91c Mon Sep 17 00:00:00 2001 From: Jacob Zhong Date: Thu, 23 Jul 2026 09:55:04 +0800 Subject: [PATCH 1/5] Address opendp-num fuzz findings (DASHU-007/008/015/022/026) Reviewed the opendp-num differential/property-fuzz findings against the current source. DASHU-020/021 are already fixed by PR #91; DASHU-023/024 (directed rounding at FBig's exponent-range extremes) are tracked with TODO comments and deferred. This commit fixes the remaining five. DASHU-015 (dashu-int, incorrect-result): `to_f64_small` reported an inexact conversion as `Approximation::Exact` at the `DoubleWord::MAX` boundary. The exactness test used a saturating `f as DoubleWord` round-trip, which clamps back to `DoubleWord::MAX` when the value rounds up to `2^BITS`. Detect that saturation case before the comparison; apply the same guard to the 16/32-bit RefLarge -> f64 fast path. DASHU-026 (dashu-float, debug panic/OOM): `round_fract`'s debug assert built `B^precision`; for a sparse sticky tail (precision == the exponent gap, e.g. from `exp_m1` of a large-magnitude input) this exhausted memory in debug builds. Replace it with a `log2_bounds`-based precondition check that allocates nothing and only fires on a proven violation. DASHU-022 (dashu-float, panic): `assert_limited_precision` fired before the exact zero/one shortcuts, so values carrying unlimited precision (precision 0) -- `FBig::try_from(0.0)`, `FBig::ONE`/`ZERO` -- panicked in exp/exp_m1/ sqrt/ln/ln_1p despite mathematically exact results. Hoist the exact shortcuts above the assertion. DASHU-008 (dashu-int, panic): GCD of two large, similarly-sized integers aborted with "internal error: not enough memory allocated". The scratchpad was reserved once from the initial operand lengths, but each euclidean step's division dispatches on the current lengths, so a later lopsided Burnikel-Ziegler step (divisor > 48 words, quotient > 32 words) was under-reserved. Size it for the worst case reachable in the loop: `mul::memory_requirement_up_to(rhs_len, rhs_len / 2)`. Applies to both the gcd and extended-gcd reservations. DASHU-007 (dashu-float, missing feature): add a correctly-rounded `FBig::log2` / `Context::log2` / `CachedFBig::log2`, computed as `ln(x)/ln(2)` at an elevated working precision. Previously only the f32-precision `log2_bounds` estimate existed, so directed log2 was wrong by many ULPs. The result magnitude tracks the division's error amplification, so a few guard digits certify the final round across the whole range. Each fix includes a regression test; the DASHU-008 test reproduces the exact memory.rs panic under the old code. Co-Authored-By: Claude --- float/CHANGELOG.md | 9 ++ float/src/exp.rs | 63 ++++++++++++- float/src/fbig_cached.rs | 1 + float/src/fbig_cached_ops.rs | 1 + float/src/log.rs | 175 ++++++++++++++++++++++++++++++++++- float/src/root.rs | 6 +- float/src/round.rs | 16 +++- integer/CHANGELOG.md | 13 +++ integer/src/convert.rs | 26 +++++- integer/src/gcd/lehmer.rs | 24 ++++- integer/tests/convert.rs | 22 +++++ integer/tests/gcd.rs | 37 ++++++++ rational/CHANGELOG.md | 7 ++ rational/tests/arith_prop.rs | 32 +++++++ 14 files changed, 415 insertions(+), 17 deletions(-) diff --git a/float/CHANGELOG.md b/float/CHANGELOG.md index 4d450cc2..29b41346 100644 --- a/float/CHANGELOG.md +++ b/float/CHANGELOG.md @@ -6,6 +6,9 @@ - `Repr::new_const`: a `const`-evaluable, normalized `Repr` constructor from a `DoubleWord` significand (the `const` counterpart of `Repr::new`). `FBig::from_parts_const` now delegates to it, and the complex literal macro uses it. +- `FBig::log2` / `Context::log2` / `CachedFBig::log2`: base-2 logarithm, correctly rounded via + `ln(x)/ln(2)` evaluated at an elevated working precision. Previously only the f32-precision + `log2_bounds` magnitude estimate was available, so directed `log2` was wrong by many ULPs. ### Fix - `to_f64`/`to_f32` now round the source once, directly to the target's precision at its own @@ -16,6 +19,12 @@ high-precision inputs: the base-changing conversion now pre-shrinks the source significand before dividing, upholding `repr_div`'s dividend-width contract (as `Context::div` already does) instead of feeding an oversized dividend into the division. +- `exp`, `exp_m1`, `sqrt`, `ln`, and `ln_1p` no longer panic on exact zero/one inputs that carry + unlimited precision (precision 0), such as `FBig::try_from(0.0)` and the `FBig::ONE`/`ZERO` + constants. The exact-result shortcuts now run before the limited-precision assertion. +- The `round_fract` debug assertion no longer materializes `B^precision`, which for a sparse sticky + tail (where `precision` is the exponent gap — e.g. `exp_m1` of a large-magnitude input) could + exhaust memory in debug builds. It now checks the precondition with `log2` bounds. ## 0.5.0 diff --git a/float/src/exp.rs b/float/src/exp.rs index 9d8152f3..5c3040ef 100644 --- a/float/src/exp.rs +++ b/float/src/exp.rs @@ -109,6 +109,11 @@ impl Context { /// Panics if the precision is unlimited and the exponent is negative. In this case, the exact /// result is likely to have infinite digits. pub fn powi(&self, base: &Repr, exp: IBig) -> FpResult> { + // TODO(DASHU-024): range handling has three known limitations at the exponent extremes: + // (1) the overflow guard below estimates the result magnitude with an f64 and misclassifies + // representable boundaries near 2^63; (2) genuine Overflow/Underflow is unwrapped mode-blindly + // to ±inf / signed zero; (3) the negative-exponent reciprocal path can panic. None affects + // ordinary inputs; fixing requires mode-aware range saturation. if base.is_infinite() { return Err(FpError::InfiniteInput); } @@ -335,11 +340,13 @@ impl Context { mut cache: Option<&mut ConstCache>, ) -> FpResult> { assert_finite(x); - assert_limited_precision(self.precision); let input_sign = x.sign(); if x.significand.is_zero() { - // exp(±0) = 1; exp_m1(±0) = ±0 (IEEE 754 §9.2.1 preserves the sign of zero) + // exp(±0) = 1; exp_m1(±0) = ±0 (IEEE 754 §9.2.1 preserves the sign of zero). + // These exact results need no rounding, so handle them before the + // limited-precision assertion: a precision-0 (unlimited) FBig such as the + // one produced by `try_from(0.0)` must still compute exp/exp_m1 exactly. return match minus_one { false => Ok(Exact(FBig::ONE)), true => { @@ -353,6 +360,8 @@ impl Context { }; } + assert_limited_precision(self.precision); + // A simple algorithm: // - let r = (x - s logB) / Bⁿ, where s = floor(x / logB), such that r < B⁻ⁿ. // - if the target precision is p digits, then there're only about p/m terms in Tyler series @@ -399,6 +408,12 @@ impl Context { Err(_) => { // |floor(x / ln B)| overflows isize — x is astronomically large, so the // result is an infinity (x → +∞) or underflows to the limit (x → −∞). + // + // TODO(DASHU-023): this branch discards the rounding mode. For a huge + // *negative* x the true exp(x) is positive but below the exponent range, so + // directed Up/Away should return the minimum positive value (and exp_m1 the + // value just above -1) rather than +0 / Exact(-1). The mode-blind saturation + // here is a known limitation. return if input_sign == Sign::Positive { Err(FpError::Overflow(Sign::Positive)) } else if minus_one { @@ -473,6 +488,50 @@ mod tests { assert_eq!(m1, -FBig::::ONE); } + #[test] + fn test_exp_m1_large_negative_no_oom() { + // Regression for DASHU-026: exp_m1 of a large-magnitude negative input subtracts + // 1 from an astronomically small exp(x), so the aligned subtraction hands + // round_fract a sparse sticky tail whose `precision` equals the exponent gap + // (~6.6e18 here). The old debug assert materialized B^precision (2^6.6e18 bits → + // certain OOM); the new check uses log2 bounds and must not allocate. + // + // x = -2^62 is large enough that the exponent gap dwarfs available memory, yet + // small enough that floor(x/ln2) still fits isize (≈6.6e18 < 9.2e18), so it does + // NOT enter the overflow branch tested above. + let ctx = Context::::new(2); + let x = Repr::new(-(IBig::from(1) << 62), 0); + let r = ctx.exp_m1::<2>(&x, None).unwrap().value(); + // exp_m1(-2^62) = -1 + ε (ε > 0 vanishingly small); under Up it rounds above -1. + assert!( + r > -FBig::::ONE, + "exp_m1(-2^62) under Up should round above -1, got {r:?}" + ); + } + + #[test] + fn test_exact_results_on_unlimited_precision() { + // Regression for DASHU-022: values carrying precision 0 (unlimited) — produced + // by `try_from(0.0)` and the `FBig::ONE`/`ZERO` constants — must still compute + // their exact-result special cases instead of panicking in + // assert_limited_precision before reaching the shortcut. + type F = FBig; + + let zero = F::try_from(0.0_f64).unwrap(); + assert_eq!(zero.exp(), F::ONE); + assert_eq!(zero.exp_m1(), F::ZERO); + assert_eq!(zero.sqrt(), F::ZERO); + assert_eq!(zero.ln_1p(), F::ZERO); + + // -0.0 preserves its sign through exp_m1 and sqrt. + let neg_zero = F::try_from(-0.0_f64).unwrap(); + assert!(neg_zero.exp_m1().repr().is_neg_zero()); + assert!(neg_zero.sqrt().repr().is_neg_zero()); + + // FBig::ONE carries unlimited precision; ln(1) = 0 is exact. + assert_eq!(F::ONE.ln(), F::ZERO); + } + #[test] fn test_powf_zero_base() { use crate::DBig; diff --git a/float/src/fbig_cached.rs b/float/src/fbig_cached.rs index 1ea2a645..87bebc33 100644 --- a/float/src/fbig_cached.rs +++ b/float/src/fbig_cached.rs @@ -539,6 +539,7 @@ mod tests { assert_eq!(x.clone().cos().into_fbig(), y.clone().cos()); assert_eq!(x.clone().exp_m1().into_fbig(), y.clone().exp_m1()); assert_eq!(x.clone().ln_1p().into_fbig(), y.clone().ln_1p()); + assert_eq!(x.clone().log2().into_fbig(), y.clone().log2()); assert_eq!(x.powf(&x.clone()).into_fbig(), y.clone().powf(&y)); } diff --git a/float/src/fbig_cached_ops.rs b/float/src/fbig_cached_ops.rs index c3f44e4c..3ac82869 100644 --- a/float/src/fbig_cached_ops.rs +++ b/float/src/fbig_cached_ops.rs @@ -666,6 +666,7 @@ macro_rules! forward_to_fbig { impl CachedFBig { forward_to_context!(ln); forward_to_context!(ln_1p); + forward_to_context!(log2); forward_to_context!(exp); forward_to_context!(exp_m1); diff --git a/float/src/log.rs b/float/src/log.rs index 02772914..c0206cb0 100644 --- a/float/src/log.rs +++ b/float/src/log.rs @@ -97,6 +97,28 @@ impl FBig { pub fn ln_1p(&self) -> Self { self.context.unwrap_fp(self.context.ln_1p(&self.repr, None)) } + + /// Calculate the base-2 logarithm function (`log2(x)`) on the float number. + /// + /// This is evaluated as `ln(x) / ln(2)` at an elevated working precision and rounded + /// once, so it is correctly rounded to this number's precision (unlike + /// [`log2_bounds`][dashu_base::EstimatedLog2::log2_bounds], which only gives an + /// f32-precision magnitude estimate). + /// + /// # Examples + /// + /// ``` + /// # use core::str::FromStr; + /// # use dashu_base::ParseError; + /// # use dashu_float::DBig; + /// let a = DBig::from_str("8.0")?; + /// assert_eq!(a.log2(), DBig::from_str("3")?); + /// # Ok::<(), ParseError>(()) + /// ``` + #[inline] + pub fn log2(&self) -> Self { + self.context.unwrap_fp(self.context.log2(&self.repr, None)) + } } impl Context { @@ -214,6 +236,48 @@ impl Context { Ok(self.ln_internal(x, false, cache)) } + /// Calculate the base-2 logarithm function (`log2(x)`) on the float number under this + /// context. + /// + /// This is evaluated as `ln(x) / ln(2)` at an elevated working precision and rounded + /// once, so it is correctly rounded to this context's precision (unlike + /// [`log2_bounds`][dashu_base::EstimatedLog2::log2_bounds], which only gives an + /// f32-precision magnitude estimate). + /// + /// # Examples + /// + /// ``` + /// # use core::str::FromStr; + /// # use dashu_base::ParseError; + /// # use dashu_float::DBig; + /// use dashu_base::Approximation::*; + /// use dashu_float::{Context, round::mode::HalfEven}; + /// + /// let context = Context::::new(8); + /// let a = DBig::from_str("10")?; + /// // log2(10) ≈ 3.3219281 + /// assert_eq!(context.log2(&a.repr(), None).unwrap().value(), DBig::from_str("3.3219281")?); + /// # Ok::<(), ParseError>(()) + /// ``` + #[inline] + pub fn log2( + &self, + x: &Repr, + cache: Option<&mut ConstCache>, + ) -> FpResult> { + if x.is_infinite() { + return Err(FpError::InfiniteInput); + } + if x.significand.is_zero() { + // log2(±0) = -inf (a value, not an error) + return Ok(Exact(FBig::new(Repr::neg_infinity(), *self))); + } + if x.sign() == Sign::Negative { + return Err(FpError::OutOfDomain); + } + Ok(self.log2_internal(x, cache)) + } + /// Calculate the natural logarithm function (`log(x+1)`) on the float number under this context. /// /// # Examples @@ -257,8 +321,10 @@ impl Context { mut cache: Option<&mut ConstCache>, ) -> Rounded> { assert_finite(x); - assert_limited_precision(self.precision); + // Exact special cases first: they need no rounding, so a precision-0 (unlimited) + // value such as `FBig::ONE` or the one from `try_from(0.0)` must still resolve + // ln/ln_1p exactly rather than tripping the limited-precision assertion below. if !one_plus && x.is_one() { return Exact(FBig::ZERO); // ln(1) = +0 } @@ -272,6 +338,8 @@ impl Context { return Exact(zero); } + assert_limited_precision(self.precision); + // A simple algorithm: // - let log(x) = log(x/2^s) + slog2 where s = floor(log2(x)) // - such that x*2^s is close to but larger than 1 (and x*2^s < 2) @@ -344,6 +412,32 @@ impl Context { }; result.with_precision(self.precision) } + + /// `log2(x) = ln(x) / ln(2)`, correctly rounded to this context's precision. + /// + /// Both `ln(x)` and `ln(2)` are evaluated at a working precision above the target and + /// the quotient is rounded once. The result magnitude `|log2(x)| ≈ |floor(log2 x)|` + /// grows with the operand, which tracks the division's error amplification, so the + /// per-ulp error is ~`B^(target - work)`: a few guard digits certify the final round + /// across the whole magnitude range (no manual scaling needed). + fn log2_internal( + &self, + x: &Repr, + mut cache: Option<&mut ConstCache>, + ) -> Rounded> { + assert_finite(x); + if x.is_one() { + return Exact(FBig::ZERO); // log2(1) = +0 (exact; also covers unlimited precision) + } + + let guard = ceil_usize(self.precision.log2_est() / B.log2_est()) + 3; + let work = Context::::new(self.precision + guard); + let ln_x = work + .ln_internal(x, false, reborrow_cache(&mut cache)) + .value(); + let ln2 = work.ln2::(reborrow_cache(&mut cache)); + (ln_x / ln2).with_precision(self.precision) + } } #[cfg(test)] @@ -359,6 +453,85 @@ mod tests { assert_eq!(r.repr().sign(), Sign::Negative); } + #[test] + fn test_log2_domain() { + let ctx = Context::::new(53); + // log2(±0) = -inf (a value, not an error) + let r = ctx.log2::<2>(&Repr::<2>::zero(), None).unwrap().value(); + assert!(r.repr().is_infinite()); + assert_eq!(r.repr().sign(), Sign::Negative); + // log2(negative) is out of domain + assert_eq!(ctx.log2::<2>(&Repr::new(IBig::from(-1), 0), None), Err(FpError::OutOfDomain)); + } + + #[test] + fn test_log2_powers_of_two() { + // log2(2^k) = k, including log2(1) = 0. This is the headline DASHU-007 case: the + // adapter previously derived directed log2 from f32-precision log2_bounds, which + // gave a small *negative* value for log2(1) and was loose by many ULPs elsewhere. + // The correctly-rounded log2 returns the exact integer value for every power of + // two (log2(1) is tagged Exact via the shortcut; other powers are value-exact). + let ctx = Context::::new(53); + for k in [0usize, 1, 2, 3, 10, 50, 100, 200, 1000] { + let x = Repr::new(IBig::from(1) << k, 0); // 2^k + let r = ctx.log2::<2>(&x, None).unwrap(); + let expected = ctx.convert_int::<2>(IBig::from(k)).value(); + assert_eq!(r.value(), expected, "log2(2^{k}) = {k}"); + } + // log2(1) = 0 is certified Exact via the special-case shortcut. + assert!(matches!(ctx.log2::<2>(&Repr::new(IBig::from(1), 0), None).unwrap(), Exact(_))); + } + + #[test] + fn test_log2_correctly_rounded() { + // log2(x) must equal the correctly-rounded ln(x)/ln(2) oracle (computed at much + // higher precision then rounded down) across magnitude ranges — tiny, normal, and + // huge inputs, in both base 2 and base 10. + let p = 30; + let ctx = Context::::new(p); + let hi = Context::::new(200); + + let cases_base2 = [ + Repr::new(IBig::from(3), 0), + Repr::new(IBig::from(10), 0), + Repr::new(IBig::from(1), -24), // near the fuzzer's 2^-24 region + Repr::new(IBig::from(1) << 50, 0), // huge + Repr::new(IBig::from(1), -50), // tiny + ]; + for x in cases_base2 { + let got = ctx.log2::<2>(&x, None).unwrap().value(); + let want = (hi.ln::<2>(&x, None).unwrap().value() / hi.ln2::<2>(None)) + .with_precision(p) + .value(); + assert_eq!(got, want, "base-2 log2 mismatch for x={x:?}"); + } + + // base 10: log2(x) = ln(x)/ln(2) with a decimal significand + let x = Repr::new(IBig::from(123456), 0); + let got = ctx.log2::<10>(&x, None).unwrap().value(); + let want = (hi.ln::<10>(&x, None).unwrap().value() / hi.ln2::<10>(None)) + .with_precision(p) + .value(); + assert_eq!(got, want, "base-10 log2 mismatch"); + } + + #[test] + fn test_log2_directed_bounds() { + // Down rounds toward -inf, Up toward +inf, so the correctly-rounded (HalfEven) + // value must lie within [Down, Up]. + use crate::round::mode::{Down, Up}; + let p = 20; + let mid = Context::::new(p); + for x in [Repr::new(IBig::from(10), 0), Repr::new(IBig::from(3), 13)] { + let down = Context::::new(p).log2::<2>(&x, None).unwrap().value(); + let up = Context::::new(p).log2::<2>(&x, None).unwrap().value(); + let half = mid.log2::<2>(&x, None).unwrap().value(); + assert!(down <= half, "down <= half failed for x={x:?}: {down:?} vs {half:?}"); + assert!(half <= up, "half <= up failed for x={x:?}: {half:?} vs {up:?}"); + assert!(down <= up, "down <= up failed for x={x:?}"); + } + } + #[test] fn test_iacoth() { let context = Context::::new(10); diff --git a/float/src/root.rs b/float/src/root.rs index 5e806285..465517f6 100644 --- a/float/src/root.rs +++ b/float/src/root.rs @@ -90,11 +90,13 @@ impl Context { if x.is_infinite() { return Err(FpError::InfiniteInput); } - assert_limited_precision(self.precision); if x.significand.is_zero() { - // sqrt(+0) = +0, sqrt(-0) = -0 (preserve the sign of zero) + // sqrt(+0) = +0, sqrt(-0) = -0 (preserve the sign of zero). Exact, so handle + // it before the limited-precision assertion: a precision-0 (unlimited) value + // such as the one from `try_from(0.0)` must still compute sqrt(0) exactly. return Ok(Approximation::Exact(FBig::new(x.clone(), *self))); } + assert_limited_precision(self.precision); if x.sign() == Sign::Negative { return Err(FpError::OutOfDomain); } diff --git a/float/src/round.rs b/float/src/round.rs index 6aaa61ff..265a2833 100644 --- a/float/src/round.rs +++ b/float/src/round.rs @@ -2,7 +2,7 @@ use core::cmp::Ordering; use core::ops::{Add, AddAssign}; -use dashu_base::{AbsOrd, Approximation, BitTest, EstimatedLog2, Sign, Signed, UnsignedAbs}; +use dashu_base::{AbsOrd, Approximation, BitTest, EstimatedLog2, Sign, Signed}; use dashu_int::{IBig, UBig, Word}; use crate::FBig; @@ -104,8 +104,18 @@ pub trait Round: Copy { /// assuming |fract| / X^precision < 1. Return the adjustment. #[inline] fn round_fract(integer: &IBig, fract: IBig, precision: usize) -> Rounding { - // this assertion is costly, so only check in debug mode - debug_assert!(fract.clone().unsigned_abs() < UBig::from_word(B).pow(precision)); + // this assertion is costly, so only check in debug mode. + // Verify the precondition |fract| < B^precision *without* materializing B^precision: + // for a sparse sticky tail produced by aligned subtraction, `precision` is the + // exponent gap and can be astronomically large, so building B^precision (the old + // check) exhausts memory in debug builds. Use log2 bounds instead, and only flag a + // *proven* violation so a valid input never trips the assertion. + debug_assert!({ + let (lb, _ub) = fract.log2_bounds(); // bounds on log2|fract| (−inf when fract == 0) + let (_b_lb, b_ub) = B.log2_bounds(); + // certain violation (|fract| >= B^precision) only when lb >= precision * b_ub + lb < b_ub * precision as f32 + }); if fract.is_zero() { return Rounding::NoOp; diff --git a/integer/CHANGELOG.md b/integer/CHANGELOG.md index a478bd03..63053385 100644 --- a/integer/CHANGELOG.md +++ b/integer/CHANGELOG.md @@ -1,5 +1,18 @@ # Changelog +## Unreleased + +### Fix +- `to_f64` of a magnitude at `DoubleWord::MAX` (e.g. `IBig::from(-(2^128 - 1)).to_f64()` on 64-bit + targets) reported the inexact conversion as `Approximation::Exact`. The exactness test used a + saturating `f as DoubleWord` round-trip, which clamps back to `DoubleWord::MAX` when the value + rounds up to `2^BITS`; that saturation case is now detected before the comparison. The same guard + is applied to the 16/32-bit `RefLarge` → f64 fast path. +- GCD of two large, similarly-sized integers no longer aborts with "internal error: not enough memory + allocated" when a later euclidean step takes the Burnikel-Ziegler path. The scratchpad was reserved + once from the *initial* operand lengths, but each step's division dispatches on the *current* + lengths; it is now sized for the worst case reachable during reduction (`mul::memory_requirement_up_to(rhs_len, rhs_len/2)`). + ## 0.5.0 ### Add diff --git a/integer/src/convert.rs b/integer/src/convert.rs index 8ea85175..fb3383a8 100644 --- a/integer/src/convert.rs +++ b/integer/src/convert.rs @@ -1085,8 +1085,16 @@ mod repr { if self.bit_len() <= 64 { let v: u64 = self.try_to_unsigned().unwrap(); let f = v as f64; - let back = f as u64; - match back.cmp(&v) { + // Same saturation guard as `to_f64_small`: `f as u64` clamps to + // `u64::MAX` when `v` rounds up to 2^64, which would otherwise + // report `Equal` (Exact) for `v == u64::MAX`. + const TWO_POW_64: f64 = (1u64 << 63) as f64 * 2.0; + let cmp = if f >= TWO_POW_64 { + Ordering::Greater + } else { + (f as u64).cmp(&v) + }; + match cmp { Ordering::Greater => Inexact(f, Sign::Positive), Ordering::Equal => Exact(f), Ordering::Less => Inexact(f, Sign::Negative), @@ -1130,9 +1138,19 @@ mod repr { fn to_f64_small(dword: DoubleWord) -> Approximation { const_assert!((DoubleWord::MAX as f64) < f64::MAX); let f = dword as f64; - let back = f as DoubleWord; - match back.partial_cmp(&dword).unwrap() { + // `f as DoubleWord` is a saturating float→int cast: when `dword` rounds up to + // 2^BITS (the f64 just above `DoubleWord::MAX`), the cast clamps back to + // `DoubleWord::MAX == dword` and the round-trip comparison reports `Equal`, + // wrongly tagging the conversion `Exact`. `dword as f64` is at most 2^BITS, so + // reaching 2^BITS marks exactly that saturation case — the value rounded up. + const TWO_POW_BITS: f64 = (1u128 << (DoubleWord::BITS - 1)) as f64 * 2.0; + let cmp = if f >= TWO_POW_BITS { + Ordering::Greater + } else { + (f as DoubleWord).partial_cmp(&dword).unwrap() + }; + match cmp { Ordering::Greater => Inexact(f, Sign::Positive), Ordering::Equal => Exact(f), Ordering::Less => Inexact(f, Sign::Negative), diff --git a/integer/src/gcd/lehmer.rs b/integer/src/gcd/lehmer.rs index 00a4bb4e..40db187a 100644 --- a/integer/src/gcd/lehmer.rs +++ b/integer/src/gcd/lehmer.rs @@ -225,9 +225,21 @@ pub(crate) fn lehmer_step(x: &mut [Word], y: &mut [Word], a: Word, b: Word, c: W /// Temporary memory required for gcd. #[inline] pub fn memory_requirement_up_to(lhs_len: usize, rhs_len: usize) -> Layout { - // Required memory: - // - temporary space for the division in the euclidean step - div::memory_requirement_exact(lhs_len, rhs_len) + // The scratchpad is reserved ONCE from the initial operand lengths, but each + // euclidean step's division dispatches on the *current* lengths. As Lehmer shrinks + // y while x stays large, a later step can divide a still-large x by a much smaller y + // and take the Burnikel-Ziegler path — even when the first (similarly-sized) + // division was "simple" and needed no scratch. `div::memory_requirement_exact` only + // bounds a *single* division with the given lengths, so it under-reserves here. + // + // Reserve the worst case reachable in the loop instead: every in-loop division has a + // divisor R ≤ rhs_len and a quotient split min(R/2, L-R) ≤ rhs_len/2, and + // mul::memory_requirement_up_to is non-decreasing in both arguments, so + // mul::memory_requirement_up_to(rhs_len, rhs_len/2) bounds them all (independent of + // lhs_len). It returns zero_layout for small operands automatically, since mul's own + // schoolbook threshold gates it. + let _ = lhs_len; + mul::memory_requirement_up_to(rhs_len, rhs_len / 2) } pub(crate) fn gcd_in_place( @@ -332,12 +344,14 @@ fn lehmer_ext_step( pub fn memory_requirement_ext_up_to(lhs_len: usize, rhs_len: usize) -> Layout { // Required memory: // - two numbers (t0 & t1) with at most the same size as lhs, add 1 buffer word - // - temporary space for a division (for euclidean step), and later a mulitplication (for coeff update) + // - temporary space for a division (for euclidean step, see memory_requirement_up_to + // for why the loop-wide bound rather than the initial lengths), and later a + // multiplication (for coeff update) let t_words = 2 * lhs_len + 2; memory::add_layout( memory::array_layout::(t_words), memory::max_layout( - div::memory_requirement_exact(lhs_len, rhs_len), // + mul::memory_requirement_up_to(rhs_len, rhs_len / 2), // worst-case in-loop division mul::memory_requirement_up_to(lhs_len, lhs_len / 2), // for coeff update ), ) diff --git a/integer/tests/convert.rs b/integer/tests/convert.rs index 568e0dc8..88c093a2 100644 --- a/integer/tests/convert.rs +++ b/integer/tests/convert.rs @@ -434,6 +434,28 @@ fn test_to_f64() { ); } +/// Regression test for DASHU-015: a magnitude at `DoubleWord::MAX` rounds up to the +/// next power of two as f64, and the saturating `f as DoubleWord` round-trip used to +/// clamp back to `DoubleWord::MAX` and report that inexact conversion as `Exact`. +#[cfg(target_pointer_width = "64")] +#[test] +fn test_to_f64_dword_max_boundary() { + // 2^128 - 1 == DoubleWord::MAX on 64-bit Word targets; it fits a single DoubleWord, + // so this exercises the `to_f64_small` fast path. + let two_pow_128: f64 = 2f64.powi(128); + let max = ubig!(340282366920938463463374607431768211455); + assert_eq!(max.to_f64(), Inexact(two_pow_128, Positive)); + let neg_max = -ibig!(340282366920938463463374607431768211455); + assert_eq!(neg_max.to_f64(), Inexact(-two_pow_128, Negative)); + + // One below the boundary still rounds up to 2^128; the guard must report it Inexact. + assert_eq!((max - ubig!(1)).to_f64(), Inexact(two_pow_128, Positive)); + + // A clearly representable value below the saturation zone is still Exact — the + // guard must not be over-eager. + assert_eq!((ubig!(1) << 127).to_f64(), Exact(2f64.powi(127))); +} + #[cfg(not(any(target_pointer_width = "16", force_bits = "16")))] #[test] fn test_from_u64_i64_const() { diff --git a/integer/tests/gcd.rs b/integer/tests/gcd.rs index e52639bd..622ac985 100644 --- a/integer/tests/gcd.rs +++ b/integer/tests/gcd.rs @@ -166,3 +166,40 @@ fn test_gcd_0() { fn test_gcd_ext_0() { let _ = ubig!(0).gcd_ext(ubig!(0)); } + +#[test] +fn test_gcd_large_lopsided_reduction() { + // Regression for DASHU-008: the gcd scratchpad is reserved once from the *initial* + // operand lengths, but each euclidean step's division dispatches on the *current* + // lengths. The pre-fix code under-reserved when the initial operands were similar + // in size, so a later lopsided euclidean step (divisor > 48 words and quotient + // > 32 words, taking the Burninkel-Ziegler + Karatsuba path) aborted with + // "internal error: not enough memory allocated" in both debug and release. + use dashu_int::UBig; + + fn big_from_seed(mut seed: u64, words: usize) -> UBig { + let mut v = UBig::from(0u64); + for _ in 0..words { + seed ^= seed << 13; + seed ^= seed >> 7; + seed ^= seed << 17; + v = (v << 64) | UBig::from(seed | 1); // |1 keeps the top word nonzero + } + v + } + + // Build a reduction that starts from similar-sized operands (so the initial division + // is "simple" and reserves zero scratch) but reaches a lopsided Burnikel-Ziegler + // step. Let R ~100 words, Q ~50 words, L = Q*R + 1 ~150 words; then + // gcd(L+R, L): step 1 (L+R)/L has quotient 1 -> next pair (L, R) + // step 2 L/R has quotient Q (>32 words) with divisor R (>48 words) + // -> Burnikel-Ziegler + Karatsuba, the under-reserved step. + let r = big_from_seed(0x9e3779b97f4a7c15, 100); + let q = big_from_seed(0xd1b54a32d192ed03, 50); + let l = &q * &r + ubig!(1); + let a = &l + &r; // L + R + let b = l; // L + + // gcd(L+R, L) = gcd(L, R) = gcd(Q*R+1, R) = gcd(1, R) = 1. + assert_eq!((&a).gcd(&b), ubig!(1)); +} diff --git a/rational/CHANGELOG.md b/rational/CHANGELOG.md index 62bcef9c..d6075d45 100644 --- a/rational/CHANGELOG.md +++ b/rational/CHANGELOG.md @@ -1,5 +1,12 @@ # Changelog +## Unreleased + +### Fix +- Constructing or reducing an `RBig` with two large, similarly-sized operands no longer aborts with + "internal error: not enough memory allocated". The underlying integer GCD scratchpad is now sized + for the worst-case division reachable during reduction (see the `dashu-int` changelog). + ## 0.5.0 ### Fix diff --git a/rational/tests/arith_prop.rs b/rational/tests/arith_prop.rs index f9192489..493c5ae7 100644 --- a/rational/tests/arith_prop.rs +++ b/rational/tests/arith_prop.rs @@ -59,3 +59,35 @@ proptest! { prop_assert_eq!(reduced.as_relaxed(), scaled.as_relaxed()); } } + +/// Regression for DASHU-008: constructing an `RBig` reduces numerator/denominator by +/// their gcd, which must not abort with "not enough memory allocated" when that gcd +/// reduction hits a lopsided Burnikel-Ziegler division. The gcd scratchpad is sized +/// from the *initial* operand lengths, but each step dispatches on the *current* +/// lengths, so two large, similarly-sized operands (zero initial scratch) that later +/// reduce through a wide quotient under-reserved and panicked. +#[test] +fn rbig_reduce_large_lopsided() { + fn big_from_seed(mut seed: u64, words: usize) -> UBig { + let mut v = UBig::from(0u64); + for _ in 0..words { + seed ^= seed << 13; + seed ^= seed >> 7; + seed ^= seed << 17; + v = (v << 64) | UBig::from(seed | 1); // |1 keeps the top word nonzero + } + v + } + + // L = Q*R + 1, then gcd(L+R, L) reduces through L/R (huge quotient Q, large divisor R). + let r = big_from_seed(0x9e3779b97f4a7c15, 100); + let q = big_from_seed(0xd1b54a32d192ed03, 50); + let l = &q * &r + UBig::from(1u64); + let num = IBig::from(&l + &r); + let den = l; + + let v = RBig::from_parts(num.clone(), den.clone()); + // gcd(L+R, L) = gcd(L, R) = gcd(Q*R+1, R) = gcd(1, R) = 1, so it is already reduced. + assert_eq!(v.numerator(), &num); + assert_eq!(v.denominator(), &den); +} From 8dd6a9affd5c3c68cc4fec85829a86fd022f538a Mon Sep 17 00:00:00 2001 From: Jacob Zhong Date: Fri, 24 Jul 2026 23:19:17 +0800 Subject: [PATCH 2/5] Drop opaque finding IDs from code/test comments The DASHU-NNN tokens reference an external fuzz-findings corpus that other developers don't have context for. Rephrase the affected comments (two production TODOs and several test docstrings) to describe the limitation or regression in self-contained terms. Co-Authored-By: Claude --- float/src/exp.rs | 17 ++++++++--------- float/src/log.rs | 10 +++++----- integer/tests/convert.rs | 6 +++--- integer/tests/gcd.rs | 10 +++++----- rational/tests/arith_prop.rs | 4 ++-- 5 files changed, 23 insertions(+), 24 deletions(-) diff --git a/float/src/exp.rs b/float/src/exp.rs index 5c3040ef..b9336d24 100644 --- a/float/src/exp.rs +++ b/float/src/exp.rs @@ -109,7 +109,7 @@ impl Context { /// Panics if the precision is unlimited and the exponent is negative. In this case, the exact /// result is likely to have infinite digits. pub fn powi(&self, base: &Repr, exp: IBig) -> FpResult> { - // TODO(DASHU-024): range handling has three known limitations at the exponent extremes: + // Known limitations in range handling at the exponent extremes: // (1) the overflow guard below estimates the result magnitude with an f64 and misclassifies // representable boundaries near 2^63; (2) genuine Overflow/Underflow is unwrapped mode-blindly // to ±inf / signed zero; (3) the negative-exponent reciprocal path can panic. None affects @@ -409,11 +409,10 @@ impl Context { // |floor(x / ln B)| overflows isize — x is astronomically large, so the // result is an infinity (x → +∞) or underflows to the limit (x → −∞). // - // TODO(DASHU-023): this branch discards the rounding mode. For a huge + // Known limitation: this branch discards the rounding mode. For a huge // *negative* x the true exp(x) is positive but below the exponent range, so // directed Up/Away should return the minimum positive value (and exp_m1 the - // value just above -1) rather than +0 / Exact(-1). The mode-blind saturation - // here is a known limitation. + // value just above -1) rather than +0 / Exact(-1). return if input_sign == Sign::Positive { Err(FpError::Overflow(Sign::Positive)) } else if minus_one { @@ -490,11 +489,11 @@ mod tests { #[test] fn test_exp_m1_large_negative_no_oom() { - // Regression for DASHU-026: exp_m1 of a large-magnitude negative input subtracts + // Regression test: exp_m1 of a large-magnitude negative input subtracts // 1 from an astronomically small exp(x), so the aligned subtraction hands // round_fract a sparse sticky tail whose `precision` equals the exponent gap - // (~6.6e18 here). The old debug assert materialized B^precision (2^6.6e18 bits → - // certain OOM); the new check uses log2 bounds and must not allocate. + // (~6.6e18 here). The debug assert in round_fract used to materialize B^precision + // (2^6.6e18 bits → certain OOM); it now uses log2 bounds and must not allocate. // // x = -2^62 is large enough that the exponent gap dwarfs available memory, yet // small enough that floor(x/ln2) still fits isize (≈6.6e18 < 9.2e18), so it does @@ -511,8 +510,8 @@ mod tests { #[test] fn test_exact_results_on_unlimited_precision() { - // Regression for DASHU-022: values carrying precision 0 (unlimited) — produced - // by `try_from(0.0)` and the `FBig::ONE`/`ZERO` constants — must still compute + // Regression test: values carrying precision 0 (unlimited) — produced by + // `try_from(0.0)` and the `FBig::ONE`/`ZERO` constants — must still compute // their exact-result special cases instead of panicking in // assert_limited_precision before reaching the shortcut. type F = FBig; diff --git a/float/src/log.rs b/float/src/log.rs index c0206cb0..325dbddf 100644 --- a/float/src/log.rs +++ b/float/src/log.rs @@ -466,11 +466,11 @@ mod tests { #[test] fn test_log2_powers_of_two() { - // log2(2^k) = k, including log2(1) = 0. This is the headline DASHU-007 case: the - // adapter previously derived directed log2 from f32-precision log2_bounds, which - // gave a small *negative* value for log2(1) and was loose by many ULPs elsewhere. - // The correctly-rounded log2 returns the exact integer value for every power of - // two (log2(1) is tagged Exact via the shortcut; other powers are value-exact). + // log2(2^k) = k, including log2(1) = 0. The only magnitude estimate available + // (log2_bounds) is f32-precision, so deriving directed log2 from it gives a small + // *negative* value for log2(1) and is loose by many ULPs elsewhere; the + // correctly-rounded log2 returns the exact integer value for every power of two + // (log2(1) is tagged Exact via the shortcut; other powers are value-exact). let ctx = Context::::new(53); for k in [0usize, 1, 2, 3, 10, 50, 100, 200, 1000] { let x = Repr::new(IBig::from(1) << k, 0); // 2^k diff --git a/integer/tests/convert.rs b/integer/tests/convert.rs index 88c093a2..1d1775a9 100644 --- a/integer/tests/convert.rs +++ b/integer/tests/convert.rs @@ -434,9 +434,9 @@ fn test_to_f64() { ); } -/// Regression test for DASHU-015: a magnitude at `DoubleWord::MAX` rounds up to the -/// next power of two as f64, and the saturating `f as DoubleWord` round-trip used to -/// clamp back to `DoubleWord::MAX` and report that inexact conversion as `Exact`. +/// Regression test: a magnitude at `DoubleWord::MAX` rounds up to the next power +/// of two as f64, and the saturating `f as DoubleWord` round-trip used to clamp +/// back to `DoubleWord::MAX` and report that inexact conversion as `Exact`. #[cfg(target_pointer_width = "64")] #[test] fn test_to_f64_dword_max_boundary() { diff --git a/integer/tests/gcd.rs b/integer/tests/gcd.rs index 622ac985..8ccbdaa1 100644 --- a/integer/tests/gcd.rs +++ b/integer/tests/gcd.rs @@ -169,11 +169,11 @@ fn test_gcd_ext_0() { #[test] fn test_gcd_large_lopsided_reduction() { - // Regression for DASHU-008: the gcd scratchpad is reserved once from the *initial* - // operand lengths, but each euclidean step's division dispatches on the *current* - // lengths. The pre-fix code under-reserved when the initial operands were similar - // in size, so a later lopsided euclidean step (divisor > 48 words and quotient - // > 32 words, taking the Burninkel-Ziegler + Karatsuba path) aborted with + // Regression test: the gcd scratchpad is reserved once from the *initial* operand + // lengths, but each euclidean step's division dispatches on the *current* lengths. + // Before the fix, two large, similarly-sized operands (zero initial scratch) that + // later reduced through a lopsided euclidean step (divisor > 48 words and quotient + // > 32 words, taking the Burnikel-Ziegler + Karatsuba path) aborted with // "internal error: not enough memory allocated" in both debug and release. use dashu_int::UBig; diff --git a/rational/tests/arith_prop.rs b/rational/tests/arith_prop.rs index 493c5ae7..c5974dde 100644 --- a/rational/tests/arith_prop.rs +++ b/rational/tests/arith_prop.rs @@ -60,8 +60,8 @@ proptest! { } } -/// Regression for DASHU-008: constructing an `RBig` reduces numerator/denominator by -/// their gcd, which must not abort with "not enough memory allocated" when that gcd +/// Regression test: constructing an `RBig` reduces numerator/denominator by their +/// gcd, which must not abort with "not enough memory allocated" when that gcd /// reduction hits a lopsided Burnikel-Ziegler division. The gcd scratchpad is sized /// from the *initial* operand lengths, but each step dispatches on the *current* /// lengths, so two large, similarly-sized operands (zero initial scratch) that later From f9d31aff833aeae25ad280fd69e892101f21c87f Mon Sep 17 00:00:00 2001 From: Jacob Zhong Date: Fri, 24 Jul 2026 23:31:08 +0800 Subject: [PATCH 3/5] Use fixed operands in GCD tests; ban random generators in unit tests MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Replace the xorshift-based `big_from_seed` in the two DASHU-008 regression tests (`test_gcd_large_lopsided_reduction`, `rbig_reduce_large_lopsided`) with fixed all-ones operands sized by word count. The tests still reproduce the exact `memory.rs:150` panic under the old code (verified) but now use deterministic inputs. Add a Code-style rule to AGENTS.md: in-crate tests must use fixed, deterministic inputs — never random/property-test generators — except when testing the `rand` integration itself. Randomized input belongs under `fuzz/`. Co-Authored-By: Claude --- AGENTS.md | 1 + integer/tests/gcd.rs | 21 ++++++--------------- rational/tests/arith_prop.rs | 17 +++++------------ 3 files changed, 12 insertions(+), 27 deletions(-) diff --git a/AGENTS.md b/AGENTS.md index d0133a64..d9439210 100644 --- a/AGENTS.md +++ b/AGENTS.md @@ -50,6 +50,7 @@ Note: always `--exclude dashu-python` when running workspace-wide commands, sinc - When borrowing an algorithm idea from GMP (or any other library), do **not** reference its function names in our docstrings or comments. Describe the algorithm in our own terms and use our own function names (e.g. write `add_mul_dword_same_len_in_place`, never `addmul_2` / `mpn_addmul_2`). External function names must not appear anywhere in the repo. - Tests for a specific algorithm/kernel belong in the same source file as the implementation, as a `#[cfg(test)] mod tests` block at the bottom — not in a separate integration test file under `tests/`. Reserve `tests/` for cross-cutting or public-API tests. - When debugging or writing test assertions, use `{:?}` (or `{:#?}` for the verbose form with digit/bit counts) to inspect arbitrary precision values. The [`Debug`] format prints a compact head‥tail representation (most significant digits `..` least significant digits) instead of dumping the entire number, making it readable even for thousand-digit integers. +- In-crate tests (`#[cfg(test)] mod tests` and `tests/`) must use **fixed, deterministic inputs** — never pseudo-random generators, fuzzing harnesses, or property-test styles that draw random data. Tests must fail or pass identically on every run. The only exception is code that tests the `rand` integration itself. For randomized / exploratory input, use the fuzz targets under `fuzz/` (the proper home for that work), not unit tests. ## Feature flags diff --git a/integer/tests/gcd.rs b/integer/tests/gcd.rs index 8ccbdaa1..03ca7f51 100644 --- a/integer/tests/gcd.rs +++ b/integer/tests/gcd.rs @@ -175,27 +175,18 @@ fn test_gcd_large_lopsided_reduction() { // later reduced through a lopsided euclidean step (divisor > 48 words and quotient // > 32 words, taking the Burnikel-Ziegler + Karatsuba path) aborted with // "internal error: not enough memory allocated" in both debug and release. - use dashu_int::UBig; - - fn big_from_seed(mut seed: u64, words: usize) -> UBig { - let mut v = UBig::from(0u64); - for _ in 0..words { - seed ^= seed << 13; - seed ^= seed >> 7; - seed ^= seed << 17; - v = (v << 64) | UBig::from(seed | 1); // |1 keeps the top word nonzero - } - v - } + use dashu_int::{UBig, Word}; // Build a reduction that starts from similar-sized operands (so the initial division // is "simple" and reserves zero scratch) but reaches a lopsided Burnikel-Ziegler - // step. Let R ~100 words, Q ~50 words, L = Q*R + 1 ~150 words; then + // step. Let R (~100 words), Q (~50 words), L = Q*R + 1 (~150 words); then // gcd(L+R, L): step 1 (L+R)/L has quotient 1 -> next pair (L, R) // step 2 L/R has quotient Q (>32 words) with divisor R (>48 words) // -> Burnikel-Ziegler + Karatsuba, the under-reserved step. - let r = big_from_seed(0x9e3779b97f4a7c15, 100); - let q = big_from_seed(0xd1b54a32d192ed03, 50); + // R and Q are fixed all-ones values chosen only for their size; the path is selected + // by operand length, not value. + let r = (UBig::from(1u8) << (100 * Word::BITS as usize)) - ubig!(1); // ~100 words + let q = (UBig::from(1u8) << (50 * Word::BITS as usize)) - ubig!(1); // ~50 words let l = &q * &r + ubig!(1); let a = &l + &r; // L + R let b = l; // L diff --git a/rational/tests/arith_prop.rs b/rational/tests/arith_prop.rs index c5974dde..acc3bf57 100644 --- a/rational/tests/arith_prop.rs +++ b/rational/tests/arith_prop.rs @@ -68,20 +68,13 @@ proptest! { /// reduce through a wide quotient under-reserved and panicked. #[test] fn rbig_reduce_large_lopsided() { - fn big_from_seed(mut seed: u64, words: usize) -> UBig { - let mut v = UBig::from(0u64); - for _ in 0..words { - seed ^= seed << 13; - seed ^= seed >> 7; - seed ^= seed << 17; - v = (v << 64) | UBig::from(seed | 1); // |1 keeps the top word nonzero - } - v - } + use dashu_int::Word; // L = Q*R + 1, then gcd(L+R, L) reduces through L/R (huge quotient Q, large divisor R). - let r = big_from_seed(0x9e3779b97f4a7c15, 100); - let q = big_from_seed(0xd1b54a32d192ed03, 50); + // R and Q are fixed all-ones values chosen only for their size; the reduction path is + // selected by operand length, not value. + let r = (UBig::from(1u8) << (100 * Word::BITS as usize)) - UBig::from(1u8); // ~100 words + let q = (UBig::from(1u8) << (50 * Word::BITS as usize)) - UBig::from(1u8); // ~50 words let l = &q * &r + UBig::from(1u64); let num = IBig::from(&l + &r); let den = l; From 8b5fa15fcefa663cfea10dc3c8295341b6ac0065 Mon Sep 17 00:00:00 2001 From: Jacob Zhong Date: Fri, 24 Jul 2026 23:36:16 +0800 Subject: [PATCH 4/5] Gate exp_m1 OOM regression test to 64-bit MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit On 32-bit targets `isize` tops out at ~2.1e9, so for x = -2^62 the value floor(x/ln2) ≈ -6.6e18 overflows isize and exp_internal takes its (deferred) overflow branch — returning Exact(-1) — instead of the round_fract path this test targets. A sharp OOM needs an exponent gap large enough that 2^gap exceeds memory, yet still fitting isize; that window only exists on 64-bit. Gate the test accordingly. The underlying round_fract fix is arch-independent. Co-Authored-By: Claude --- float/src/exp.rs | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/float/src/exp.rs b/float/src/exp.rs index b9336d24..ed387158 100644 --- a/float/src/exp.rs +++ b/float/src/exp.rs @@ -487,6 +487,13 @@ mod tests { assert_eq!(m1, -FBig::::ONE); } + // A sharp OOM regression needs an exponent gap large enough that 2^gap exceeds any + // memory (gap ≳ 1e11), yet with floor(x/ln2) still fitting isize so the overflow + // branch is not taken. That window only exists where isize is 64-bit: on 32-bit, + // isize tops out at ~2.1e9 — below any OOM-inducing gap — so the overflow branch + // always intervenes first. The fix itself (log2_bounds in round_fract) is + // arch-independent; only this dedicated sharp test is 64-bit-only. + #[cfg(target_pointer_width = "64")] #[test] fn test_exp_m1_large_negative_no_oom() { // Regression test: exp_m1 of a large-magnitude negative input subtracts @@ -495,9 +502,8 @@ mod tests { // (~6.6e18 here). The debug assert in round_fract used to materialize B^precision // (2^6.6e18 bits → certain OOM); it now uses log2 bounds and must not allocate. // - // x = -2^62 is large enough that the exponent gap dwarfs available memory, yet - // small enough that floor(x/ln2) still fits isize (≈6.6e18 < 9.2e18), so it does - // NOT enter the overflow branch tested above. + // x = -2^62: the gap (~6.6e18) dwarfs memory, yet floor(x/ln2) ≈ -6.6e18 still + // fits 64-bit isize (≈9.2e18), so this does NOT enter the overflow branch. let ctx = Context::::new(2); let x = Repr::new(-(IBig::from(1) << 62), 0); let r = ctx.exp_m1::<2>(&x, None).unwrap().value(); From 3ecb03a9b55f444951a903f2dc37cf6fb6acae02 Mon Sep 17 00:00:00 2001 From: Jacob Zhong Date: Sat, 25 Jul 2026 00:00:33 +0800 Subject: [PATCH 5/5] Restore TODO markers on the deferred powi/exp limitation comments The earlier ID cleanup dropped the TODO tag along with the issue ID. Keep `TODO:` so the limitations still surface in a TODO grep. Co-Authored-By: Claude --- float/src/exp.rs | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/float/src/exp.rs b/float/src/exp.rs index ed387158..54948119 100644 --- a/float/src/exp.rs +++ b/float/src/exp.rs @@ -109,7 +109,7 @@ impl Context { /// Panics if the precision is unlimited and the exponent is negative. In this case, the exact /// result is likely to have infinite digits. pub fn powi(&self, base: &Repr, exp: IBig) -> FpResult> { - // Known limitations in range handling at the exponent extremes: + // TODO: range handling has three known limitations at the exponent extremes: // (1) the overflow guard below estimates the result magnitude with an f64 and misclassifies // representable boundaries near 2^63; (2) genuine Overflow/Underflow is unwrapped mode-blindly // to ±inf / signed zero; (3) the negative-exponent reciprocal path can panic. None affects @@ -409,7 +409,7 @@ impl Context { // |floor(x / ln B)| overflows isize — x is astronomically large, so the // result is an infinity (x → +∞) or underflows to the limit (x → −∞). // - // Known limitation: this branch discards the rounding mode. For a huge + // TODO: this branch discards the rounding mode. For a huge // *negative* x the true exp(x) is positive but below the exponent range, so // directed Up/Away should return the minimum positive value (and exp_m1 the // value just above -1) rather than +0 / Exact(-1).