diff --git a/src/uucore/src/lib/features/parser/num_parser.rs b/src/uucore/src/lib/features/parser/num_parser.rs index 726e89f67..33acdbee7 100644 --- a/src/uucore/src/lib/features/parser/num_parser.rs +++ b/src/uucore/src/lib/features/parser/num_parser.rs @@ -8,7 +8,7 @@ // spell-checker:ignore powf copysign prec inity infinit bigdecimal extendedbigdecimal biguint underflowed use bigdecimal::{ - BigDecimal, + BigDecimal, Context, num_bigint::{BigInt, BigUint, Sign}, }; use num_traits::Signed; @@ -286,6 +286,32 @@ fn make_error<'a>(overflow: bool, negative: bool) -> ExtendedParserError<'a, Ext } } +/// Compute bd**exp using exponentiation by squaring algorithm, while maintaining the +/// precision specified in ctx (the number of digits would otherwise explode). +// TODO: We do lose a little bit of precision, and the last digits may not be correct. +// TODO: Upstream this to bigdecimal-rs. +fn pow_with_context(bd: BigDecimal, exp: u32, ctx: &bigdecimal::Context) -> BigDecimal { + if exp == 0 { + return 1.into(); + } + + fn trim_precision(bd: BigDecimal, ctx: &bigdecimal::Context) -> BigDecimal { + if bd.digits() > ctx.precision().get() { + bd.with_precision_round(ctx.precision(), ctx.rounding_mode()) + } else { + bd + } + } + + let bd = trim_precision(bd, ctx); + let ret = if exp % 2 == 0 { + pow_with_context(bd.square(), exp / 2, ctx) + } else { + &bd * pow_with_context(bd.square(), (exp - 1) / 2, ctx) + }; + trim_precision(ret, ctx) +} + // Construct an ExtendedBigDecimal based on parsed data fn construct_extended_big_decimal<'a>( digits: BigUint, @@ -339,13 +365,14 @@ fn construct_extended_big_decimal<'a>( // Confusingly, exponent is in base 2 for hex floating point numbers. // Note: We cannot overflow/underflow BigDecimal here, as we will not be able to reach the // maximum/minimum scale (i64 range). - let pow2 = BigDecimal::from_bigint(BigInt::from(2).pow(abs_exponent.to_u32().unwrap()), 0); - - if !exponent.is_negative() { - bd * pow2 + let base: BigDecimal = if !exponent.is_negative() { + 2.into() } else { - bd / pow2 - } + BigDecimal::from(2).inverse() + }; + let pow2 = pow_with_context(base, abs_exponent.to_u32().unwrap(), &Context::default()); + + bd * pow2 } else { // scale != 0, which means that integral_only is not set, so only base 10 and 16 are allowed. unreachable!(); @@ -771,6 +798,26 @@ mod tests { ExtendedBigDecimal::extended_parse("0xf.fffffffffffffffffffff") ); + // Test very large exponents (they used to take forever as we kept all digits in the past) + // Wolfram Alpha can get us (close to?) these values with a bit of log trickery: + // 2**3000000000 = 10**log_10(2**3000000000) = 10**(3000000000 * log_10(2)) + // TODO: We do lose a little bit of precision, and the last digits are not be correct. + assert_eq!( + Ok(ExtendedBigDecimal::BigDecimal( + // Wolfram Alpha says 9.8162042336235053508313854078782835648991393286913072670026492205522618203568834202759669215027003865... × 10^903089986 + BigDecimal::from_str("9.816204233623505350831385407878283564899139328691307267002649220552261820356883420275966921514831318e+903089986").unwrap() + )), + ExtendedBigDecimal::extended_parse("0x1p3000000000") + ); + assert_eq!( + Ok(ExtendedBigDecimal::BigDecimal( + // Wolfram Alpha says 1.3492131462369983551036088935544888715959511045742395978049631768570509541390540646442193112226520316... × 10^-9030900 + BigDecimal::from_str("1.349213146236998355103608893554488871595951104574239597804963176857050954139054064644219311222656999e-9030900").unwrap() + )), + // Couldn't get a answer from Wolfram Alpha for smaller negative exponents + ExtendedBigDecimal::extended_parse("0x1p-30000000") + ); + // ExtendedBigDecimal overflow/underflow assert!(matches!( ExtendedBigDecimal::extended_parse(&format!("0x1p{}", u32::MAX as u64 + 1)),