Fix rounding in scm_i_divide2double for negative arguments.

* libguile/numbers.c (INUM_LOSSLESSLY_CONVERTIBLE_TO_DOUBLE):
  New macro.
  (scm_i_divide2double): Use INUM_LOSSLESSLY_CONVERTIBLE_TO_DOUBLE to
  determine if our fast path is safe.  Previously, negative arguments
  were not checked properly.

* test-suite/tests/numbers.test (exact->inexact): Add tests.

diff --git a/libguile/numbers.c b/libguile/numbers.c
index 1f4b9a8..0abcb25 100644 (file)
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -100,6 +100,13 @@ typedef scm_t_signed_bits scm_t_inum;
 #define DOUBLE_IS_POSITIVE_INFINITY(x) (isinf(x) && ((x) > 0))
 #define DOUBLE_IS_NEGATIVE_INFINITY(x) (isinf(x) && ((x) < 0))
 
+/* Test an inum to see if it can be converted to a double without loss
+   of precision.  Note that this will sometimes return 0 even when 1
+   could have been returned, e.g. for large powers of 2.  It is designed
+   to be a fast check to optimize common cases. */
+#define INUM_LOSSLESSLY_CONVERTIBLE_TO_DOUBLE(n)                        \
+  (SCM_I_FIXNUM_BIT-1 <= DBL_MANT_DIG                                   \
+   || ((n) ^ ((n) >> (SCM_I_FIXNUM_BIT-1))) < (1L << DBL_MANT_DIG))
 
 #if ! HAVE_DECL_MPZ_INITS
 
@@ -506,10 +513,10 @@ scm_i_divide2double (SCM n, SCM d)
 
   if (SCM_LIKELY (SCM_I_INUMP (d)))
     {
-      if (SCM_LIKELY (SCM_I_INUMP (n)
-                      && (SCM_I_FIXNUM_BIT-1 <= DBL_MANT_DIG
-                          || (SCM_I_INUM (n) < (1L << DBL_MANT_DIG)
-                              && SCM_I_INUM (d) < (1L << DBL_MANT_DIG)))))
+      if (SCM_LIKELY
+          (SCM_I_INUMP (n)
+           && INUM_LOSSLESSLY_CONVERTIBLE_TO_DOUBLE (SCM_I_INUM (n))
+           && INUM_LOSSLESSLY_CONVERTIBLE_TO_DOUBLE (SCM_I_INUM (d))))
         /* If both N and D can be losslessly converted to doubles, then
            we can rely on IEEE floating point to do proper rounding much
            faster than we can. */

diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index eca4536..3d25e6a 100644 (file)
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -4021,6 +4021,19 @@
     (let ((big (ash 1 4096)))
       (= 1.0 (exact->inexact (/ (1+ big) big)))))
 
+  ;; In guile 2.0.9, 'exact->inexact' guaranteed proper rounding when
+  ;; applied to non-negative fractions, but on 64-bit systems would
+  ;; sometimes double-round when applied to negative fractions,
+  ;; specifically when the numerator was a fixnum not exactly
+  ;; representable as a double.
+  (with-test-prefix "frac inum/inum, numerator not exactly representable as a double"
+    (let ((n (+ 1 (expt 2 dbl-mant-dig))))
+      (for-each (lambda (d)
+                  (test (/ n d)
+                        (/ n d)
+                        (exact->inexact (/ n d))))
+                '(3 5 6 7 9 11 13 17 19 23 0.0 -0.0 +nan.0 +inf.0 -inf.0))))
+
   (test "round up to odd"
         ;; =====================================================
         ;; 11111111111111111111111111111111111111111111111111000101 ->