Correct float rounding algorithm

Lib/test/test_float.py

@@ -832,8 +832,6 @@ class RoundTestCase(unittest.TestCase):
         self.assertRaises(TypeError, round, NAN, "ceci n'est pas un integer")
         self.assertRaises(TypeError, round, -0.0, 1j)
 
-    # TODO: RUSTPYTHON
-    @unittest.expectedFailure
     def test_large_n(self):
         for n in [324, 325, 400, 2**31-1, 2**31, 2**32, 2**100]:
             self.assertEqual(round(123.456, n), 123.456)
@@ -846,8 +844,6 @@ class RoundTestCase(unittest.TestCase):
         self.assertEqual(round(1e150, 309), 1e150)
         self.assertEqual(round(1.4e-315, 315), 1e-315)
 
-    # TODO: RUSTPYTHON
-    @unittest.expectedFailure
     def test_small_n(self):
         for n in [-308, -309, -400, 1-2**31, -2**31, -2**31-1, -2**100]:
             self.assertEqual(round(123.456, n), 0.0)
@@ -855,8 +851,6 @@ class RoundTestCase(unittest.TestCase):
             self.assertEqual(round(1e300, n), 0.0)
             self.assertEqual(round(1e-320, n), 0.0)
 
-    # TODO: RUSTPYTHON
-    @unittest.expectedFailure
     def test_overflow(self):
         self.assertRaises(OverflowError, round, 1.6e308, -308)
         self.assertRaises(OverflowError, round, -1.7e308, -308)

common/src/float_ops.rs

@@ -1,5 +1,6 @@
 use num_bigint::{BigInt, ToBigInt};
 use num_traits::{Float, Signed, ToPrimitive, Zero};
+use std::f64;
 
 pub fn ufrexp(value: f64) -> (f64, i32) {
     if 0.0 == value {
@@ -366,6 +367,69 @@ pub fn ulp(x: f64) -> f64 {
     }
 }
 
+pub fn round_float_digits(x: f64, ndigits: i32) -> Option<f64> {
+    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()
+        }
+    } 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)
+        } 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
+        }
+    };
+    Some(float)
+}
+
 #[test]
 fn test_to_hex() {
     use rand::Rng;

vm/src/builtins/float.rs

@@ -1,8 +1,7 @@
 use num_bigint::{BigInt, ToBigInt};
 use num_complex::Complex64;
 use num_rational::Ratio;
-use num_traits::{pow, Signed, ToPrimitive, Zero};
-use std::f64;
+use num_traits::{Signed, ToPrimitive, Zero};
 
 use super::bytes::PyBytes;
 use super::int::{self, PyInt, PyIntRef};
@@ -395,43 +394,13 @@ impl PyFloat {
         let ndigits = ndigits.flatten();
         let value = if let Some(ndigits) = ndigits {
             let ndigits = ndigits.borrow_value();
-            let float = if ndigits.is_zero() {
-                let fract = self.value.fract();
-                if (fract.abs() - 0.5).abs() < f64::EPSILON {
-                    if self.value.trunc() % 2.0 == 0.0 {
-                        self.value - fract
-                    } else {
-                        self.value + fract
-                    }
-                } else {
-                    self.value.round()
-                }
-            } else {
-                let ndigits = match ndigits.to_isize() {
-                    Some(n) => n,
-                    None if ndigits.is_positive() => isize::MAX,
-                    None => isize::MIN,
-                };
-                const NDIGITS_MAX: isize = ((f64::MANTISSA_DIGITS as i32 - f64::MIN_EXP) as f64
-                    * f64::consts::LOG10_2) as isize;
-                const NDIGITS_MIN: isize =
-                    -(((f64::MAX_EXP + 1) as f64 * f64::consts::LOG10_2) as isize);
-                if ndigits > NDIGITS_MAX {
-                    self.value
-                } else if ndigits > NDIGITS_MIN {
-                    0.0f64.copysign(self.value)
-                } else if ndigits >= 0 {
-                    (self.value * pow(10.0, ndigits as usize)).round() / pow(10.0, ndigits as usize)
-                } else {
-                    let result = (self.value / pow(10.0, (-ndigits) as usize)).round()
-                        * pow(10.0, (-ndigits) as usize);
-                    if result.is_nan() {
-                        0.0
-                    } else {
-                        result
-                    }
-                }
+            let ndigits = match ndigits.to_i32() {
+                Some(n) => n,
+                None if ndigits.is_positive() => i32::MAX,
+                None => i32::MIN,
             };
+            let float = float_ops::round_float_digits(self.value, ndigits)
+                .ok_or_else(|| vm.new_overflow_error("overflow ocurred during round".to_owned()))?;
             vm.ctx.new_float(float)
         } else {
             let fract = self.value.fract();