Bug #14635: Float#round(n) returns a wrong result when n is big

• Status: Closed
• Assignee: -
• Target version: -
• Backport: 2.3: UNKNOWN, 2.4: UNKNOWN, 2.5: UNKNOWN
• Tags: core, patch

[mame (Yusuke Endoh)]

First of all, don't confuse that this is a usual floating-point error issue.

The following looks inconsistent:

3.0e-31           #=> 3.0e-31
3.0e-31.round(31) #=> 3.0000000000000003e-31

What it should be

A Float value is actually a range.

3.0e-31 represents a range of 0.299999999999999959315060e-30 .. 0.300000000000000003105637e-30 (the bounds are approximate). I call this range "A".
3.0000000000000003e-31 represents a range of 0.300000000000000003105637e-30 .. 0.300000000000000046896214e-30. I call this range "B".

x.round(31) should (1) multiple x with 10**31, (2) round it as an integer, and (3) divide it with 10**31.

In this case:

(1) 3.0e-31 * 10**31 is a range of 2.99999999999999959315060 .. 3.00000000000000003105637.
(2) The rounded result is 3, whichever value is chosen from the range above.
(3) 3.0 / 10**31 is within the range "A", not within the range "B", so the result should be 3.0e-31, not 3.0000000000000003e-31.

How the bug occurs

The reason why 3.0e-31.round(31) returns 3.0000000000000003e-31, is the implementation issue of Float#round. It does the following:

(1) f = pow(10, b)
(2) x = round(x * f) as an integer
(3) return x / f

However, a double variable f cannot represent pow(10, 31) precisely. In other words, the 10**31 must be handled as an integer, but the implementation handles it as an inexact floating-point value. This is the issue.

How to fix

The issue is simple, but it might be very difficult to fix. strtod handles a string "3.0e-31" correctly. So, by doing the same as strtod, this issue would be fixed. However, the strtod implementation looks very difficult, at least to me. Contribution from mathematician is welcome.
(Honestly, I don't want to see such a complication in the source code. Another simpler approach would be more preferable.)

References

This issue has been already reported in #5273 by marcandre. But the status of the ticket looks unclear; I cannot see how many issues remains. So, I created this ticket for just one bug that I could confirm.

[mame (Yusuke Endoh)]

• Related to Bug #5273: Float#round returns the wrong floats for higher precision added

[mame (Yusuke Endoh)]

• Subject changed from Float#round sometimes returns a wrong result to Float#round(n) returns a wrong result when n is big

[mame (Yusuke Endoh)]

I've found a much simpler solution: when n is big, it should first translate the float to a rational, then call Rational#round, and finally translate the resulting rational to a float. It is slow, only when n >= 23 for Float#round(n).

Currently, Rational#to_f has the same inaccuracy issue, which can be fixed by #14637. The following patch includes the hunk for #14637.

diff --git a/bignum.c b/bignum.c
index b4c7560034..fd5f385cac 100644
--- a/bignum.c
+++ b/bignum.c
@@ -6178,9 +6178,7 @@ rb_big_fdiv_double(VALUE x, VALUE y)
 	    return big_fdiv_int(x, rb_int2big(FIX2LONG(y)));
     }
     else if (RB_BIGNUM_TYPE_P(y)) {
-	dy = rb_big2dbl(y);
-	if (isinf(dx) || isinf(dy))
-	    return big_fdiv_int(x, y);
+	return big_fdiv_int(x, y);
     }
     else if (RB_FLOAT_TYPE_P(y)) {
 	dy = RFLOAT_VALUE(y);
diff --git a/internal.h b/internal.h
index 9b6a213151..0bf20b19b0 100644
--- a/internal.h
+++ b/internal.h
@@ -1689,6 +1689,7 @@ VALUE rb_cstr_to_rat(const char *, int);
 VALUE rb_rational_abs(VALUE self);
 VALUE rb_rational_cmp(VALUE self, VALUE other);
 VALUE rb_numeric_quo(VALUE x, VALUE y);
+VALUE rb_flo_round_by_rational(int argc, VALUE *argv, VALUE num);
 
 /* re.c */
 VALUE rb_reg_compile(VALUE str, int options, const char *sourcefile, int sourceline);
diff --git a/numeric.c b/numeric.c
index 01856c7f20..15b27e9132 100644
--- a/numeric.c
+++ b/numeric.c
@@ -2239,6 +2239,10 @@ flo_round(int argc, VALUE *argv, VALUE num)
 	frexp(number, &binexp);
 	if (float_round_overflow(ndigits, binexp)) return num;
 	if (float_round_underflow(ndigits, binexp)) return DBL2NUM(0);
+	if (ndigits > DBL_MANT_DIG * log(2.0) / log(5.0)) {
+	    /* In this case, pow(10, ndigits) cannot be accurate. */
+	    return rb_flo_round_by_rational(argc, argv, num);
+	}
 	f = pow(10, ndigits);
 	x = ROUND_CALL(mode, round, (number, f));
 	return DBL2NUM(x / f);
diff --git a/rational.c b/rational.c
index d88f50f886..01bb88d1ae 100644
--- a/rational.c
+++ b/rational.c
@@ -1533,6 +1533,13 @@ nurat_round_n(int argc, VALUE *argv, VALUE self)
     return f_round_common(argc, argv, self, round_func);
 }
 
+static VALUE float_to_r(VALUE self);
+VALUE
+rb_flo_round_by_rational(int argc, VALUE *argv, VALUE num)
+{
+    return nurat_to_f(nurat_round_n(argc, argv, float_to_r(num)));
+}
+
 static double
 nurat_to_double(VALUE self)
 {
@@ -2016,7 +2023,6 @@ integer_denominator(VALUE self)
     return INT2FIX(1);
 }
 
-static VALUE float_to_r(VALUE self);
 /*
  * call-seq:
  *    flo.numerator  ->  integer

[marcandre (Marc-Andre Lafortune)]

• Has duplicate Bug #17183: Float round working weirdly added

[marcandre (Marc-Andre Lafortune)]

• Has duplicate deleted (Bug #17183: Float round working weirdly)

[marcandre (Marc-Andre Lafortune)]

• Related to Bug #17183: Float round working weirdly added

[jeremyevans0 (Jeremy Evans)]

This issue still exists in the master branch. I've submitted a pull request that uses a modified version of @mame's patch to fix this issue: https://github.com/ruby/ruby/pull/4682

• Status changed from Feedback to Closed

Applied in changeset git|d16b68cb2204eeb5af8bd39149202b630374c67f.

---

Use Rational for Float#round with ndigits > 14

ndigits higher than 14 can result in values that are slightly too
large due to floating point limitations. Converting to rational
for the calculation and then back to float fixes these issues.

Fixes [Bug #14635]
Fixes [Bug #17183]

Co-authored by: Yusuke Endoh [e-mail]