[oct-cvs] Oct commit: oct qd-class.lisp qd-format.lisp qd-fun.lisp qd-io.lisp qd-methods.lisp qd-package.lisp qd-test.lisp rt-tests.lisp tests.lisp timing.lisp

rtoy rtoy at common-lisp.net
Mon Oct 15 18:21:47 UTC 2007


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

Modified Files:
	qd-class.lisp qd-format.lisp qd-fun.lisp qd-io.lisp 
	qd-methods.lisp qd-package.lisp qd-test.lisp rt-tests.lisp 
	tests.lisp timing.lisp 
Log Message:
o Rename QUAD-DOUBLE-INTERNAL package to OCT-INTERNAL, with nickname
  OCTI instead of QDI.
o Rename OCT package to NET.COMMON-LISP.OCT, with a nickname of OCT
o Remove nickname of QD.  (Conflicts with other packages dealing with
  quad-doubles.)
o Update all uses of QDI: to OCTI:

qd-fun.lisp:
o Add REM-PI/2 to do a simpler computation if the arg is small
  enough.  Otherwise, use the accurate but expensive rem operation.
o Renamed ACCURATE-SIN-QD to SIN-QD, etc.
o Update SIN-QD etc to use REM-PI/2.


--- /project/oct/cvsroot/oct/qd-class.lisp	2007/10/13 02:14:43	1.27
+++ /project/oct/cvsroot/oct/qd-class.lisp	2007/10/15 18:21:46	1.28
@@ -47,15 +47,15 @@
 
 #-cmu
 (defmethod print-object ((qd qd-real) stream)
-  (format stream "~/qdi::qd-format/" (qd-value qd)))
+  (format stream "~/octi::qd-format/" (qd-value qd)))
 
 #+cmu
 (defun print-qd (q stream)
   (declare (type %quad-double q))
   (if (or (ext:float-infinity-p (qd-0 q))
 	  (ext:float-nan-p (qd-0 q)))
-      (format stream "~/qdi::qd-format/" q)
-      (format stream "#q~/qdi::qd-format/" q)))
+      (format stream "~/octi::qd-format/" q)
+      (format stream "#q~/octi::qd-format/" q)))
 #+cmu
 (defmethod print-object ((qd qd-real) stream)
   (print-qd (qd-value qd) stream))
@@ -67,7 +67,7 @@
   (make-instance 'qd-real :value (qd-value x)))
 
 (defmethod print-object ((qd qd-complex) stream)
-  (format stream "#q(~<#q~/qdi::qd-format/ #q~/qdi::qd-format/~:@>)"
+  (format stream "#q(~<#q~/octi::qd-format/ #q~/octi::qd-format/~:@>)"
 	  (list (qd-real qd)
 		(qd-imag qd))))
 
--- /project/oct/cvsroot/oct/qd-format.lisp	2007/10/13 02:14:43	1.7
+++ /project/oct/cvsroot/oct/qd-format.lisp	2007/10/15 18:21:46	1.8
@@ -70,7 +70,7 @@
 				       1 0))
 				nil)))
 	    (multiple-value-bind (fstr flen lpoint tpoint)
-		(qdi::qd-to-string (qd-value num) spaceleft fdig k fmin)
+		(octi::qd-to-string (qd-value num) spaceleft fdig k fmin)
 	      (when (and d (zerop d)) (setq tpoint nil))
 	      (when w 
 		(cl:decf spaceleft flen)
--- /project/oct/cvsroot/oct/qd-fun.lisp	2007/10/15 03:26:38	1.83
+++ /project/oct/cvsroot/oct/qd-fun.lisp	2007/10/15 18:21:46	1.84
@@ -646,7 +646,52 @@
 	  (values mod f)
 	  (values (mod (1+ mod) 4) (sub-qd f +qd-pi/2+))))))
 
-(defun accurate-sin-qd (a)
+(defun rem-pi/2-int (qd)
+  ;; Compute qd rem pi/2 = k*pi/2+y.  So we compute k + y*2/pi =
+  ;; qd*2/pi.
+  ;;
+  ;; First convert qd to 2^e*I.  We already have 2/pi in the form
+  ;; 2^-1584*J.  Then qd*2/pi = 2^(e-1584)*I*J.  Extract out the
+  ;; integer and fractional parts of this.  For the integer part we
+  ;; only need it modulo 4, because of the periodicity.  For the
+  ;; fractional part, we only need 212 (or so bits of fraction).
+  ;;
+  ;; FIXME: But we don't really need to compute all the bits of I*J.
+  ;; In the product, we really only need the 2 bits to the left of the
+  ;; binary point, and then 212 bits to the right.  This doesn't
+  ;; require doing the full multiplication.
+  (multiple-value-bind (i e s)
+      (integer-decode-qd qd)
+    (let* ((exp (- e (integer-length +2/pi-bits+)))
+	   (prod (* (* s i) +2/pi-bits+))
+	   (mod (ldb (byte 2 (- exp)) prod))
+	   ;; A quad-double has about 212 bits, but we add another 53
+	   ;; (5 doubles) for some extra accuracty.
+	   (qd-bits 265)
+	   (frac (ldb (byte qd-bits (- (- exp) qd-bits)) prod))
+	   (f (mul-qd (scale-float-qd (rational-to-qd frac) (- qd-bits))
+		      +qd-pi/2+)))
+      ;; We want the remainder part to be <= pi/4 because the trig
+      ;; functions want that.  So if the fraction is too big, adjust
+      ;; it, and mod value.
+      (if (<= (abs (qd-0 f)) (/ pi 4))
+	  (values mod f)
+	  (values (mod (1+ mod) 4) (sub-qd f +qd-pi/2+))))))
+
+(defun rem-pi/2 (a)
+  ;; If the number is small enough, we don't need to use the full
+  ;; precision algorithm to compute the remainder.  The value of 1024
+  ;; here is rather arbitrary.  We should do an analysis to figure
+  ;; where the breakpoint should be.
+  (cond ((<= (abs (qd-0 a)) 256)
+	 (let ((quot (truncate (qd-0 (nint-qd (div-qd a +qd-pi/2+))))))
+	   (values (mod quot 4)
+		   (sub-qd a (mul-qd-d +qd-pi/2+ (float quot 1d0))))))
+	(t
+	 (rem-pi/2-int a))))
+      
+
+(defun sin-qd (a)
   "Sin(a)"
   (declare (type %quad-double a))
   ;; To compute sin(x), choose integers a, b so that
@@ -666,10 +711,10 @@
   ;; cos(s+k*pi/1024) = cos(s)*cos(k*pi/1024)
   ;;                     - sin(s)*sin(k*pi/1024)
   (when (zerop-qd a)
-    (return-from accurate-sin-qd a))
+    (return-from sin-qd a))
 
   (multiple-value-bind (j r)
-      (rem-pi/2-int a)
+      (rem-pi/2 a)
     (multiple-value-bind (k tmp)
 	(divrem-qd r +qd-pi/1024+)
       (let* ((k (truncate (qd-0 k)))
@@ -723,7 +768,7 @@
 		   ;; cos(j*pi/2) = 0, sin(j*pi/2) = -1
 		   (neg-qd c)))))))))
 
-(defun accurate-cos-qd (a)
+(defun cos-qd (a)
   "Cos(a)"
   ;; Just like sin-qd, but for cos.
   (declare (type %quad-double a))
@@ -744,10 +789,10 @@
   ;; cos(s+k*pi/1024) = cos(s)*cos(k*pi/1024)
   ;;                     - sin(s)*sin(k*pi/1024)
   (when (zerop-qd a)
-    (return-from accurate-cos-qd +qd-one+))
+    (return-from cos-qd +qd-one+))
 
   (multiple-value-bind (j r)
-      (rem-pi/2-int a)
+      (rem-pi/2 a)
     (multiple-value-bind (k tmp)
 	(divrem-qd r +qd-pi/1024+)
       (let* ((k (truncate (qd-0 k)))
@@ -811,7 +856,7 @@
 	      +qd-one+)))
 
   (multiple-value-bind (j r)
-      (rem-pi/2-int a)
+      (rem-pi/2 a)
     (multiple-value-bind (k tmp)
 	(divrem-qd tmp +qd-pi/1024+)
       (let* ((k (truncate (qd-0 k)))
@@ -930,9 +975,9 @@
 	 (yy (div-qd y r)))
     #+nil
     (progn
-      (format t "r  = ~/qdi::qd-format/~%" r)
-      (format t "xx = ~/qdi::qd-format/~%" xx)
-      (format t "yy = ~/qdi::qd-format/~%" yy))
+      (format t "r  = ~/octi::qd-format/~%" r)
+      (format t "xx = ~/octi::qd-format/~%" xx)
+      (format t "yy = ~/octi::qd-format/~%" yy))
     
     ;; Compute double-precision approximation to atan
     (let ((z (make-qd-d (atan (qd-0 y) (qd-0 x))))
--- /project/oct/cvsroot/oct/qd-io.lisp	2007/10/10 15:24:07	1.18
+++ /project/oct/cvsroot/oct/qd-io.lisp	2007/10/15 18:21:46	1.19
@@ -376,7 +376,7 @@
 			 (scale-float (float q1 1d0)
 				      (cl:- len (* 2 53))))
 		 #+(or)
-		 (format t "~/qdi::qd-format/~%" (mul-qd xx yy)))
+		 (format t "~/octi::qd-format/~%" (mul-qd xx yy)))
 	       (if (minusp sign)
 		   (neg-qd (mul-qd xx yy))
 		   (mul-qd xx yy))))))))
--- /project/oct/cvsroot/oct/qd-methods.lisp	2007/10/13 02:14:43	1.63
+++ /project/oct/cvsroot/oct/qd-methods.lisp	2007/10/15 18:21:47	1.64
@@ -26,23 +26,23 @@
 (in-package #:oct)
 
 (defconstant +pi+
-  (make-instance 'qd-real :value qdi:+qd-pi+)
+  (make-instance 'qd-real :value octi:+qd-pi+)
   "Quad-double value of pi")
 
 (defconstant +pi/2+
-  (make-instance 'qd-real :value qdi:+qd-pi/2+)
+  (make-instance 'qd-real :value octi:+qd-pi/2+)
   "Quad-double value of pi/2")
 
 (defconstant +pi/4+
-  (make-instance 'qd-real :value qdi:+qd-pi/4+)
+  (make-instance 'qd-real :value octi:+qd-pi/4+)
   "Quad-double value of pi/4")
 
 (defconstant +2pi+
-  (make-instance 'qd-real :value qdi:+qd-2pi+)
+  (make-instance 'qd-real :value octi:+qd-2pi+)
   "Quad-double value of 2*pi")
 
 (defconstant +log2+
-  (make-instance 'qd-real :value qdi:+qd-log2+)
+  (make-instance 'qd-real :value octi:+qd-log2+)
   "Quad-double value of log(2), natural log of 2")
 
 #+cmu
@@ -57,7 +57,7 @@
 
 (defconstant +most-positive-quad-double-float+
   (make-instance 'qd-real
-		 :value (qdi::%make-qd-d most-positive-double-float
+		 :value (octi::%make-qd-d most-positive-double-float
 					 (cl:scale-float most-positive-double-float (cl:* 1 -53))
 					 (cl:scale-float most-positive-double-float (cl:* 2 -53))
 					 (cl:scale-float most-positive-double-float (cl:* 3 -53)))))
@@ -79,13 +79,7 @@
 
 (defmethod make-qd ((x cl:rational))
   ;; We should do something better than this.
-  (let ((top (numerator x))
-	(bot (denominator x)))
-    (if (= bot 1)
-	(make-instance 'qd-real :value (qdi::make-float (signum top) (abs top) 0 0 0))
-	(make-instance 'qd-real
-		       :value (div-qd (qdi::make-float (signum top) (abs top) 0 0 0)
-				      (qdi::make-float (signum bot) (abs bot) 0 0 0))))))
+  (make-instance 'qd-real :value (rational-to-qd x)))
 
 
 (defmethod add1 ((a number))
@@ -424,7 +418,7 @@
 (declaim (inline qd-cssqs))
 (defun qd-cssqs (z)
   (multiple-value-bind (rho k)
-      (qdi::hypot-aux-qd (qd-value (realpart z))
+      (octi::hypot-aux-qd (qd-value (realpart z))
 			 (qd-value (imagpart z)))
     (values (make-instance 'qd-real :value rho)
 	    k)))
@@ -811,7 +805,7 @@
 
 (defmethod random ((x qd-real) &optional (state *random-state*))
   (* x (make-instance 'qd-real
-		      :value (qdi:random-qd state))))
+		      :value (octi:random-qd state))))
 
 (defmethod float-digits ((x cl:real))
   (cl:float-digits x))
--- /project/oct/cvsroot/oct/qd-package.lisp	2007/10/10 15:24:07	1.42
+++ /project/oct/cvsroot/oct/qd-package.lisp	2007/10/15 18:21:47	1.43
@@ -23,9 +23,9 @@
 ;;;; FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
 ;;;; OTHER DEALINGS IN THE SOFTWARE.
 
-(defpackage #:quad-double-internal
+(defpackage #:oct-internal
   (:use #:cl #+cmu #:extensions)
-  (:nicknames #:qdi)
+  (:nicknames #:octi)
   (:export #:%quad-double
 	   #:read-qd
 	   #:add-qd
@@ -102,9 +102,9 @@
 		#:two-prod
 		#:two-sqr))
 
-(defpackage #:quad-double
-  (:use #:cl #:quad-double-internal)
-  (:nicknames #:oct #:qd)
+(defpackage #:net.common-lisp.oct
+  (:use #:cl #:oct-internal)
+  (:nicknames #:oct)
   (:shadow #:+
 	   #:-
 	   #:*
@@ -240,6 +240,7 @@
 	   #:float-digits
 	   #:rational
 	   #:rationalize
+	   #:make-qd
 	   )
   ;; Constants
   (:export #:+pi+
--- /project/oct/cvsroot/oct/qd-test.lisp	2007/09/18 03:05:56	1.19
+++ /project/oct/cvsroot/oct/qd-test.lisp	2007/10/15 18:21:47	1.20
@@ -38,8 +38,8 @@
 	(cl:- (log err 2d0)))))
 
 (defun print-result (est true)
-  (format t "est: ~/qdi::qd-format/~%" est)
-  (format t "tru: ~/qdi::qd-format/~%" true)
+  (format t "est: ~/octi::qd-format/~%" est)
+  (format t "tru: ~/octi::qd-format/~%" true)
   (format t "err: ~A~%" (qd-0 (sub-qd est true)))
   (format t "bits: ~,1f~%" (bit-accuracy est true)))
   
@@ -218,24 +218,24 @@
   (let* ((arg (div-qd +qd-one+ (sqrt-qd (make-qd-d 3d0))))
 	 (y (div-qd (funcall fun arg) +qd-pi+))
 	 (true (div-qd +qd-one+ (make-qd-d 6d0))))
-    (format t "atan(1/sqrt(3))/pi = ~/qdi::qd-format/~%" y)
-    (format t "1/6                = ~/qdi::qd-format/~%" true)
+    (format t "atan(1/sqrt(3))/pi = ~/octi::qd-format/~%" y)
+    (format t "1/6                = ~/octi::qd-format/~%" true)
     (format t "bits               = ~,1f~%"
 	    (bit-accuracy y true)))
   ;; atan(sqrt(3)) = %pi/3
   (let* ((arg (sqrt-qd (make-qd-d 3d0)))
 	 (y (div-qd (funcall fun arg) +qd-pi+))
 	 (true (div-qd +qd-one+ (make-qd-d 3d0))))
-    (format t "atan(sqrt(3))/pi   = ~/qdi::qd-format/~%" y)
-    (format t "1/3                = ~/qdi::qd-format/~%" true)
+    (format t "atan(sqrt(3))/pi   = ~/octi::qd-format/~%" y)
+    (format t "1/3                = ~/octi::qd-format/~%" true)
     (format t "bits               = ~,1f~%"
 	    (bit-accuracy y true)))
   ;; atan(1) = %pi/4
   (let* ((arg +qd-one+)
 	 (y (div-qd (funcall fun arg) +qd-pi+))
 	 (true (div-qd +qd-one+ (make-qd-d 4d0))))
-    (format t "atan(1)/pi         = ~/qdi::qd-format/~%" y)
-    (format t "1/4                = ~/qdi::qd-format/~%" true)
+    (format t "atan(1)/pi         = ~/octi::qd-format/~%" y)
+    (format t "1/4                = ~/octi::qd-format/~%" true)
     (format t "bits               = ~,1f~%"
 	    (bit-accuracy y true))))
 
@@ -244,22 +244,22 @@
   (let* ((arg (div-qd +qd-pi+ (make-qd-d 6d0)))
 	 (y (sin-qd arg))
 	 (true (make-qd-d 0.5d0)))
-    (format t "sin(pi/6)      = ~/qdi::qd-format/~%" y)
-    (format t "1/2            = ~/qdi::qd-format/~%" true)
+    (format t "sin(pi/6)      = ~/octi::qd-format/~%" y)
+    (format t "1/2            = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (div-qd +qd-pi+ (make-qd-d 4d0)))
 	 (y (sin-qd arg))
 	 (true (sqrt-qd (make-qd-d 0.5d0))))
-    (format t "sin(pi/4)      = ~/qdi::qd-format/~%" y)
-    (format t "1/sqrt(2)      = ~/qdi::qd-format/~%" true)
+    (format t "sin(pi/4)      = ~/octi::qd-format/~%" y)
+    (format t "1/sqrt(2)      = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (div-qd +qd-pi+ (make-qd-d 3d0)))
 	 (y (sin-qd arg))
 	 (true (div-qd (sqrt-qd (make-qd-d 3d0)) (make-qd-d 2d0))))
-    (format t "sin(pi/3)      = ~/qdi::qd-format/~%" y)
-    (format t "sqrt(3)/2      = ~/qdi::qd-format/~%" true)
+    (format t "sin(pi/3)      = ~/octi::qd-format/~%" y)
+    (format t "sqrt(3)/2      = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   )
@@ -269,22 +269,22 @@
   (let* ((arg (div-qd +qd-pi+ (make-qd-d 6d0)))
 	 (y (funcall f arg))
 	 (true (div-qd +qd-one+ (sqrt-qd (make-qd-d 3d0)))))
-    (format t "tan(pi/6)      = ~/qdi::qd-format/~%" y)
-    (format t "1/sqrt(3)      = ~/qdi::qd-format/~%" true)
+    (format t "tan(pi/6)      = ~/octi::qd-format/~%" y)
+    (format t "1/sqrt(3)      = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (div-qd +qd-pi+ (make-qd-d 4d0)))
 	 (y (funcall f arg))
 	 (true +qd-one+))
-    (format t "tan(pi/4)      = ~/qdi::qd-format/~%" y)
-    (format t "1              = ~/qdi::qd-format/~%" true)
+    (format t "tan(pi/4)      = ~/octi::qd-format/~%" y)
+    (format t "1              = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (div-qd +qd-pi+ (make-qd-d 3d0)))
 	 (y (funcall f arg))
 	 (true (sqrt-qd (make-qd-d 3d0))))
-    (format t "tan(pi/3)      = ~/qdi::qd-format/~%" y)
-    (format t "sqrt(3)        = ~/qdi::qd-format/~%" true)
+    (format t "tan(pi/3)      = ~/octi::qd-format/~%" y)
+    (format t "sqrt(3)        = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   )
@@ -294,22 +294,22 @@
   (let* ((arg (make-qd-d 0.5d0))
 	 (y (asin-qd arg))
 	 (true (div-qd +qd-pi+ (make-qd-d 6d0))))
-    (format t "asin(1/2)      = ~/qdi::qd-format/~%" y)
-    (format t "pi/6           = ~/qdi::qd-format/~%" true)
+    (format t "asin(1/2)      = ~/octi::qd-format/~%" y)
+    (format t "pi/6           = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (sqrt-qd (make-qd-d 0.5d0)))
 	 (y (asin-qd arg))
 	 (true (div-qd +qd-pi+ (make-qd-d 4d0))))
-    (format t "asin(1/sqrt(2))= ~/qdi::qd-format/~%" y)
-    (format t "pi/4           = ~/qdi::qd-format/~%" true)
+    (format t "asin(1/sqrt(2))= ~/octi::qd-format/~%" y)
+    (format t "pi/4           = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (div-qd (sqrt-qd (make-qd-d 3d0)) (make-qd-d 2d0)))
 	 (y (asin-qd arg))
 	 (true (div-qd +qd-pi+ (make-qd-d 3d0))))
-    (format t "asin(sqrt(3)/2)= ~/qdi::qd-format/~%" y)
-    (format t "pi/3           = ~/qdi::qd-format/~%" true)
+    (format t "asin(sqrt(3)/2)= ~/octi::qd-format/~%" y)
+    (format t "pi/3           = ~/octi::qd-format/~%" true)
     (format t "bits           = ~,1f~%"
 	    (bit-accuracy y true)))
   )
@@ -319,22 +319,22 @@
   (let* ((arg (make-qd-d 2d0))
 	 (y (funcall f arg))
 	 (true +qd-log2+))
-    (format t "log(2)      = ~/qdi::qd-format/~%" y)
-    (format t "log(2)      = ~/qdi::qd-format/~%" true)
+    (format t "log(2)      = ~/octi::qd-format/~%" y)
+    (format t "log(2)      = ~/octi::qd-format/~%" true)
     (format t "bits        = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (make-qd-d 10d0))
 	 (y (funcall f arg))
 	 (true +qd-log10+))
-    (format t "log(10)     = ~/qdi::qd-format/~%" y)
-    (format t "log(10)     = ~/qdi::qd-format/~%" true)
+    (format t "log(10)     = ~/octi::qd-format/~%" y)
+    (format t "log(10)     = ~/octi::qd-format/~%" true)
     (format t "bits        = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (add-d-qd 1d0 (scale-float-qd (make-qd-d 1d0) -80)))
 	 (y (funcall f arg))
 	 (true (qd-from-string "8.2718061255302767487140834995607996176476940491239977084112840149578911975528492q-25")))
-    (format t "log(1+2^-80) = ~/qdi::qd-format/~%" y)
-    (format t "log(1+2^-80) = ~/qdi::qd-format/~%" true)
+    (format t "log(1+2^-80) = ~/octi::qd-format/~%" y)
+    (format t "log(1+2^-80) = ~/octi::qd-format/~%" true)
     (format t "bits         = ~,1f~%"
 	    (bit-accuracy y true)))
   )
@@ -344,22 +344,22 @@
   (let* ((arg (make-qd-d 1d0))
 	 (y (funcall f arg))
 	 (true +qd-log2+))
-    (format t "log1p(1)    = ~/qdi::qd-format/~%" y)
-    (format t "log(2)      = ~/qdi::qd-format/~%" true)
+    (format t "log1p(1)    = ~/octi::qd-format/~%" y)
+    (format t "log(2)      = ~/octi::qd-format/~%" true)
     (format t "bits        = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (make-qd-d 9d0))
 	 (y (funcall f arg))
 	 (true (qd-from-string "2.3025850929940456840179914546843642076011014886287729760333279009675726096773525q0")))
-    (format t "log1p(9)    = ~/qdi::qd-format/~%" y)
-    (format t "log(10)     = ~/qdi::qd-format/~%" true)
+    (format t "log1p(9)    = ~/octi::qd-format/~%" y)
+    (format t "log(10)     = ~/octi::qd-format/~%" true)
     (format t "bits        = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg (scale-float-qd (make-qd-d 1d0) -80))
 	 (y (funcall f arg))
 	 (true (qd-from-string "8.2718061255302767487140834995607996176476940491239977084112840149578911975528492q-25")))
-    (format t "log1p(2^-80) = ~/qdi::qd-format/~%" y)
-    (format t "log(1+2^-80) = ~/qdi::qd-format/~%" true)
+    (format t "log1p(2^-80) = ~/octi::qd-format/~%" y)
+    (format t "log(1+2^-80) = ~/octi::qd-format/~%" true)
     (format t "bits         = ~,1f~%"
 	    (bit-accuracy y true)))
   )
@@ -369,23 +369,23 @@
   (let* ((arg +qd-zero+)
 	 (y (funcall f arg))
 	 (true +qd-zero+))
-    (format t "exp(0)-1    = ~/qdi::qd-format/~%" y)
-    (format t "0           = ~/qdi::qd-format/~%" true)
+    (format t "exp(0)-1    = ~/octi::qd-format/~%" y)
+    (format t "0           = ~/octi::qd-format/~%" true)
     (format t "bits        = ~,1f~%"
 	    (bit-accuracy y true)))
   (let* ((arg +qd-one+)
 	 (y (funcall f arg))
 	 (true (qd-from-string "1.7182818284590452353602874713526624977572470936999595749669676277240766303535475945713821785251664274274663919320030599218174135966290435729003342952q0")))
-    (format t "exp(1)-1    = ~/qdi::qd-format/~%" y)
-    (format t "e-1         = ~/qdi::qd-format/~%" true)
+    (format t "exp(1)-1    = ~/octi::qd-format/~%" y)
+    (format t "e-1         = ~/octi::qd-format/~%" true)
     (format t "bits        = ~,1f~%"
 	    (bit-accuracy y true)))
 
   (let* ((arg (scale-float-qd +qd-one+ -100))
 	 (y (funcall f arg))
 	 (true (qd-from-string "7.888609052210118054117285652830973804370994921943802079729680186943164342372119432861876389514693341738324702996270767390039172777809233288470357147q-31")))
-    (format t "exp(2^-80)-1 = ~/qdi::qd-format/~%" y)
-    (format t "exp(2^-80)-1 = ~/qdi::qd-format/~%" true)
+    (format t "exp(2^-80)-1 = ~/octi::qd-format/~%" y)
+    (format t "exp(2^-80)-1 = ~/octi::qd-format/~%" true)
     (format t "bits         = ~,1f~%"
 	    (bit-accuracy y true)))
 
--- /project/oct/cvsroot/oct/rt-tests.lisp	2007/10/15 15:45:33	1.4
+++ /project/oct/cvsroot/oct/rt-tests.lisp	2007/10/15 18:21:47	1.5
@@ -50,7 +50,7 @@
 ;; Pi via Machin's formula
 (rt:deftest oct.pi.machin
     (let* ((*standard-output* *null*)
-	   (val (make-instance 'qd-real :value (qdi::test2 nil)))
+	   (val (make-instance 'qd-real :value (octi::test2 nil)))
 	   (true oct:+pi+))
       (check-accuracy 213 val true))
   nil)
@@ -58,7 +58,7 @@
 ;; Pi via Salamin-Brent algorithm
 (rt:deftest oct.pi.salamin-brent
     (let* ((*standard-output* *null*)
-	   (val (make-instance 'qd-real :value (qdi::test3 nil)))
+	   (val (make-instance 'qd-real :value (octi::test3 nil)))
 	   (true oct:+pi+))
       (check-accuracy 202 val true))
   nil)
@@ -66,7 +66,7 @@
 ;; Pi via Borweign's Quartic formula
 (rt:deftest oct.pi.borweign
     (let* ((*standard-output* *null*)
-	   (val (make-instance 'qd-real :value (qdi::test4 nil)))
+	   (val (make-instance 'qd-real :value (octi::test4 nil)))
 	   (true oct:+pi+))
       (check-accuracy 211 val true))
   nil)
@@ -74,16 +74,16 @@
 ;; e via Taylor series
 (rt:deftest oct.e.taylor
     (let* ((*standard-output* *null*)
-	   (val (make-instance 'qd-real :value (qdi::test5 nil)))
-	   (true (make-instance 'qd-real :value qdi::+qd-e+)))
+	   (val (make-instance 'qd-real :value (octi::test5 nil)))
+	   (true (make-instance 'qd-real :value octi::+qd-e+)))
       (check-accuracy 212 val true))
   nil)
 
 ;; log(2) via Taylor series
 (rt:deftest oct.log2.taylor
     (let* ((*standard-output* *null*)
-	   (val (make-instance 'qd-real :value (qdi::test6 nil)))
-	   (true (make-instance 'qd-real :value qdi::+qd-log2+)))
+	   (val (make-instance 'qd-real :value (octi::test6 nil)))
+	   (true (make-instance 'qd-real :value octi::+qd-log2+)))
       (check-accuracy 212 val true))
   nil)
 
@@ -126,7 +126,7 @@
 
 (defun atan-qd/duplication (arg)
   (make-instance 'qd-real
-		 :value (qdi::atan-qd/duplication (qd-value arg))))
+		 :value (octi::atan-qd/duplication (qd-value arg))))
 
 ;;; Tests of atan where we know the analytical result.  Same tests,
 ;;; but using the atan duplication formula.
@@ -139,15 +139,15 @@
 
 (rt:deftest oct.atan/dup.2
     (let* ((arg (sqrt #q3))
-	 (y (/ (atan-qd/duplication arg) +pi+))
-	 (true (/ #q3)))
+	   (y (/ (atan-qd/duplication arg) +pi+))
+	   (true (/ #q3)))
       (check-accuracy 212 y true))
   nil)
 
 (rt:deftest oct.atan/dup.3
     (let* ((arg #q1)
-	 (y (/ (atan-qd/duplication arg) +pi+))
-	 (true (/ #q4)))
+	   (y (/ (atan-qd/duplication arg) +pi+))
+	   (true (/ #q4)))
     (check-accuracy 212 y true))
   nil)
 
@@ -160,8 +160,8 @@
 
 (rt:deftest oct.atan/dup.5
     (let* ((arg #q-1q100)
-	 (y (/ (atan-qd/duplication arg) +pi+))
-	 (true #q-.5))
+	   (y (/ (atan-qd/duplication arg) +pi+))
+	   (true #q-.5))
     (check-accuracy 212 y true))
   nil)
 
@@ -169,7 +169,7 @@
 ;;; but using a CORDIC implementation.
 (defun atan-qd/cordic (arg)
   (make-instance 'qd-real
-		 :value (qdi::atan-qd/cordic (qd-value arg))))
+		 :value (octi::atan-qd/cordic (qd-value arg))))
 
 (rt:deftest oct.atan/cordic.1
     (let* ((arg (/ (sqrt #q3)))
@@ -180,15 +180,15 @@
 
 (rt:deftest oct.atan/cordic.2
     (let* ((arg (sqrt #q3))
-	 (y (/ (atan-qd/cordic arg) +pi+))
-	 (true (/ #q3)))
+	   (y (/ (atan-qd/cordic arg) +pi+))
+	   (true (/ #q3)))
       (check-accuracy 212 y true))
   nil)
 
 (rt:deftest oct.atan/cordic.3
     (let* ((arg #q1)
-	 (y (/ (atan-qd/cordic arg) +pi+))
-	 (true (/ #q4)))
+	   (y (/ (atan-qd/cordic arg) +pi+))
+	   (true (/ #q4)))
     (check-accuracy 212 y true))
   nil)
 
@@ -201,8 +201,8 @@
 
 (rt:deftest oct.atan/cordic.5
     (let* ((arg #q-1q100)
-	 (y (/ (atan-qd/cordic arg) +pi+))
-	 (true #q-.5))
+	   (y (/ (atan-qd/cordic arg) +pi+))
+	   (true #q-.5))
     (check-accuracy 212 y true))
   nil)
 
@@ -210,25 +210,39 @@
 ;;; Tests of sin where we know the analytical result.
 (rt:deftest oct.sin.1
     (let* ((arg (/ +pi+ 6))
-	 (y (sin arg))
-	 (true #q.5))
+	   (y (sin arg))
+	   (true #q.5))
     (check-accuracy 212 y true))
   nil)
 
 (rt:deftest oct.sin.2
     (let* ((arg (/ +pi+ 4))
-	 (y (sin arg))
-	 (true (sqrt #q.5)))
+	   (y (sin arg))
+	   (true (sqrt #q.5)))
     (check-accuracy 212 y true))
   nil)
 
 (rt:deftest oct.sin.3
     (let* ((arg (/ +pi+ 3))
-	 (y (sin arg))
-	 (true (/ (sqrt #q3) 2)))
+	   (y (sin arg))
+	   (true (/ (sqrt #q3) 2)))
     (check-accuracy 212 y true))
   nil)
 
+(rt:deftest oct.big-sin.1
+    (let* ((arg (oct:make-qd (ash 1 120)))
+	   (y (sin arg))
+	   (true #q3.778201093607520226555484700569229919605866976512306642257987199414885q-1))
+      (check-accuracy 205 y true))
+  nil)
+
+(rt:deftest oct.big-sin.2
+    (let* ((arg (oct:make-qd (ash 1 1023)))
+	   (y (sin arg))
+	   (true #q5.631277798508840134529434079444683477103854907361251399182750155357133q-1))
+      (check-accuracy 205 y true))
+  nil)
+
 ;;; Tests of tan where we know the analytical result.
 (rt:deftest oct.tan.1
     (let* ((arg (/ +pi+ 6))
@@ -256,7 +270,7 @@
 
 (defun tan/cordic (arg)
   (make-instance 'qd-real
-		 :value (qdi::tan-qd/cordic (qd-value arg))))
+		 :value (octi::tan-qd/cordic (qd-value arg))))
 
 (rt:deftest oct.tan/cordic.1
     (let* ((arg (/ +pi+ 6))
@@ -307,14 +321,14 @@
 (rt:deftest oct.log.1
     (let* ((arg #q2)
 	   (y (log arg))
-	   (true (make-instance 'qd-real :value qdi::+qd-log2+)))
+	   (true (make-instance 'qd-real :value octi::+qd-log2+)))
       (check-accuracy 212 y true))
   nil)
 
 (rt:deftest oct.log.2
     (let* ((arg #q10)
 	 (y (log arg))
-	 (true (make-instance 'qd-real :value qdi::+qd-log10+)))
+	 (true (make-instance 'qd-real :value octi::+qd-log10+)))
     (check-accuracy 207 y true))
   nil)
 
@@ -329,19 +343,19 @@
 
 (defun log/newton (arg)
   (make-instance 'qd-real
-		 :value (qdi::log-qd/newton (qd-value arg))))
+		 :value (octi::log-qd/newton (qd-value arg))))
 
 (rt:deftest oct.log/newton.1
     (let* ((arg #q2)
 	   (y (log/newton arg))
-	   (true (make-instance 'qd-real :value qdi::+qd-log2+)))
+	   (true (make-instance 'qd-real :value octi::+qd-log2+)))
       (check-accuracy 212 y true))
   nil)
 
 (rt:deftest oct.log/newton.2
     (let* ((arg #q10)
 	 (y (log/newton arg))
-	 (true (make-instance 'qd-real :value qdi::+qd-log10+)))
+	 (true (make-instance 'qd-real :value octi::+qd-log10+)))
     (check-accuracy 207 y true))
   nil)
 
@@ -356,19 +370,19 @@
 
 (defun log/agm (arg)
   (make-instance 'qd-real
-		 :value (qdi::log-qd/agm (qd-value arg))))
+		 :value (octi::log-qd/agm (qd-value arg))))
 
 (rt:deftest oct.log/agm.1
     (let* ((arg #q2)
 	   (y (log/agm arg))
-	   (true (make-instance 'qd-real :value qdi::+qd-log2+)))
+	   (true (make-instance 'qd-real :value octi::+qd-log2+)))
       (check-accuracy 203 y true))
   nil)
 
 (rt:deftest oct.log/agm.2
     (let* ((arg #q10)
 	 (y (log/agm arg))
-	 (true (make-instance 'qd-real :value qdi::+qd-log10+)))
+	 (true (make-instance 'qd-real :value octi::+qd-log10+)))
     (check-accuracy 205 y true))
   nil)
 
@@ -383,19 +397,19 @@
 
 (defun log/agm2 (arg)
   (make-instance 'qd-real
-		 :value (qdi::log-qd/agm2 (qd-value arg))))
+		 :value (octi::log-qd/agm2 (qd-value arg))))
 
 (rt:deftest oct.log/agm2.1
     (let* ((arg #q2)
 	   (y (log/agm2 arg))
-	   (true (make-instance 'qd-real :value qdi::+qd-log2+)))
+	   (true (make-instance 'qd-real :value octi::+qd-log2+)))
       (check-accuracy 203 y true))
   nil)
 
 (rt:deftest oct.log/agm2.2
     (let* ((arg #q10)
 	 (y (log/agm2 arg))
-	 (true (make-instance 'qd-real :value qdi::+qd-log10+)))
+	 (true (make-instance 'qd-real :value octi::+qd-log10+)))
     (check-accuracy 205 y true))
   nil)
 
@@ -409,19 +423,19 @@
 ;;; Tests of log using AGM3, a faster variation of AGM2.
 (defun log/agm3 (arg)
   (make-instance 'qd-real
-		 :value (qdi::log-qd/agm3 (qd-value arg))))
+		 :value (octi::log-qd/agm3 (qd-value arg))))
 
 (rt:deftest oct.log/agm3.1
     (let* ((arg #q2)
 	   (y (log/agm3 arg))
-	   (true (make-instance 'qd-real :value qdi::+qd-log2+)))
+	   (true (make-instance 'qd-real :value octi::+qd-log2+)))
       (check-accuracy 203 y true))
   nil)
 
 (rt:deftest oct.log/agm3.2
     (let* ((arg #q10)
 	 (y (log/agm3 arg))
-	 (true (make-instance 'qd-real :value qdi::+qd-log10+)))
+	 (true (make-instance 'qd-real :value octi::+qd-log10+)))
     (check-accuracy 205 y true))
   nil)
 
@@ -474,7 +488,7 @@
 
 (defun log1p/dup (arg)
   (make-instance 'qd-real
-		 :value (qdi::log1p-qd/duplication (qd-value arg))))
+		 :value (octi::log1p-qd/duplication (qd-value arg))))
 
 (rt:deftest oct.log1p.1
     (let* ((arg #q9)
@@ -495,7 +509,7 @@
 
 (defun expm1/series (arg)
   (make-instance 'qd-real
-		 :value (qdi::expm1-qd/series (qd-value arg))))
+		 :value (octi::expm1-qd/series (qd-value arg))))
 
 (rt:deftest oct.expm1/series.1
   (let* ((arg #q0)
@@ -522,7 +536,7 @@
 
 (defun expm1/dup (arg)
   (make-instance 'qd-real
-		 :value (qdi::expm1-qd/duplication (qd-value arg))))
+		 :value (octi::expm1-qd/duplication (qd-value arg))))
 
 
 (rt:deftest oct.expm1/dup.1
@@ -551,7 +565,7 @@
 ;; thing.
 (rt:deftest oct.integer-decode.1
     (multiple-value-bind (frac exp s)
-	(qdi:integer-decode-qd (qdi::%make-qd-d -0.03980126756814893d0
+	(octi:integer-decode-qd (octi::%make-qd-d -0.03980126756814893d0
 						-2.7419792323327893d-18
 						0d0 0d0))
       (unless (and (eql frac 103329998279901916046530991816704)
--- /project/oct/cvsroot/oct/tests.lisp	2007/10/13 02:14:43	1.2
+++ /project/oct/cvsroot/oct/tests.lisp	2007/10/15 18:21:47	1.3
@@ -43,10 +43,10 @@
   (format t "bits: ~,1f~%" (bit-accuracy est true)))
 
 (defconstant +e+
-  (make-instance 'qd-real :value qdi::+qd-e+))
+  (make-instance 'qd-real :value octi::+qd-e+))
 
 (defconstant +log2+
-  (make-instance 'qd-real :value qdi::+qd-log2+))
+  (make-instance 'qd-real :value octi::+qd-log2+))
   
 (defun test2 ()
   ;; pi/4 = 4 * arctan(1/5) - arctan(1/239)
@@ -268,7 +268,7 @@
     (print-result y true))
   (let* ((arg #q10)
 	 (y (log arg))
-	 (true (make-instance 'qd-real :value qdi::+qd-log10+)))
+	 (true (make-instance 'qd-real :value octi::+qd-log10+)))
     (format t "log(10)~%")
     (print-result y true))
   (let* ((arg (+ 1 (scale-float #q1 -80)))
@@ -287,7 +287,7 @@
     (destructuring-bind (arg true)
 	f
       (let ((y (sqrt arg)))
-	(format t "sqrt(~/qdi::qd-format/)~%" (qd-value arg))
+	(format t "sqrt(~/octi::qd-format/)~%" (qd-value arg))
 	(print-result y true)))))
   
 (defun all-tests ()
--- /project/oct/cvsroot/oct/timing.lisp	2007/08/27 17:49:19	1.2
+++ /project/oct/cvsroot/oct/timing.lisp	2007/10/15 18:21:47	1.3
@@ -38,9 +38,9 @@
 	       (setf sum (cl:+ sum 1d0)))
 	     sum))
 	 (sum-%qd ()
-	   (let ((sum (qdi::make-qd-d 0d0))
-		 (one (qdi::make-qd-d 1d0)))
-	     (declare (type qdi::%quad-double sum)
+	   (let ((sum (octi::make-qd-d 0d0))
+		 (one (octi::make-qd-d 1d0)))
+	     (declare (type octi::%quad-double sum)
 		      (optimize (speed 3)))
 	     (dotimes (k n)
 	       (declare (fixnum k))
@@ -76,9 +76,9 @@
 	       (setf sum (cl:* sum 1d0)))
 	     sum))
 	 (mul-%qd ()
-	   (let ((sum (qdi::make-qd-d 0d0))
-		 (one (qdi::make-qd-d 1d0)))
-	     (declare (type qdi::%quad-double sum)
+	   (let ((sum (octi::make-qd-d 0d0))
+		 (one (octi::make-qd-d 1d0)))
+	     (declare (type octi::%quad-double sum)
 		      (optimize (speed 3)))
 	     (dotimes (k n)
 	       (declare (fixnum k))
@@ -113,9 +113,9 @@
 	       (setf sum (cl:/ sum 1d0)))
 	     sum))
 	 (div-%qd ()
-	   (let ((sum (qdi::make-qd-d 7d0))
-		 (one (qdi::make-qd-d 1d0)))
-	     (declare (type qdi::%quad-double sum)
+	   (let ((sum (octi::make-qd-d 7d0))
+		 (one (octi::make-qd-d 1d0)))
+	     (declare (type octi::%quad-double sum)
 		      (optimize (speed 3)))
 	     (dotimes (k n)
 	       (declare (fixnum k))
@@ -150,8 +150,8 @@
 	       (setf sum (cl:sqrt sum)))
 	     sum))
 	 (sqrt-%qd ()
-	   (let ((sum (qdi::make-qd-d 7d0)))
-	     (declare (type qdi::%quad-double sum)
+	   (let ((sum (octi::make-qd-d 7d0)))
+	     (declare (type octi::%quad-double sum)
 		      (optimize (speed 3)))
 	     (dotimes (k n)
 	       (declare (fixnum k))




More information about the oct-cvs mailing list