[commit: ghc] master: Fix integer overflow when encoding doubles (Trac #15271) (8e1a559)

[git at git.haskell.org]

Repository : ssh://git@git.haskell.org/ghc

On branch  : master
Link       : http://ghc.haskell.org/trac/ghc/changeset/8e1a5593f940b9d2882a5e81abf028e7b99aff3a/ghc

>---------------------------------------------------------------

commit 8e1a5593f940b9d2882a5e81abf028e7b99aff3a
Author: Fangyi Zhou <fangyi.zhou15 at imperial.ac.uk>
Date:   Sun Oct 28 12:28:53 2018 -0400

    Fix integer overflow when encoding doubles (Trac #15271)
    
    Summary:
    Ticket #15271 reports a case where 1e1000000000 is incorrectly
    converted to 0.0. After some investigation, I discovered the number is
    converted to rational correctly, but converting the ratio into a double
    introduced an error.
    
    Tracking down to how the conversion is done, I found the rts float
    implementation uses `ldexp`, whose signature is
    `double ldexp (double x, int exp);`
    The callsite passes an `I_` to the second argument, which is
    platform-dependent. On machines where `I_` is 64 bits and `int` is 32 bits, we
    observe integer overflow behaviour.
    
    Here is a mapping from rational to exponent with observations
    1e646457008  -> 2147483645 (result = infinity, positive in int32)
    1e646457009  -> 2147483648 (result = 0.0, overflow to negative in int32)
    1e1000000000 -> 3321928042 (result = infinity, overflow to positive in int32)
    1e1555550000 -> 5167425196 (result = 0.0, overflow to negative in int32)
    
    We fix this issue by comparing STG_INT_MIN/MAX and INT_MIN/MAX and bound the
    value appropriately.
    
    Test Plan: New test cases
    
    Reviewers: bgamari, erikd, simonmar
    
    Reviewed By: bgamari
    
    Subscribers: rwbarton, carter
    
    GHC Trac Issues: #15271
    
    Differential Revision: https://phabricator.haskell.org/D5271


>---------------------------------------------------------------

8e1a5593f940b9d2882a5e81abf028e7b99aff3a
 rts/StgPrimFloat.c                                 | 26 ++++++++++++++++++----
 testsuite/tests/numeric/should_run/T15271.hs       |  5 +++++
 .../tests/numeric/should_run/T15271.stdout         |  2 --
 3 files changed, 27 insertions(+), 6 deletions(-)

diff --git a/rts/StgPrimFloat.c b/rts/StgPrimFloat.c
index 277ae66..f1f6736 100644
--- a/rts/StgPrimFloat.c
+++ b/rts/StgPrimFloat.c
@@ -14,6 +14,7 @@
 
 #include <math.h>
 #include <float.h>
+#include <limits.h>
 
 #define IEEE_FLOATING_POINT 1
 
@@ -47,6 +48,23 @@
 
 #define __abs(a)                (( (a) >= 0 ) ? (a) : (-(a)))
 
+/** Trac #15271: Some large ratios are converted into double incorrectly.
+  * This occurs when StgInt has 64 bits and C int has 32 bits, where wrapping
+  * occurs and an incorrect signed value is passed into ldexp */
+STATIC_INLINE int
+truncExponent(I_ e)
+{
+#if INT_MAX < STG_INT_MAX
+  if (RTS_UNLIKELY(e > INT_MAX))
+    e = INT_MAX;
+#endif
+#if INT_MIN > STG_INT_MIN
+  if (RTS_UNLIKELY(e < INT_MIN))
+    e = INT_MIN;
+#endif
+  return e;
+}
+
 /* Special version for words */
 StgDouble
 __word_encodeDouble (W_ j, I_ e)
@@ -57,7 +75,7 @@ __word_encodeDouble (W_ j, I_ e)
 
   /* Now raise to the exponent */
   if ( r != 0.0 ) /* Lennart suggests this avoids a bug in MIPS's ldexp */
-    r = ldexp(r, e);
+    r = ldexp(r, truncExponent(e));
 
   return r;
 }
@@ -72,7 +90,7 @@ __int_encodeDouble (I_ j, I_ e)
 
   /* Now raise to the exponent */
   if ( r != 0.0 ) /* Lennart suggests this avoids a bug in MIPS's ldexp */
-    r = ldexp(r, e);
+    r = ldexp(r, truncExponent(e));
 
   /* sign is encoded in the size */
   if (j < 0)
@@ -91,7 +109,7 @@ __int_encodeFloat (I_ j, I_ e)
 
   /* Now raise to the exponent */
   if ( r != 0.0 ) /* Lennart suggests this avoids a bug in MIPS's ldexp */
-    r = ldexp(r, e);
+    r = ldexp(r, truncExponent(e));
 
   /* sign is encoded in the size */
   if (j < 0)
@@ -110,7 +128,7 @@ __word_encodeFloat (W_ j, I_ e)
 
   /* Now raise to the exponent */
   if ( r != 0.0 ) /* Lennart suggests this avoids a bug in MIPS's ldexp */
-    r = ldexp(r, e);
+    r = ldexp(r, truncExponent(e));
 
   return r;
 }
diff --git a/testsuite/tests/numeric/should_run/T15271.hs b/testsuite/tests/numeric/should_run/T15271.hs
new file mode 100644
index 0000000..c8f2c95
--- /dev/null
+++ b/testsuite/tests/numeric/should_run/T15271.hs
@@ -0,0 +1,5 @@
+main = do
+  print 1e646457008
+  print 1e646457009 -- T15271: This incorrectly printed 0.0
+  print 1e1555550000 -- This is still infinity
+  print 1e1000000000 -- T15271: This incorrectly printed 0.0
diff --git a/libraries/base/tests/T7034.stdout b/testsuite/tests/numeric/should_run/T15271.stdout
similarity index 66%
copy from libraries/base/tests/T7034.stdout
copy to testsuite/tests/numeric/should_run/T15271.stdout
index 2675153..6e89f24 100644
--- a/libraries/base/tests/T7034.stdout
+++ b/testsuite/tests/numeric/should_run/T15271.stdout
@@ -2,5 +2,3 @@ Infinity
 Infinity
 Infinity
 Infinity
-Infinity
-Infinity