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

[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  d17228b86105b8a6ac652459b1100266ad0311b3 (commit)
      from  46e7193fb5b2fad35e67366ff375d9ebbb524c5c (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 d17228b86105b8a6ac652459b1100266ad0311b3
Author: Raymond Toy <toy.raymond at gmail.com>
Date:   Thu Mar 17 22:14:25 2011 -0400

    Add Fresnel integrals; fix issue in incomplete-gamma for large args.
    
    o INCOMPLETE-GAMMA was returning bad values for large (complex)
      arguments.  Fix this by using incomplete gamma tail function since
      the incomplete gamma function approaches gamma for large arguments.
    
    o Implement Fresnel S and C functions.

diff --git a/qd-gamma.lisp b/qd-gamma.lisp
index a25bc7d..dfc87a6 100644
--- a/qd-gamma.lisp
+++ b/qd-gamma.lisp
@@ -292,7 +292,9 @@
 	(if (< z (- a 1))
 	    (cf-incomplete-gamma a z)
 	    (- (gamma a) (cf-incomplete-gamma-tail a z)))
-	(cf-incomplete-gamma a z))))
+	(if (< (abs z) (abs a))
+	    (cf-incomplete-gamma a z)
+	    (- (gamma a) (cf-incomplete-gamma-tail a z))))))
 
 (defun erf (z)
   "Error function:
@@ -341,3 +343,48 @@
   ;; Wolfram gives E(v,z) = z^(v-1)*gamma_incomplete_tail(1-v,z)
   (* (expt z (- v 1))
      (incomplete-gamma-tail (- 1 v) z)))
+
+(defun fresnel-s (z)
+  "Fresnel S:
+
+   S(z) = integrate(sin(%pi*t^2/2), t, 0, z) "
+  (let ((sqrt-pi (sqrt (float-pi z))))
+    (flet ((fs (z)
+	     ;; Wolfram gives
+	     ;;
+	     ;;  S(z) = (1+%i)/4*(erf(c*z) - %i*erf(conjugate(c)*z))
+	     ;;
+	     ;; where c = sqrt(%pi)/2*(1+%i).
+	     (* #c(1/4 1/4)
+		(- (erf (* #c(1/2 1/2) sqrt-pi z))
+		   (* #c(0 1)
+		      (erf (* #c(1/2 -1/2) sqrt-pi z)))))))
+      (if (realp z)
+	  ;; FresnelS is real for a real argument. And it is odd.
+	  (if (minusp z)
+	      (- (realpart (fs (- z))))
+	      (realpart (fs z)))
+	  (fs z)))))
+
+(defun fresnel-c (z)
+  "Fresnel C:
+
+   C(z) = integrate(cos(%pi*t^2/2), t, 0, z) "
+  (let ((sqrt-pi (sqrt (float-pi z))))
+    (flet ((fs (z)
+	     ;; Wolfram gives
+	     ;;
+	     ;;  C(z) = (1-%i)/4*(erf(c*z) + %i*erf(conjugate(c)*z))
+	     ;;
+	     ;; where c = sqrt(%pi)/2*(1+%i).
+	     (* #c(1/4 -1/4)
+		(+ (erf (* #c(1/2 1/2) sqrt-pi z))
+		   (* #c(0 1)
+		      (erf (* #c(1/2 -1/2) sqrt-pi z)))))))
+      (if (realp z)
+	  ;; FresnelS is real for a real argument. And it is odd.
+	  (if (minusp z)
+	      (- (realpart (fs (- z))))
+	      (realpart (fs z)))
+	  (fs z)))))
+     
\ No newline at end of file

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

Summary of changes:
 qd-gamma.lisp |   49 ++++++++++++++++++++++++++++++++++++++++++++++++-
 1 files changed, 48 insertions(+), 1 deletions(-)


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