From 26627f70af3ede3610e4092175965e0250456a4d Mon Sep 17 00:00:00 2001
From: Nikodemus Siivola <nikodemus@random-state.net>
Date: Tue, 7 Jun 2011 14:14:21 +0300
Subject: [PATCH] more conservative bounds in FP interval arithmetic

 Make BOUND-BINOP return closed intervals when floating
 point rounding can cause an open bound to close.

 Consider:

  ;; Can return zero.
  (defun fii (x)
    (declare (type (single-float (0.0)) x))
    (/ x 2.0))

  ;; Can return 2.0.
  (defun fii (x y)
    (declare (type (single-float 2.0) x)
             (type (single-float (0.0)) y))
    (+ x y))
---
 NEWS                      |    2 ++
 src/compiler/srctran.lisp |   37 ++++++++++++++++++++++++++++++++-----
 tests/compiler.pure.lisp  |   34 ++++++++++++++++++++++++++++++++++
 3 files changed, 68 insertions(+), 5 deletions(-)

diff --git a/NEWS b/NEWS
index 4cf6c16..98f1151 100644
--- a/NEWS
+++ b/NEWS
@@ -2,6 +2,8 @@
 changes relative to sbcl-1.0.49:
   * enhancement: errors from FD handlers now provide a restart to remove
     the offending handler.
+  * bug fix: bound derivation for floating point operations is now more
+    careful about rounding possibly closing open bounds. (lp#793771)
 
 changes in sbcl-1.0.49 relative to sbcl-1.0.48:
   * minor incompatible change: WITH-LOCKED-HASH-TABLE no longer disables
diff --git a/src/compiler/srctran.lisp b/src/compiler/srctran.lisp
index f8f40b4..aec88cb 100644
--- a/src/compiler/srctran.lisp
+++ b/src/compiler/srctran.lisp
@@ -419,11 +419,38 @@
      (t (,op ,x ,y))))
 
 (defmacro bound-binop (op x y)
-  `(and ,x ,y
-       (with-float-traps-masked (:underflow :overflow :inexact :divide-by-zero)
-         (set-bound (safely-binop ,op (type-bound-number ,x)
-                                  (type-bound-number ,y))
-                    (or (consp ,x) (consp ,y))))))
+  (with-unique-names (xb yb res)
+    `(and ,x ,y
+          (with-float-traps-masked (:underflow :overflow :inexact :divide-by-zero)
+            (let* ((,xb (type-bound-number ,x))
+                   (,yb (type-bound-number ,y))
+                   (,res (safely-binop ,op ,xb ,yb)))
+              (set-bound ,res
+                         (and (or (consp ,x) (consp ,y))
+                              ;; Open bounds can very easily be messed up
+                              ;; by FP rounding, so take care here.
+                              ,(case op
+                                 (*
+                                  ;; Multiplying a greater-than-zero with
+                                  ;; less than one can round to zero.
+                                  `(or (not (fp-zero-p ,res))
+                                       (cond ((and (consp ,x) (fp-zero-p ,xb))
+                                              (>= (abs ,yb) 1))
+                                             ((and (consp ,y) (fp-zero-p ,yb))
+                                              (>= (abs ,xb) 1)))))
+                                 (/
+                                  ;; Dividing a greater-than-zero with
+                                  ;; greater than one can round to zero.
+                                  `(or (not (fp-zero-p ,res))
+                                       (cond ((and (consp ,x) (fp-zero-p ,xb))
+                                              (<= (abs ,yb) 1))
+                                             ((and (consp ,y) (fp-zero-p ,yb))
+                                              (<= (abs ,xb) 1)))))
+                                 ((+ -)
+                                  ;; Adding or subtracting greater-than-zero
+                                  ;; can end up with identity.
+                                  `(and (not (fp-zero-p ,xb))
+                                        (not (fp-zero-p ,yb))))))))))))
 
 (defun coerce-for-bound (val type)
   (if (consp val)
diff --git a/tests/compiler.pure.lisp b/tests/compiler.pure.lisp
index b8643a1..7b2e915 100644
--- a/tests/compiler.pure.lisp
+++ b/tests/compiler.pure.lisp
@@ -3892,3 +3892,37 @@
     (let* ((cell (cons t t)))
       (funcall f cell :ok)
       (assert (equal '(:ok . t) cell)))))
+
+(with-test (:name (:bug-793771 +))
+  (let ((f (compile nil `(lambda (x y)
+                            (declare (type (single-float 2.0) x)
+                                     (type (single-float (0.0)) y))
+                           (+ x y)))))
+    (assert (equal `(function ((single-float 2.0) (single-float (0.0)))
+                              (values (single-float 2.0) &optional))
+                   (sb-kernel:%simple-fun-type f)))))
+
+(with-test (:name (:bug-793771 -))
+  (let ((f (compile nil `(lambda (x y)
+                            (declare (type (single-float * 2.0) x)
+                                     (type (single-float (0.0)) y))
+                           (- x y)))))
+    (assert (equal `(function ((single-float * 2.0) (single-float (0.0)))
+                              (values (single-float * 2.0) &optional))
+                   (sb-kernel:%simple-fun-type f)))))
+
+(with-test (:name (:bug-793771 *))
+  (let ((f (compile nil `(lambda (x)
+                            (declare (type (single-float (0.0)) x))
+                           (* x 0.1)))))
+    (assert (equal `(function ((single-float (0.0)))
+                              (values (or (member 0.0) (single-float (0.0))) &optional))
+                   (sb-kernel:%simple-fun-type f)))))
+
+(with-test (:name (:bug-793771 /))
+  (let ((f (compile nil `(lambda (x)
+                            (declare (type (single-float (0.0)) x))
+                           (/ x 3.0)))))
+    (assert (equal `(function ((single-float (0.0)))
+                              (values (or (member 0.0) (single-float (0.0))) &optional))
+                   (sb-kernel:%simple-fun-type f)))))
-- 
1.7.10.4