Update of /project/oct/cvsroot/oct In directory clnet:/tmp/cvs-serv31151
Modified Files: qd-fun.lisp qd-io.lisp qd-rep.lisp qd.lisp Log Message: o Add default implementation of float-infinity-p, float-nan-p, float-trapping-nan-p. These return NIL by default, unless the Lisp implementation has a suitable version. o Remove CMU conditionalization for float-infinity-p, float-nan-p, float-trapping-nan-p.
--- /project/oct/cvsroot/oct/qd-fun.lisp 2007/10/16 02:39:22 1.86 +++ /project/oct/cvsroot/oct/qd-fun.lisp 2007/10/16 17:05:13 1.87 @@ -62,7 +62,6 @@ ;; 2^k*sqrt(f), and sqrt(f) doesn't have round-off problems. (when (zerop-qd a) (return-from sqrt-qd a)) - #+cmu (when (float-infinity-p (qd-0 a)) (return-from sqrt-qd a))
@@ -215,7 +214,7 @@ (defun expm1-qd (a) "exp(a) - 1, done accurately" (declare (type %quad-double a)) - #+cmu + (when (float-infinity-p (qd-0 a)) (return-from expm1-qd (if (minusp (float-sign (qd-0 a))) @@ -310,7 +309,7 @@ "log1p(x) = log(1+x), done more accurately than just evaluating log(1+x)" (declare (type %quad-double x)) - #+cmu + (when (float-infinity-p (qd-0 x)) x) (log1p-qd/duplication x)) @@ -323,7 +322,6 @@ ((and (zerop-qd a) (plusp (float-sign (qd-0 a)))) (%make-qd-d (/ -1d0 (qd-0 a)) 0d0 0d0 0d0)) - #+cmu ((float-infinity-p (qd-0 a)) a) ((minusp (float-sign (qd-0 a))) @@ -1055,12 +1053,10 @@ ;; where D(x) = exp(x) - 1. This helps for x near 0. (cond ((zerop a) a) - #+cmu ((float-infinity-p (qd-0 a)) a) (t (let ((d (expm1-qd a))) - #+cmu (when (float-infinity-p (qd-0 d)) (return-from sinh-qd d)) (scale-float-qd (add-qd d @@ -1072,7 +1068,6 @@ (declare (type %quad-double a)) ;; cosh(x) = 1/2*(exp(x)+exp(-x)) (let ((e (exp-qd a))) - #+cmu (when (float-infinity-p (qd-0 e)) (return-from cosh-qd e)) (scale-float-qd (add-qd e (div-qd +qd-one+ e)) @@ -1172,7 +1167,6 @@ (add-qd (scale-float-qd (log1p-qd a^2) -1) (log1p-qd (div-qd a (sqrt-qd (add-qd-d a^2 1d0))))))) - #+cmu ((float-infinity-p (qd-0 a)) a) (t @@ -1210,7 +1204,6 @@ (cond ((< (abs (qd-0 a)) (sqrt most-positive-double-float)) (let ((y (sub-qd-d a 1d0))) (log1p-qd (add-qd y (sqrt-qd (mul-qd y (add-qd-d y 2d0))))))) - #+cmu ((float-infinity-p (qd-0 a)) a) (t --- /project/oct/cvsroot/oct/qd-io.lisp 2007/10/15 18:53:44 1.20 +++ /project/oct/cvsroot/oct/qd-io.lisp 2007/10/16 17:05:13 1.21 @@ -298,7 +298,7 @@ (stream stream)) ;; We should do something with colon-p and at-sign-p (declare (ignore colon-p at-sign-p par)) - (cond ((ext:float-infinity-p (qd-0 arg)) + (cond ((float-infinity-p (qd-0 arg)) (qd-output-infinity arg stream)) ((ext:float-nan-p (qd-0 arg)) (qd-output-nan arg stream)) --- /project/oct/cvsroot/oct/qd-rep.lisp 2007/10/15 18:53:44 1.8 +++ /project/oct/cvsroot/oct/qd-rep.lisp 2007/10/16 17:05:13 1.9 @@ -189,3 +189,32 @@ (,a3 (qd-3 ,q))) ,@body)))
+ +;; Some simple support for infinity and NaN. For CMUCL, we can import +;; the desired functions from the EXTENSIONS package. + +;; Implementation for Allegro +#+allegro +(progn +(defmacro float-infinity-p (x) + (= (abs ,x) #.excl::*infinity-double*)) + +(defun float-nan-p (x) + (excl::nan-p x)) + +(defun float-trapping-nan-p (x) + nil) +) ; end progn + + +;; Default implementation. Assume we can't recognize any of these. + +#-(or cmu allegro) +(progn +(defun float-infinity-p (x) + nil) +(defun float-nan-p (x) + nil) +(defun float-trapping-nan-p (x) + nil) +) ; end progn --- /project/oct/cvsroot/oct/qd.lisp 2007/10/16 02:39:22 1.58 +++ /project/oct/cvsroot/oct/qd.lisp 2007/10/16 17:05:13 1.59 @@ -298,8 +298,7 @@ (double-float b) (optimize (speed 3) (space 0)) - #+cmu - (inline ext:float-infinity-p)) + (inline float-infinity-p)) (let* ((c0 0d0) (e c0) (c1 c0) @@ -307,7 +306,7 @@ (c3 c0)) (declare (double-float e c0 c1 c2 c3)) (two-sum c0 e (qd-0 a) b) - #+cmu + (when (float-infinity-p c0) (return-from add-qd-d (%make-qd-d c0 0d0 0d0 0d0))) (two-sum c1 e (qd-1 a) e) @@ -404,9 +403,8 @@ (s2 (cl:+ a2 b2)) (s3 (cl:+ a3 b3))) (declare (double-float s0 s1 s2 s3) - #+cmu - (inline ext:float-infinity-p)) - #+cmu + (inline float-infinity-p)) + (when (float-infinity-p s0) (return-from add-qd (%make-qd-d s0 0d0 0d0 0d0))) (let ((v0 (cl:- s0 a0)) @@ -487,11 +485,10 @@ (double-float b) (optimize (speed 3) (space 0)) - #+cmu - (inline ext:float-infinity-p)) + (inline float-infinity-p)) (multiple-value-bind (p0 q0) (two-prod (qd-0 a) b) - #+cmu + (when (float-infinity-p p0) (return-from mul-qd-d (%make-qd-d p0 0d0 0d0 0d0))) (multiple-value-bind (p1 q1) @@ -608,8 +605,7 @@ (declare (type %quad-double a b) (optimize (speed 3) (space 0)) - #+cmu - (inline ext:float-infinity-p)) + (inline float-infinity-p)) (with-qd-parts (a0 a1 a2 a3) a (declare (double-float a0 a1 a2 a3)) @@ -820,14 +816,13 @@ (declare (type %quad-double a b) (optimize (speed 3) (space 0)) - #+cmu - (inline ext:float-infinity-p)) + (inline float-infinity-p)) (let ((a0 (qd-0 a)) (b0 (qd-0 b))) (let* ((q0 (cl:/ a0 b0)) (r (sub-qd a (mul-qd-d b q0))) (q1 (cl:/ (qd-0 r) b0))) - #+cmu + (when (float-infinity-p q0) (return-from div-qd (%make-qd-d q0 0d0 0d0 0d0))) (setf r (sub-qd r (mul-qd-d b q1))) @@ -863,13 +858,12 @@ (double-float b) (optimize (speed 3) (space 0)) - #+cmu - (inline ext:float-infinity-p)) + (inline float-infinity-p)) ;; Compute approximate quotient using high order doubles, then ;; correct it 3 times using the remainder. Analogous to long ;; division. (let ((q0 (cl:/ (qd-0 a) b))) - #+cmu + (when (float-infinity-p q0) (return-from div-qd-d (%make-qd-d q0 0d0 0d0 0d0)))
@@ -901,8 +895,7 @@ (double-double-float b) (optimize (speed 3) (space 0)) - #+cmu - (inline ext:float-infinity-p)) + (inline float-infinity-p)) (let* ((q0 (cl:/ (qd-0 a) (kernel:double-double-hi b))) (r (sub-qd-dd a (cl:* b q0)))) (when (float-infinity-p q0)