[oct-cvs] Oct commit: oct qd.lisp

[rtoy]

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

Modified Files:
	qd.lisp 
Log Message:
o Replace THREE-SUM with a macro to make sure it's inlined
  everywhere. 
o Update code for new THREE-SUM macro.


--- /project/oct/cvsroot/oct/qd.lisp	2007/09/16 14:23:24	1.49
+++ /project/oct/cvsroot/oct/qd.lisp	2007/09/17 03:08:25	1.50
@@ -53,17 +53,25 @@
     (c::two-sum ,x ,y)))
 
 
-(defun three-sum (a b c)
-  (declare (double-float a b c)
-	   (optimize (speed 3)))
-  (let* ((t1 0d0)
-	 (t2 t1)
-	 (t3 t1))
-    (declare (double-float t1 t2 t3))
-    (two-sum t1 t2 a b)
-    (two-sum a t3 c t1)
-    (two-sum b c t2 t3)
-    (values a b c)))
+(defmacro three-sum (s0 s1 s2 a b c)
+  (let ((aa (gensym))
+	(bb (gensym))
+	(cc (gensym))
+	(t1 (gensym))
+	(t2 (gensym))
+	(t3 (gensym)))
+    `(let* ((,aa ,a)
+	    (,bb ,b)
+	    (,cc ,c)
+	    (,t1 0d0)
+	    (,t2 ,t1)
+	    (,t3 ,t1))
+      (declare (double-float ,t1 ,t2 ,t3))
+      (two-sum ,t1 ,t2 ,aa ,bb)
+      (two-sum ,s0 ,t3 ,cc ,t1)
+      (two-sum ,s1 ,s2 ,t2 ,t3))))
+      
+	  
 
 (defun three-sum2 (a b c)
   (declare (double-float a b c)
@@ -310,17 +318,20 @@
 	 (t0 s0)
 	 (s1 s0)
 	 (t1 s0)
+	 (s2 s0)
 	 (s3 s0))
-    (declare (double-float s0 s1 s3 t0 t1))
+    (declare (double-float s0 s1 s3 t0 t1 s2))
     (two-sum s0 t1 (qd-0 a) (kernel:double-double-hi b))
     (two-sum s1 t1 (qd-1 a) (kernel:double-double-lo b))
     (two-sum s1 t0 s1 t0)
-    (multiple-value-bind (s2 t0 t1)
-	(three-sum (qd-2 a) t0 t1)
-      (two-sum s3 t0 t0 (qd-3 a))
-      (let ((t0 (cl:+ t0 t1)))
-	(multiple-value-call #'%make-qd-d
-	  (renorm-5 s0 s1 s2 s3 t0))))))
+    (three-sum s2 t0 t1 t0 t0 (qd-3 a))
+    (two-sum s3 t0 t0 (qd-3 a))
+    (let ((t0 (cl:+ t0 t1)))
+      (declare (double-float t0))
+      (multiple-value-bind (a0 a1 a2 a3)
+	  (renorm-5 s0 s1 s2 s3 t0)
+	(declare (double-float a0 a1 a2 a3))
+	(%make-qd-d a0 a1 a2 a3)))))
 
 #+cmu
 (defun add-dd-qd (a b)
@@ -400,18 +411,17 @@
 		      (t2 (cl:+ w2 u2))
 		      (t3 (cl:+ w3 u3)))
 		  (two-sum s1 t0 s1 t0)
-		  (multiple-value-bind (s2 t0 t1)
-		      (three-sum s2 t0 t1)
-		    (multiple-value-bind (s3 t0)
-			(three-sum2 s3 t0 t2)
-		      (declare (double-float t0))
-		      (setf t0 (cl:+ t0 t1 t3))
-		      ;; Renormalize
-		      (multiple-value-setq (s0 s1 s2 s3)
-			(renorm-5 s0 s1 s2 s3 t0))
-		      (if (and (zerop a0) (zerop b0))
-			  (%make-qd-d (+ a0 b0) 0d0 0d0 0d0)
-			  (%make-qd-d s0 s1 s2 s3)))))))))))))
+		  (three-sum s2 t0 t1 s2 t0 t1)
+		  (multiple-value-bind (s3 t0)
+		      (three-sum2 s3 t0 t2)
+		    (declare (double-float t0))
+		    (setf t0 (cl:+ t0 t1 t3))
+		    ;; Renormalize
+		    (multiple-value-setq (s0 s1 s2 s3)
+		      (renorm-5 s0 s1 s2 s3 t0))
+		    (if (and (zerop a0) (zerop b0))
+			(%make-qd-d (+ a0 b0) 0d0 0d0 0d0)
+			(%make-qd-d s0 s1 s2 s3))))))))))))
 
 (defun neg-qd (a)
   (declare (type %quad-double a))
@@ -471,17 +481,16 @@
 	       (s2 p0))
 	  (declare (double-float s0 s1 s2 p3))
 	  (two-sum s1 s2 q0 p1)
-	  (multiple-value-bind (s2 q1 p2)
-	      (three-sum s2 q1 p2)
-	    (multiple-value-bind (q1 q2)
-		(three-sum2 q1 q2 p3)
-	      (let ((s3 q1)
-		    (s4 (cl:+ q2 p2)))
-		(multiple-value-bind (s0 s1 s2 s3)
-		    (renorm-5 s0 s1 s2 s3 s4)
-		  (if (zerop s0)
-		      (%make-qd-d (float-sign p0 0d0) 0d0 0d0 0d0)
-		      (%make-qd-d s0 s1 s2 s3)))))))))))
+	  (three-sum s2 q1 p2 s2 q1 p2)
+	  (multiple-value-bind (q1 q2)
+	      (three-sum2 q1 q2 p3)
+	    (let ((s3 q1)
+		  (s4 (cl:+ q2 p2)))
+	      (multiple-value-bind (s0 s1 s2 s3)
+		  (renorm-5 s0 s1 s2 s3 s4)
+		(if (zerop s0)
+		    (%make-qd-d (float-sign p0 0d0) 0d0 0d0 0d0)
+		    (%make-qd-d s0 s1 s2 s3))))))))))
 
 ;; a0 * b0                        0
 ;;      a0 * b1                   1
@@ -597,14 +606,11 @@
 		(multiple-value-bind (p5 q5)
 		    (two-prod a2 b0)
 		  ;; Start accumulation
-		  (multiple-value-setq (p1 p2 q0)
-		    (three-sum p1 p2 q0))
+		  (three-sum p1 p2 q0 p1 p2 q0)
 
 		  ;; six-three-sum of p2, q1, q2, p3, p4, p5
-		  (multiple-value-setq (p2 q1 q2)
-		    (three-sum p2 q1 q2))
-		  (multiple-value-setq (p3 p4 p5)
-		    (three-sum p3 p4 p5))
+		  (three-sum p2 q1 q2 p2 q1 q2)
+		  (three-sum p3 p4 p5 p3 p4 p5)
 		  ;; Compute (s0,s1,s2) = (p2,q1,q2) + (p3,p4,p5)
 		  (let* ((s0 0d0)
 			 (s1 s0)
@@ -779,8 +785,9 @@
 	      (two-sum p3 p4 p3 t0)
 	      (incf p4 (cl:+ q0 t1))
 
-	      (multiple-value-call #'%make-qd-d
-		(renorm-5 p0 p1 p2 p3 p4)))))))))
+	      (multiple-value-bind (a0 a1 a2 a3)
+		  (renorm-5 p0 p1 p2 p3 p4)
+		(%make-qd-d a0 a1 a2 a3)))))))))
 	      
 
 (defun div-qd (a b)