From c2910c06f3ff1c7e26d83240fa4137b2ce6478a3 Mon Sep 17 00:00:00 2001 From: Changjoon Date: Sun, 3 May 2026 00:16:27 +0900 Subject: [PATCH] Round float at the decimal level to match CPython's _Py_dg_dtoa (#7761) * Round float at the decimal level to match CPython's _Py_dg_dtoa CPython's `float.__round__` (Objects/floatobject.c) routes through `_Py_dg_dtoa` and rounds at the decimal level. The previous `round_float_digits` multiplied by 10**ndigits and rounded at the IEEE 754 binary level, which diverges for values that aren't exactly representable. For example, 2.675 stores as 2.67499...; dtoa correctly rounds it down to 2.67, but `(2.675 * 100.0).round() / 100.0` lands on 2.68 because the multiplication produces a phantom 267.5 tie that round-half-to-even snaps up. Rust's `{:.*}` float formatting uses dtoa-style algorithms (Grisu3 + Dragon4 fallback) and matches CPython's `_Py_dg_dtoa` byte-for-byte. Replace the multiply-then-round step with `format!` + `parse` for ndigits >= 0. The ndigits < 0 path is unchanged because dividing typical inputs by 10**|ndigits| produces genuine ties rather than synthesizing them. Verified byte-identical with CPython 3.14.4 over a 108-case random fuzz plus targeted half-tie probes. Unmasks `test_float.RoundTestCase.test_matches_float_format` and `test_previous_round_bugs`. * Use #[expect] with reason for float_cmp suppression Co-authored-by: ShaharNaveh <50263213+ShaharNaveh@users.noreply.github.com> --------- Co-authored-by: ShaharNaveh <50263213+ShaharNaveh@users.noreply.github.com> --- Lib/test/test_float.py | 2 - crates/common/src/float_ops.rs | 103 ++++++++++++-------------- extra_tests/snippets/builtin_round.py | 50 +++++++++++++ 3 files changed, 96 insertions(+), 59 deletions(-) diff --git a/Lib/test/test_float.py b/Lib/test/test_float.py index 53f454d18..0938c89cb 100644 --- a/Lib/test/test_float.py +++ b/Lib/test/test_float.py @@ -965,7 +965,6 @@ class RoundTestCase(unittest.TestCase, FloatsAreIdenticalMixin): self.assertRaises(OverflowError, round, 1.6e308, -308) self.assertRaises(OverflowError, round, -1.7e308, -308) - @unittest.expectedFailure # TODO: RUSTPYTHON; AssertionError: 56294995342131.51 != 56294995342131.5 @unittest.skipUnless(getattr(sys, 'float_repr_style', '') == 'short', "applies only when using short float repr style") def test_previous_round_bugs(self): @@ -984,7 +983,6 @@ class RoundTestCase(unittest.TestCase, FloatsAreIdenticalMixin): self.assertEqual(round(85.0, -1), 80.0) self.assertEqual(round(95.0, -1), 100.0) - @unittest.expectedFailure # TODO: RUSTPYTHON; AssertionError: 0.01 != 0.0 @unittest.skipUnless(getattr(sys, 'float_repr_style', '') == 'short', "applies only when using short float repr style") def test_matches_float_format(self): diff --git a/crates/common/src/float_ops.rs b/crates/common/src/float_ops.rs index c63da1077..c5e2e6064 100644 --- a/crates/common/src/float_ops.rs +++ b/crates/common/src/float_ops.rs @@ -1,6 +1,6 @@ use core::f64; use malachite_bigint::{BigInt, ToBigInt}; -use num_traits::{Float, Signed, ToPrimitive, Zero}; +use num_traits::{Signed, ToPrimitive}; pub const fn decompose_float(value: f64) -> (f64, i32) { if 0.0 == value { @@ -208,64 +208,53 @@ pub fn ulp(x: f64) -> f64 { } pub fn round_float_digits(x: f64, ndigits: i32) -> Option { - let float = if ndigits.is_zero() { - let fract = x.fract(); - if (fract.abs() - 0.5).abs() < f64::EPSILON { - if x.trunc() % 2.0 == 0.0 { - x - fract - } else { - x + fract - } - } else { - x.round() - } + // Mirror CPython's `float.__round__` (Objects/floatobject.c), which uses + // `_Py_dg_dtoa` to round at the decimal level. Multiplying by 10**ndigits + // and rounding at the IEEE 754 binary level diverges for values that + // aren't exactly representable: 2.675 stores as 2.67499..., which dtoa + // correctly rounds down to 2.67, but `(2.675 * 100.0).round() / 100.0` + // lands on 2.68 because the multiplication produces a phantom 267.5 tie. + // Rust's `{:.*}` float formatting uses dtoa-style algorithms and matches + // CPython's `_Py_dg_dtoa` byte-for-byte. + if !x.is_finite() { + return Some(x); + } + + const NDIGITS_MAX: i32 = + ((f64::MANTISSA_DIGITS as i32 - f64::MIN_EXP) as f64 * f64::consts::LOG10_2) as i32; + const NDIGITS_MIN: i32 = -(((f64::MAX_EXP + 1) as f64 * f64::consts::LOG10_2) as i32); + + if ndigits > NDIGITS_MAX { + return Some(x); + } + if ndigits < NDIGITS_MIN { + return Some(0.0f64.copysign(x)); + } + + let result: f64 = if ndigits >= 0 { + let s = format!("{:.*}", ndigits as usize, x); + s.parse().ok()? } else { - const NDIGITS_MAX: i32 = - ((f64::MANTISSA_DIGITS as i32 - f64::MIN_EXP) as f64 * f64::consts::LOG10_2) as i32; - const NDIGITS_MIN: i32 = -(((f64::MAX_EXP + 1) as f64 * f64::consts::LOG10_2) as i32); - if ndigits > NDIGITS_MAX { - x - } else if ndigits < NDIGITS_MIN { - 0.0f64.copysign(x) + // ndigits < 0: divide-then-round avoids the phantom-tie problem + // because dividing typical inputs by 10**|ndigits| produces genuine + // half-integer ties rather than synthesizing them. + let pow1 = 10.0f64.powi(-ndigits); + let y = x / pow1; + let z = y.round(); + #[expect( + clippy::float_cmp, + reason = "exact half-tie detection for banker's rounding correction" + )] + let z = if (y - z).abs() == 0.5 { + 2.0 * (y / 2.0).round() } else { - let (y, pow1, pow2) = if ndigits >= 0 { - // according to cpython: pow1 and pow2 are each safe from overflow, but - // pow1*pow2 ~= pow(10.0, ndigits) might overflow - let (pow1, pow2) = if ndigits > 22 { - (10.0.powf((ndigits - 22) as f64), 1e22) - } else { - (10.0.powf(ndigits as f64), 1.0) - }; - let y = (x * pow1) * pow2; - if !y.is_finite() { - return Some(x); - } - (y, pow1, Some(pow2)) - } else { - let pow1 = 10.0.powf((-ndigits) as f64); - (x / pow1, pow1, None) - }; - let z = y.round(); - #[allow(clippy::float_cmp)] - let z = if (y - z).abs() == 0.5 { - 2.0 * (y / 2.0).round() - } else { - z - }; - let z = if let Some(pow2) = pow2 { - // ndigits >= 0 - (z / pow2) / pow1 - } else { - z * pow1 - }; - - if !z.is_finite() { - // overflow - return None; - } - z - } + }; + z * pow1 }; - Some(float) + + if !result.is_finite() { + return None; + } + Some(result) } diff --git a/extra_tests/snippets/builtin_round.py b/extra_tests/snippets/builtin_round.py index 99c4ed1d2..e94d97542 100644 --- a/extra_tests/snippets/builtin_round.py +++ b/extra_tests/snippets/builtin_round.py @@ -1,3 +1,5 @@ +import math + from testutils import assert_raises assert round(1.2) == 1 @@ -43,3 +45,51 @@ class X: assert round(X(), 1) == 1.1 assert round(X(), None) == 1.1 assert round(X()) == 1.1 + + +# Banker's rounding at the decimal level: values like 2.675 store as +# 2.67499... in IEEE 754, so a multiply-then-round implementation creates +# a phantom 267.5 tie that round-half-to-even snaps up to 2.68. CPython's +# `_Py_dg_dtoa` rounds at the decimal level and returns 2.67. +assert round(2.675, 2) == 2.67 +assert round(2.685, 2) == 2.69 +assert round(-2.675, 2) == -2.67 +assert round(0.05, 1) == 0.1 +assert round(0.15, 1) == 0.1 +assert round(0.35, 1) == 0.3 +assert round(0.45, 1) == 0.5 +assert round(0.65, 1) == 0.7 +assert round(0.85, 1) == 0.8 +assert round(0.95, 1) == 0.9 +assert round(0.645, 2) == 0.65 +assert round(0.665, 2) == 0.67 +assert round(0.685, 2) == 0.69 +assert round(0.695, 2) == 0.69 +assert round(1.685, 2) == 1.69 +assert round(3.745, 2) == 3.75 + +# Exact-halfway ties at integer level use banker's (round-half-to-even). +assert round(0.5, 0) == 0.0 +assert round(1.5, 0) == 2.0 +assert round(2.5, 0) == 2.0 + +# Negative ndigits uses divide-then-round; banker's still applies. +assert round(1235, -1) == 1240 +assert round(1245, -1) == 1240 +assert round(1234.5, -1) == 1230.0 +assert round(150, -2) == 200 +assert round(250, -2) == 200 + +# NaN and infinities round to themselves with no error when ndigits is given. +assert math.isnan(round(float("nan"), 2)) +assert round(float("inf"), 2) == float("inf") +assert round(float("-inf"), 2) == float("-inf") + +# Signed zero is preserved through rounding. +assert math.copysign(1.0, round(-0.0, 2)) == -1.0 +assert math.copysign(1.0, round(0.0, 2)) == 1.0 + +# Out-of-range ndigits short-circuits without overflow. +assert round(1.0, 1000) == 1.0 +assert round(1.0, -1000) == 0.0 +assert round(1.7976931348623157e308, 0) == 1.7976931348623157e308