Changeset cb1a5d41baf9d4db12a2563a230d0a1e55c8adea

Timestamp:

04/09/12 08:52:08 (14 months ago)

Author:

Raymond Toy <rtoy@…>

Children:

f8943af6bff60e23d679089db5207b4834aa83ff

Parents:

6cfb0ac4b6bcc1a25bc119e87fd2b57bfa1f4355

git-committer:

Raymond Toy <rtoy@…> (04/09/12 08:52:08)

Message:

Fix some mistakes in exp-integral-e when changing algorithm to use
series.

Files:

• qd-gamma.lisp (modified) (2 diffs)

qd-gamma.lisp

@@ -533 +533 @@
   ;;
   ;;
-  (let* ((prec (float-contagion v z))
-         (v (apply-contagion v prec))
-         (z (apply-contagion z prec)))
+  (with-floating-point-contagion (v z)
     (cond ((and (realp v) (minusp v))
            ;; E(-v, z) = z^(-v-1)*incomplete_gamma_tail(v+1,z)
@@ -541 +539 @@
              (* (expt z (- v 1))
                 (incomplete-gamma-tail (+ -v 1) z))))
-          ((< (abs z) 1)
-           ;; Use series for small z
+          ((or (< (abs z) 1) (>= (abs (phase z)) 3.1))
+           ;; Use series for small z or if z is near the negative real
+           ;; axis because the continued fraction does not converge on
+           ;; the negative axis and converges slowly near the negative
+           ;; axis.
            (s-exp-integral-e v z))
-          ((>= (abs (phase z)) 3.1)
-           ;; The continued fraction doesn't converge on the negative
-           ;; real axis, and converges very slowly near the negative
-           ;; real axis, so use the incomplete-gamma-tail function in
-           ;; this region.  "Closeness" to the negative real axis is
-           ;; teken to mean that z is in a sector near the axis.
-           ;;
-           ;; E(v,z) = z^(v-1)*incomplete_gamma_tail(1-v,z)
-           (* (expt z (- v 1))
-              (incomplete-gamma-tail (+ -v 1) z))))
-        ((or (< (abs z) 1) (>= (abs (phase z)) 3.1))
-         ;; Use series for small z or if z is near the negative real
-         ;; axis because the continued fraction does not converge on
-         ;; the negative axis and converges slowly near the negative
-         ;; axis.
-         (s-exp-integral-e v z))
-        (t
-         ;; Use continued fraction for everything else.
-         (cf-exp-integral-e v z))))
+          (t
+           ;; Use continued fraction for everything else.
+           (cf-exp-integral-e v z)))))
 
 ;; Series for Fresnel S