Changeset 6cd96249b94dc29962cecda996a5799d8ac27271
 Timestamp:
 04/08/12 09:35:35 (2 years ago)
 Author:
 Raymond Toy <toy.raymond@…>
 Children:
 dbc1e376add1d492dc35c37b3098c2ffb796d6f1
 Parents:
 c86217a5f7191f1a29566d6ee24cb57344dd546d
 gitcommitter:
 Raymond Toy <toy.raymond@…> (04/08/12 09:35:35)
 Message:

Define macro WITHFLOATINGPOINTCONTAGION.
qdmethods.lisp:
o Define the macro
qdgamma.lisp:
o Use it.
 Files:

Legend:
 Unmodified
 Added
 Removed

rc86217

r6cd962


109  109  ;; Or 
110  110  ;; log(gamma(z)) = log(pi)log(z)log(sin(pi*z))log(gamma(z)) 
111   ( (applycontagion (log pi) precision) 
112   (log ( z)) 
113   (applycontagion (log (sin (* pi z))) precision) 
114   (loggamma ( z)))) 
 111  (let ((p (floatpi z))) 
 112  ( (log p) 
 113  (log ( z)) 
 114  (log (sin (* p z))) 
 115  (loggamma ( z))))) 
115  116  (t 
116  117  (let ((absz (abs z))) 
… 
… 

378  379  
379  380  integrate(t^(a1)*exp(t), t, z, inf)" 
380   (let* ((prec (floatcontagion a z)) 
381   (a (applycontagion a prec)) 
382   (z (applycontagion z prec))) 
 381  (withfloatingpointcontagion (a z) 
383  382  (if (zerop a) 
384  383  ;; incomplete_gamma_tail(0, z) = exp_integral_e(1,z) 
… 
… 

410  409  
411  410  integrate(t^(a1)*exp(t), t, 0, z)" 
412   (let* ((prec (floatcontagion a z)) 
413   (a (applycontagion a prec)) 
414   (z (applycontagion z prec))) 
 411  (withfloatingpointcontagion (a z) 
415  412  (if (and (< (abs a) 1) (< (abs z) 1)) 
416  413  (sincompletegamma a z) 
… 
… 

540  537  (* (expt z ( v 1)) 
541  538  (incompletegammatail (+ v 1) z)))) 
542   ((< (abs z) 1) 
543   ;; Use series for small z 
 539  ((or (< (abs z) 1) (>= (abs (phase z)) 3.1)) 
 540  ;; Use series for small z or if z is near the negative real 
 541  ;; axis because the continued fraction does not converge on 
 542  ;; the negative axis and converges slowly near the negative 
 543  ;; axis. 
544  544  (sexpintegrale v z)) 
545   ((>= (abs (phase z)) 3.1) 
546   ;; The continued fraction doesn't converge on the negative 
547   ;; real axis, and converges very slowly near the negative 
548   ;; real axis, so use the incompletegammatail function in 
549   ;; this region. "Closeness" to the negative real axis is 
550   ;; teken to mean that z is in a sector near the axis. 
551   ;; 
552   ;; E(v,z) = z^(v1)*incomplete_gamma_tail(1v,z) 
553   (* (expt z ( v 1)) 
554   (incompletegammatail ( 1 v) z))) 
555  545  (t 
556  546  ;; Use continued fraction for everything else. 

r06f8a5

r6cd962


80  80  (coerce (imagpart number) precision))))) 
81  81  
 82  ;; WITHFLOATINGPOINTCONTAGION  macro 
 83  ;; 
 84  ;; Determines the highest precision of the variables in VARLIST and 
 85  ;; converts each of the values to that precision. 
 86  (defmacro withfloatingpointcontagion (varlist &body body) 
 87  (let ((precision (gensym "PRECISION"))) 
 88  `(let ((,precision (floatcontagion ,@varlist))) 
 89  (let (,@(mapcar #'(lambda (v) 
 90  `(,v (applycontagion ,v ,precision))) 
 91  varlist)) 
 92  ,@body)))) 
82  93  
83  94  (defmethod add1 ((a number)) 