Show
Ignore:
Timestamp:
04/11/12 09:18:31 (2 years ago)
Author:
Raymond Toy <rtoy@…>
Children:
c0e12ddf6b61f571555d46c6f168e6bebabc80b1, 25129063aaa6168763b447db6f033359a19b60e0
Parents:
bba9f8940c9f904bf14adc405d795a38ac333c24
git-committer:
Raymond Toy <rtoy@…> (04/11/12 09:18:31)
Message:

Correct some comments, remove unused code.

Files:
1 modified

Legend:

Unmodified
Added
Removed
  • qd-bessel.lisp

    re9cd1a r7c5a31  
    3838;; B[k](p) = 1/2^(k+3/2)*integrate(exp(-p*u)*u^(k-1/2),u,0,1) 
    3939;;         = 1/2^(k+3/2)/p^(k+1/2)*integrate(t^(k-1/2)*exp(-t),t,0,p) 
    40 ;;         = 1/2^(k+3/2)/p^(k+1/2) * g(k+1/2, p) 
     40;;         = 1/2^(k+3/2)/p^(k+1/2) * G(k+1/2, p) 
    4141;; 
    4242;; where G(a,z) is the lower incomplete gamma function. 
     
    110110;; 
    111111;; Hence I(-%i*z, v) = conj(I(%i*z, v)) when both z and v are real. 
     112;; 
     113;; Also note that when v is an integer of the form (2*m+1)/2, then 
     114;;   r[2*k+1](-2*%i*v) = r[2*k+1](-%i*(2*m+1)) 
     115;;                     = -%i*(2*m+1)*product(-(2*m+1)^2+(2*j-1)^2, j, 1, k) 
     116;; so the product is zero when k >= m and the series I(p, q) is 
     117;; finite. 
    112118(defun exp-arc-i (p q) 
    113119  (let* ((sqrt2 (sqrt (float 2 (realpart p)))) 
     
    145151 
    146152(defun exp-arc-i-2 (p q) 
    147   (let* ((sqrt2 (sqrt (float 2 (realpart p)))) 
    148          (exp/p/sqrt2 (/ (exp (- p)) p sqrt2)) 
    149          (v (* #c(0 -2) q)) 
     153  (let* ((v (* #c(0 -2) q)) 
    150154         (v2 (expt v 2)) 
    151155         (eps (epsilon (realpart p)))) 
    152     (when *debug-exparc* 
    153       (format t "sqrt2 = ~S~%" sqrt2) 
    154       (format t "exp/p/sqrt2 = ~S~%" exp/p/sqrt2)) 
    155156    (do* ((k 0 (1+ k)) 
    156157          (bk (bk 0 p) 
     
    163164          (sum term (+ sum term))) 
    164165         ((< (abs term) (* (abs sum) eps)) 
     166          (when *debug-exparc* 
     167            (format t "Final k= ~D~%" k) 
     168            (format t " bk    = ~S~%" bk) 
     169            (format t " ratio = ~S~%" ratio) 
     170            (format t " term  = ~S~%" term) 
     171            (format t " sum   - ~S~%" sum)) 
    165172          (* sum #c(0 2) (/ (exp p) q))) 
    166173      (when *debug-exparc*