* Use the g[k] formula instead of r[2*k+1] because we fail to
authorRaymond Toy <toy.raymond@gmail.com>
Thu, 12 Apr 2012 16:30:10 +0000 (09:30 -0700)
committerRaymond Toy <rtoy@google.com>
Thu, 12 Apr 2012 16:30:10 +0000 (09:30 -0700)
   converge when v = 0.
 * Clear hash tables in bessel-j.

qd-bessel.lisp

index 588b5af..fb02f75 100644 (file)
     (do* ((k 0 (1+ k))
          (bk (bk 0 p)
              (bk k p))
-         (ratio v
+         ;; Compute g[k](p)/(2*k)!, not r[2*k+1](p)/(2*k)!
+         (ratio 1
                 (* ratio (/ (+ v2 (expt (1- (* 2 k)) 2))
                             (* 2 k (1- (* 2 k))))))
          (term (* ratio bk)
            (format t " ratio = ~S~%" ratio)
            (format t " term  = ~S~%" term)
            (format t " sum   - ~S~%" sum))
-         (* sum #c(0 2) (/ (exp p) q)))
+         (* sum 4 (exp p)))
       (when *debug-exparc*
        (format t "k      = ~D~%" k)
        (format t " bk    = ~S~%" bk)
 ;; 
 (defun bessel-j (v z)
   (let ((vv (ftruncate v)))
+    ;; Clear the caches for now.
+    (an-clrhash)
+    (%big-a-clrhash)
     (cond ((= vv v)
           ;; v is an integer
           (integer-bessel-j-exp-arc v z))