Match CPython's _float_div_mod, fixing divmod and % zero-handling (#7745)

float_ops::divmod, mod_, and floordiv each carried their own conversion
from Rust's dividend-sign `%` to CPython's divisor-sign convention. Both
divmod and mod_ mishandled the zero-remainder case where the dividend
is a non-zero exact multiple of the divisor (e.g. divmod(6.0, -3.0),
6.0 % -3.0): the sign-correction branch fired on a zero remainder and
produced (-3.0, -3.0) and -3.0 respectively, violating the magnitude
invariant 0 <= abs(r) < abs(b). divmod also leaked the wrong signed-
zero quotient when the true quotient was zero (divmod(-1.0, -2.0)
returned (-0.0, -1.0) instead of (+0.0, -1.0)).

These are independent bugs in two functions, but both come from the
same root cause: zero-remainder needs a separate path from the sign-
correction branch.

Mirror CPython's `_float_div_mod` (Objects/floatobject.c) by making
divmod the canonical implementation and turning mod_ and floordiv into
thin wrappers. divmod(a, b) == (a // b, a % b) now holds by
construction.

Closes #7722

crates/common/src/float_ops.rs

@@ -68,38 +68,42 @@ pub const fn div(v1: f64, v2: f64) -> Option<f64> {
 }
 
 pub fn mod_(v1: f64, v2: f64) -> Option<f64> {
-    if v2 != 0.0 {
-        let val = v1 % v2;
-        match (v1.signum() as i32, v2.signum() as i32) {
-            (1, 1) | (-1, -1) => Some(val),
-            _ if (v1 == 0.0) || (v1.abs() == v2.abs()) => Some(val.copysign(v2)),
-            _ => Some((val + v2).copysign(v2)),
-        }
-    } else {
-        None
-    }
+    divmod(v1, v2).map(|(_, m)| m)
 }
 
 pub fn floordiv(v1: f64, v2: f64) -> Option<f64> {
-    if v2 != 0.0 {
-        Some((v1 / v2).floor())
-    } else {
-        None
-    }
+    divmod(v1, v2).map(|(d, _)| d)
 }
 
+// Canonical (floordiv, mod) for floats matching CPython's _float_div_mod
+// (Objects/floatobject.c). `mod_` and `floordiv` delegate here so that
+// `divmod(a, b) == (a // b, a % b)` holds by construction.
 pub fn divmod(v1: f64, v2: f64) -> Option<(f64, f64)> {
-    if v2 != 0.0 {
-        let mut m = v1 % v2;
-        let mut d = (v1 - m) / v2;
+    if v2 == 0.0 {
+        return None;
+    }
+    let mut m = v1 % v2;
+    let mut d = (v1 - m) / v2;
+    if m != 0.0 {
+        // Non-zero remainder must have the sign of the divisor.
         if v2.is_sign_negative() != m.is_sign_negative() {
             m += v2;
             d -= 1.0;
         }
-        Some((d, m))
     } else {
-        None
+        // Zero remainder: sign matches divisor (IEEE 754 / CPython contract).
+        m = (0.0_f64).copysign(v2);
     }
+    let d = if d != 0.0 {
+        let f = d.floor();
+        // Snap up if (v1 - m) / v2 undershot the true integer quotient by
+        // more than half an ULP (mirrors CPython's `if (div - *floordiv > 0.5)`).
+        if d - f > 0.5 { f + 1.0 } else { f }
+    } else {
+        // Zero quotient: take the sign of the true quotient v1 / v2.
+        (0.0_f64).copysign(v1 / v2)
+    };
+    Some((d, m))
 }
 
 // nextafter algorithm based off of https://gitlab.com/bronsonbdevost/next_afterf

extra_tests/snippets/builtin_divmod.py

@@ -1,3 +1,5 @@
+import math
+
 from testutils import assert_raises
 
 assert divmod(11, 3) == (3, 2)
@@ -7,3 +9,74 @@ assert divmod(-86340, 86400) == (-1, 60)
 
 assert_raises(ZeroDivisionError, divmod, 5, 0, _msg="divmod by zero")
 assert_raises(ZeroDivisionError, divmod, 5.0, 0.0, _msg="divmod by zero")
+
+
+def signbit(x):
+    return math.copysign(1.0, x) < 0
+
+
+# Zero remainder with opposite-sign divisor — quotient and remainder must
+# both be zero with the divisor's sign, not propagate through the
+# sign-correction branch.
+q, r = divmod(0.0, -1.0)
+assert q == 0.0 and signbit(q)
+assert r == 0.0 and signbit(r)
+
+q, r = divmod(6.0, -3.0)
+assert q == -2.0
+assert r == 0.0 and signbit(r)
+
+q, r = divmod(-100.0, 10.0)
+assert q == -10.0
+assert r == 0.0 and not signbit(r)
+
+# Zero quotient — sign matches the true quotient v1 / v2, not the sign that
+# leaks from the (v1 - m) / v2 intermediate calculation.
+q, r = divmod(-1.0, -2.0)
+assert q == 0.0 and not signbit(q)
+assert r == -1.0
+
+q, r = divmod(-0.0, 1.0)
+assert q == 0.0 and signbit(q)
+assert r == 0.0 and not signbit(r)
+
+# Spec invariant: divmod(a, b) == (a // b, a % b), including signed zero.
+for a, b in [
+    (0.0, -1.0),
+    (6.0, -3.0),
+    (-6.0, 3.0),
+    (100.0, -10.0),
+    (-100.0, 10.0),
+    (-1.0, -2.0),
+    (-0.0, 1.0),
+    (-0.0, -1.0),
+    (7.0, 3.0),
+    (-7.0, 3.0),
+    (7.0, -3.0),
+    (-7.0, -3.0),
+    (3.7, 1.5),
+    (-3.7, 1.5),
+]:
+    dm = divmod(a, b)
+    assert dm[0] == a // b
+    assert dm[1] == a % b
+    assert signbit(dm[0]) == signbit(a // b)
+    assert signbit(dm[1]) == signbit(a % b)
+
+# Spec invariants for float divmod:
+#   q * b + r == a, r == 0 or sign(r) == sign(b), 0 <= abs(r) < abs(b).
+for a, b in [
+    (7.0, 3.0),
+    (7.0, -3.0),
+    (-7.0, 3.0),
+    (-7.0, -3.0),
+    (6.0, -3.0),
+    (100.0, -10.0),
+    (3.7, 1.5),
+    (5.5, 2.0),
+    (-5.5, 2.0),
+]:
+    q, r = divmod(a, b)
+    assert q * b + r == a
+    assert r == 0.0 or (r < 0.0) == (b < 0.0)
+    assert abs(r) < abs(b)