Update of /project/oct/cvsroot/oct In directory clnet:/tmp/cvs-serv8007

Modified Files: qd.lisp Log Message: New SCALE-FLOAT-QD implementation that shouldn't suffer from premature overflow/underflow. (Still issue if the exponent is very large or very small, though, but not if the exponent is < 2000 or so.)

--- /project/oct/cvsroot/oct/qd.lisp 2007/09/17 17:15:04 1.52 +++ /project/oct/cvsroot/oct/qd.lisp 2007/09/17 19:04:23 1.53 @@ -1085,6 +1085,7 @@ q0-sign))))))

(declaim (inline scale-float-qd)) +#+(or) (defun scale-float-qd (qd k) (declare (type %quad-double qd) (type fixnum k) @@ -1112,6 +1113,29 @@ (cl:* (qd-2 qd) scale) (cl:* (qd-3 qd) scale))))

+(defun scale-float-qd (qd k) + (declare (type %quad-double qd) + (type (integer -2000 2000) k) + (optimize (speed 3) (space 0))) + ;; Split the exponent in half and do the scaling in two parts. + ;; Requires 2 multiplications, but should not prematurely return 0, + ;; and should be faster than the original version above. + (let* ((k1 (floor k 2)) + (k2 (- k k1)) + (s1 (scale-float 1d0 k1)) + (s2 (scale-float 1d0 k2))) + (multiple-value-bind (q0 q1 q2 q3) + (qd-parts qd) + (%make-qd-d (cl:* (cl:* q0 s1) s2) + (cl:* (cl:* q1 s1) s2) + (cl:* (cl:* q2 s1) s2) + (cl:* (cl:* q3 s1) s2))))) + + + + + + (defun decode-float-qd (q) (declare (type %quad-double q)) (multiple-value-bind (frac exp sign)