[Maxima-commits] Maxima, A Computer Algebra System branch master updated. 5.28.0-128-gdd921f0

[Raymond T. <rt...@us...>]

This is an automated email from the git hooks/post-receive script. It was
generated because a ref change was pushed to the repository containing
the project "Maxima, A Computer Algebra System".

The branch, master has been updated
       via  dd921f09514a8f86b700df2b821341967e17dc7a (commit)
       via  801a2951d5ee373e738ffc505952d1c51bcf1875 (commit)
      from  ec732db182bc9ecf208c1f2906dfbb6af46f3553 (commit)

Those revisions listed above that are new to this repository have
not appeared on any other notification email; so we list those
revisions in full, below.

- Log -----------------------------------------------------------------
commit dd921f09514a8f86b700df2b821341967e17dc7a
Author: Raymond Toy <toy...@gm...>
Date:   Fri Nov 16 21:21:53 2012 -0800

    o Fix warning about unused EPS variable
    o Use mformat instead of format to print out debugging info.

diff --git a/src/gamma.lisp b/src/gamma.lisp
index 05d1ff4..ab2ec17 100755
--- a/src/gamma.lisp
+++ b/src/gamma.lisp
@@ -946,8 +946,8 @@
 		(merror (intl:gettext "gamma_incomplete: continued fractions failed for gamma_incomplete(~:M, ~:M).") a x))
 	    (when *debug-gamma* 
 	      (format t "~&in coninued fractions:~%")
-	      (format t "~&   : i = ~A~%" i)
-	      (format t "~&   : h = ~A~%" h))
+	      (mformat t "~&   : i = ~M~%" i)
+	      (mformat t "~&   : h = ~M~%" h))
 	    (setq d (add (mul an d) b))
 	    (when (eq ($sign (sub (simplify (list '(mabs) d)) gm-min)) '$neg)
 	      (setq d gm-min))
@@ -997,15 +997,15 @@
              (merror (intl:gettext "gamma_incomplete: series expansion failed for gamma_incomplete(~:M, ~:M).") a x))
          (when *debug-gamma* 
            (format t "~&GAMMA-INCOMPLETE in series:~%")
-           (format t "~&   : i    = ~A~%" i)
-           (format t "~&   : ap   = ~A~%" ap)
-           (format t "~&   : x/ap = ~A~%" (div x ap))
-           (format t "~&   : del  = ~A~%" del)
-           (format t "~&   : sum  = ~A~%" sum))
+           (mformat t "~&   : i    = ~M~%" i)
+           (mformat t "~&   : ap   = ~M~%" ap)
+           (mformat t "~&   : x/ap = ~M~%" (div x ap))
+           (mformat t "~&   : del  = ~M~%" del)
+           (mformat t "~&   : sum  = ~M~%" sum))
          (when (eq ($sign (sub (simplify (list '(mabs) del)) 
                                (mul (simplify (list '(mabs) sum)) gm-eps)))
                    '$neg)
-           (when *debug-gamma* (format t "~&Series converged to ~A.~%" sum))
+           (when *debug-gamma* (mformat t "~&Series converged to ~M.~%" sum))
            (return 
              (sub (simplify (list '(%gamma) a))
                   ($rectform
@@ -2018,8 +2018,7 @@
 	 ;; The threshold is approximately erf(x) = 0.5.  (Doesn't
 	 ;; have to be super accurate.)  This gives max accuracy when
 	 ;; using the identity  erf(x) = 1 - erfc(x).
-	 (let ((z^2 (* z z))
-	       (eps (epsilon z)))
+	 (let ((z^2 (* z z)))
 	   (/ (* 2 z (sum-power-series z^2
 				       #'(lambda (k)
 					   (let ((2k (+ k k)))

commit 801a2951d5ee373e738ffc505952d1c51bcf1875
Author: Raymond Toy <toy...@gm...>
Date:   Fri Nov 16 21:07:23 2012 -0800

    o EXPT-EXTRA-BITS Needs to take into account both the base and the
      power.
    o Clean up implementation a bit.

diff --git a/src/numeric.lisp b/src/numeric.lisp
index bcb985d..a90a952 100755
--- a/src/numeric.lisp
+++ b/src/numeric.lisp
@@ -154,38 +154,65 @@
 (defmethod maxima::to ((z t))
   z)
 
-
-(defun expt-extra-bits (x a)
-  ;; When x^a using exp(a*log(x)), we need extra bits because the
-  ;; integer part of a*log(x) doesn't contribute to the accuracy of
-  ;; the result.  If x = 2^n*f, where |f| < 1, log2(x) = n + log(f).
-  ;; If a = 2^m*g, where |g| < 1, then a*log2(x) = a*n + a*log2(f) =
-  ;; a*n + 2^m*g*log2(f).
-  (declare (ignore x))
-  (let ((rp (realpart a)))
-    (if (typep rp '(or bigfloat complex-bigfloat))
-	(max 1 (third (slot-value rp 'real)))
-	(integer-length (truncate rp)))))
-
-;; Executes the body BODY with extra precision.  The precision is
-;; increased by EXTRA, and the list of variables given in VARLIST have
-;; the precision increased.  The precision of the first value of the
-;; body is then reduced back to the normal precision.
+;; MAX-EXPONENT roughly computes the log2(|x|).  If x is real and x =
+;; 2^n*f, with |f| < 1, MAX-EXPONENT returns |n|.  For complex
+;; numbers, we return one more than the max of the exponent of the
+;; real and imaginary parts.
+(defmethod max-exponent ((x bigfloat))
+  (cl:abs (third (slot-value x 'real))))
+
+(defmethod max-exponent ((x complex-bigfloat))
+  (cl:1+ (cl:max (cl:abs (third (slot-value x 'real)))
+		 (cl:abs (third (slot-value x 'imag))))))
+
+(defmethod max-exponent ((x cl:float))
+  (cl:abs (nth-value 1 (cl:decode-float x))))
+
+(defmethod max-exponent ((x cl:rational))
+  (cl:ceiling (cl:log (cl:abs x) 2)))
+
+(defmethod max-exponent ((x cl:complex))
+  (cl:1+ (cl:max (max-exponent (cl:realpart x))
+		 (max-exponent (cl:imagpart x)))))
+
+;; When computing x^a using exp(a*log(x)), we need extra bits because
+;; the integer part of a*log(x) doesn't contribute to the accuracy of
+;; the result.  The number of extra bits needed is basically the
+;; "size" of a plus the number of bits for ceiling(log(x)).  We need
+;; ceiling(log(x)) extra bits because that's how many bits are taken
+;; up by the log(x).  The "size" of a is, basically, the exponent of
+;; a. If a = 2^n*f where |f| < 1, then the size is abs(n) because
+;; that's how many extra bits are added to the integer part of
+;; a*log(x).
+(defmethod expt-extra-bits ((x t) (a t))
+  (max 1 (+ (integer-length (max-exponent x))
+	    (max-exponent a))))
+
+;;; WITH-EXTRA-PRECISION - Internal
+;;;
+;;;   Executes the body BODY with extra precision.  The precision is
+;;; increased by EXTRA, and the list of variables given in VARLIST have
+;;; the precision increased.  The precision of the first value of the
+;;; body is then reduced back to the normal precision.
 (defmacro with-extra-precision ((extra (&rest varlist)) &body body)
-  (let ((result (gensym)))
+  (let ((result (gensym))
+	(old-fpprec (gensym)))
     `(let ((,result
-	     (let ((maxima::fpprec (cl:+ maxima::fpprec ,extra)))
-	       (let ,(mapcar #'(lambda (v)
-				 ;; Could probably do this in a faster
-				 ;; way, but conversion to a maxima
-				 ;; form autoamatically increases the
-				 ;; precision of the bigfloat to the
-				 ;; new precision.  Conversion of that
-				 ;; to a bigfloat object preserves the
-				 ;; precision.
-				 `(,v (bigfloat:to (maxima::to ,v))))
-		      varlist)
-		 ,@body))))
+	     (let ((,old-fpprec maxima::fpprec))
+	       (unwind-protect
+		    (let ((maxima::fpprec (cl:+ maxima::fpprec ,extra)))
+		      (let ,(mapcar #'(lambda (v)
+					;; Could probably do this in a faster
+					;; way, but conversion to a maxima
+					;; form automatically increases the
+					;; precision of the bigfloat to the
+					;; new precision.  Conversion of that
+					;; to a bigfloat object preserves the
+					;; precision.
+					`(,v (bigfloat:to (maxima::to ,v))))
+			     varlist)
+			,@body))
+		 (setf maxima::fpprec ,old-fpprec)))))
        ;; Conversion of the result to a maxima number adjusts the
        ;; precision appropriately.
        (bigfloat:to (maxima::to ,result)))))

-----------------------------------------------------------------------

Summary of changes:
 src/gamma.lisp   |   19 +++++------
 src/numeric.lisp |   87 +++++++++++++++++++++++++++++++++++------------------
 2 files changed, 66 insertions(+), 40 deletions(-)


hooks/post-receive
-- 
Maxima, A Computer Algebra System