Changeset 7c5a31


Ignore:
Timestamp:
04/11/12 16:18:31 (3 years ago)
Author:
Raymond Toy <rtoy@…>
Branches:
master
Children:
c0e12d, b1d9be
Parents:
bba9f8
Message:

Correct some comments, remove unused code.

File:
1 edited

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*
Note: See TracChangeset for help on using the changeset viewer.