Tighter floating-point type constraints in some cases
[sbcl.git] / src / compiler / srctran.lisp
index f8f40b4..cc7cb91 100644 (file)
 (defun bound-func (f x)
   (declare (type function f))
   (and x
-       (with-float-traps-masked (:underflow :overflow :inexact :divide-by-zero)
-         ;; With these traps masked, we might get things like infinity
-         ;; or negative infinity returned. Check for this and return
-         ;; NIL to indicate unbounded.
-         (let ((y (funcall f (type-bound-number x))))
-           (if (and (floatp y)
-                    (float-infinity-p y))
-               nil
-               (set-bound y (consp x)))))))
+       (handler-case
+         (with-float-traps-masked (:underflow :overflow :inexact :divide-by-zero)
+           ;; With these traps masked, we might get things like infinity
+           ;; or negative infinity returned. Check for this and return
+           ;; NIL to indicate unbounded.
+           (let ((y (funcall f (type-bound-number x))))
+             (if (and (floatp y)
+                      (float-infinity-p y))
+                 nil
+                 (set-bound y (consp x)))))
+         ;; Some numerical operations will signal SIMPLE-TYPE-ERROR, e.g.
+         ;; in the course of converting a bignum to a float.  Default to
+         ;; NIL in that case.
+         (simple-type-error ()))))
 
 (defun safe-double-coercion-p (x)
   (or (typep x 'double-float)
      (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)
         `(- (ash x ,len))
         `(ash x ,len))))
 
+;;; These must come before the ones below, so that they are tried
+;;; first. Since %FLOOR and %CEILING are inlined, this allows
+;;; the general case to be handled by TRUNCATE transforms.
+(deftransform floor ((x y))
+  `(%floor x y))
+
+(deftransform ceiling ((x y))
+  `(%ceiling x y))
+
 ;;; If arg is a constant power of two, turn FLOOR into a shift and
 ;;; mask. If CEILING, add in (1- (ABS Y)), do FLOOR and correct a
 ;;; remainder.
       `(if (minusp x)
            (- (logand (- x) ,mask))
            (logand x ,mask)))))
+
+;;; Return an expression to calculate the integer quotient of X and
+;;; constant Y, using multiplication, shift and add/sub instead of
+;;; division. Both arguments must be unsigned, fit in a machine word and
+;;; Y must neither be zero nor a power of two. The quotient is rounded
+;;; towards zero.
+;;; The algorithm is taken from the paper "Division by Invariant
+;;; Integers using Multiplication", 1994 by Torbj\"{o}rn Granlund and
+;;; Peter L. Montgomery, Figures 4.2 and 6.2, modified to exclude the
+;;; case of division by powers of two.
+;;; The algorithm includes an adaptive precision argument.  Use it, since
+;;; we often have sub-word value ranges.  Careful, in this case, we need
+;;; p s.t 2^p > n, not the ceiling of the binary log.
+;;; Also, for some reason, the paper prefers shifting to masking.  Mask
+;;; instead.  Masking is equivalent to shifting right, then left again;
+;;; all the intermediate values are still words, so we just have to shift
+;;; right a bit more to compensate, at the end.
+;;;
+;;; The following two examples show an average case and the worst case
+;;; with respect to the complexity of the generated expression, under
+;;; a word size of 64 bits:
+;;;
+;;; (UNSIGNED-DIV-TRANSFORMER 10 MOST-POSITIVE-WORD) ->
+;;; (ASH (%MULTIPLY (LOGANDC2 X 0) 14757395258967641293) -3)
+;;;
+;;; (UNSIGNED-DIV-TRANSFORMER 7 MOST-POSITIVE-WORD) ->
+;;; (LET* ((NUM X)
+;;;        (T1 (%MULTIPLY NUM 2635249153387078803)))
+;;;   (ASH (LDB (BYTE 64 0)
+;;;             (+ T1 (ASH (LDB (BYTE 64 0)
+;;;                             (- NUM T1))
+;;;                        -1)))
+;;;        -2))
+;;;
+(defun gen-unsigned-div-by-constant-expr (y max-x)
+  (declare (type (integer 3 #.most-positive-word) y)
+           (type word max-x))
+  (aver (not (zerop (logand y (1- y)))))
+  (labels ((ld (x)
+             ;; the floor of the binary logarithm of (positive) X
+             (integer-length (1- x)))
+           (choose-multiplier (y precision)
+             (do* ((l (ld y))
+                   (shift l (1- shift))
+                   (expt-2-n+l (expt 2 (+ sb!vm:n-word-bits l)))
+                   (m-low (truncate expt-2-n+l y) (ash m-low -1))
+                   (m-high (truncate (+ expt-2-n+l
+                                        (ash expt-2-n+l (- precision)))
+                                     y)
+                           (ash m-high -1)))
+                  ((not (and (< (ash m-low -1) (ash m-high -1))
+                             (> shift 0)))
+                   (values m-high shift)))))
+    (let ((n (expt 2 sb!vm:n-word-bits))
+          (precision (integer-length max-x))
+          (shift1 0))
+      (multiple-value-bind (m shift2)
+          (choose-multiplier y precision)
+        (when (and (>= m n) (evenp y))
+          (setq shift1 (ld (logand y (- y))))
+          (multiple-value-setq (m shift2)
+            (choose-multiplier (/ y (ash 1 shift1))
+                               (- precision shift1))))
+        (cond ((>= m n)
+               (flet ((word (x)
+                        `(truly-the word ,x)))
+                 `(let* ((num x)
+                         (t1 (%multiply-high num ,(- m n))))
+                    (ash ,(word `(+ t1 (ash ,(word `(- num t1))
+                                            -1)))
+                         ,(- 1 shift2)))))
+              ((and (zerop shift1) (zerop shift2))
+               (let ((max (truncate max-x y)))
+                 ;; Explicit TRULY-THE needed to get the FIXNUM=>FIXNUM
+                 ;; VOP.
+                 `(truly-the (integer 0 ,max)
+                             (%multiply-high x ,m))))
+              (t
+               `(ash (%multiply-high (logandc2 x ,(1- (ash 1 shift1))) ,m)
+                     ,(- (+ shift1 shift2)))))))))
+
+;;; If the divisor is constant and both args are positive and fit in a
+;;; machine word, replace the division by a multiplication and possibly
+;;; some shifts and an addition. Calculate the remainder by a second
+;;; multiplication and a subtraction. Dead code elimination will
+;;; suppress the latter part if only the quotient is needed. If the type
+;;; of the dividend allows to derive that the quotient will always have
+;;; the same value, emit much simpler code to handle that. (This case
+;;; may be rare but it's easy to detect and the compiler doesn't find
+;;; this optimization on its own.)
+(deftransform truncate ((x y) (word (constant-arg word))
+                        *
+                        :policy (and (> speed compilation-speed)
+                                     (> speed space)))
+  "convert integer division to multiplication"
+  (let* ((y      (lvar-value y))
+         (x-type (lvar-type x))
+         (max-x  (or (and (numeric-type-p x-type)
+                          (numeric-type-high x-type))
+                     most-positive-word)))
+    ;; Division by zero, one or powers of two is handled elsewhere.
+    (when (zerop (logand y (1- y)))
+      (give-up-ir1-transform))
+    `(let* ((quot ,(gen-unsigned-div-by-constant-expr y max-x))
+            (rem (ldb (byte #.sb!vm:n-word-bits 0)
+                      (- x (* quot ,y)))))
+       (values quot rem))))
 \f
 ;;;; arithmetic and logical identity operation elimination
 
                    (consp (arg-info-default info))
                    (not (lambda-var-specvar var))
                    (not (lambda-var-sets var))
-                   (every #'ref-good-for-more-context-p (lambda-var-refs var)))))
+                   (every #'ref-good-for-more-context-p (lambda-var-refs var))
+                   (policy node (= 3 rest-conversion)))))
         (cond (context-ok
                (destructuring-bind (context count &optional used) (arg-info-default info)
                  (declare (ignore used))