[oct-scm] [oct-git]OCT: A portable Lisp implementation for quad-double precision floats branch master updated. f203389e5e78d1f001f68447ac2a9dd86dcfbbf6

[Raymond Toy]

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 "OCT:  A portable Lisp implementation for quad-double precision floats".

The branch, master has been updated
       via  f203389e5e78d1f001f68447ac2a9dd86dcfbbf6 (commit)
      from  f8943af6bff60e23d679089db5207b4834aa83ff (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 f203389e5e78d1f001f68447ac2a9dd86dcfbbf6
Author: Raymond Toy <toy.raymond at gmail.com>
Date:   Tue Apr 10 00:45:13 2012 -0700

     * Implement sum-big-ia.
     * Add series for Bessel J.  (Not working yet.)

diff --git a/qd-bessel.lisp b/qd-bessel.lisp
index fbe9760..9759044 100644
--- a/qd-bessel.lisp
+++ b/qd-bessel.lisp
@@ -324,8 +324,41 @@
        (incomplete-gamma-tail a (* theta z)))))
 
 (defun sum-big-ia (big-n v z)
-  )
+  (let ((big-n-1/2 (+ big-n 1/2))
+	(eps (epsilon z)))
+    (do* ((n 0 (1+ n))
+	  (term (* (big-a 0 v)
+		   (big-i 0 big-n-1/2 z v))
+		(* (big-a n v)
+		   (big-i n big-n-1/2 z v)))
+	  (sum term (+ sum term)))
+	 ((<= (abs term) (* eps (abs sum)))
+	  sum)
+      #+nil
+      (progn
+	(format t "n = ~D~%" n)
+	(format t " term = ~S~%" term)
+	(format t " sum  = ~S~%" sum)))))
 
+;; Series for bessel J:
+;;
+;; (z/2)^v*sum((-1)^k/Gamma(k+v+1)/k!*(z^2//4)^k, k, 0, inf)
+(defun s-bessel-j (v z)
+  (with-floating-point-contagion (v z)
+    (let ((z2/4 (* z z 1/4))
+	  (eps (epsilon z)))
+      (do* ((k 0 (+ 1 k))
+	    (f (gamma (+ v 1))
+	       (* f (* k (+ v k))))
+	    (term (/ f)
+		  (/ (* (- term) z2/4) f))
+	    (sum term (+ sum term)))
+	   ((<= (abs term) (* eps (abs sum)))
+	    (* sum (expt (* z 1/2) v)))
+	(format t "k = ~D~%" k)
+	(format t " term = ~S~%" term)
+	(format t " sum  = ~S~%" sum)))))
+  
 (defun bessel-j (v z)
   (let ((vv (ftruncate v)))
     (cond ((= vv v)
@@ -338,7 +371,8 @@
 		(* z
 		   (/ (sin vpi) vpi)
 		   (+ (/ -1 z)
-		      (sum-ab big-n v z)))))))))
+		      (sum-ab big-n v z)
+		      (sum-big-ia big-n v z)))))))))
 
 (defun paris-series (v z n)
   (labels ((pochhammer (a k)

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

Summary of changes:
 qd-bessel.lisp |   38 ++++++++++++++++++++++++++++++++++++--
 1 files changed, 36 insertions(+), 2 deletions(-)


hooks/post-receive
-- 
OCT:  A portable Lisp implementation for quad-double precision floats