Zsh Mailing List Archive
Messages sorted by: Reverse Date, Date, Thread, Author

[PATCH] Fix gamma() in mathfunc module.



Use tgamma() if available. Otherwise use signgam*exp(lgamma).
A simple test is added. A typo in mathfunc.yo is also fixed.
---

This patch makes gamma() in mathfunc module always return Gamma(x), not
log(abs(Gamma(x))).

In the module zsh/mathfunc, there is a function gamma() to calculate the
Gamma function. I guess little people would use Gamma function in zsh,
but anyway its behavior is different between BSD and non-BSD systems.

On Linux (and most of non-BSD systems?), C functions gamma(x) and
lgamma(x) (in libm.a) both return log(abs(Gamma(x))). As a result, in zsh,
both 'echo $(( gamma(2) ))' and 'echo $(( lgamma(2) ))' return '0.'
(Gamma(2) = 1, log(1) = 0). I guess this is not the intended result.

On Mac OS X (and most of BSD-based systems?), gamma(x) returns Gamma(x)
while lgamma(x) returns log(abs(Gamma(x))). So 'echo $(( gamma(2) ))'
gives '1.' ( = Gamma(2) ) as expected.

On both Linux and OS X, the recommended way to calculate Gamma(x)
is to use C function tgamma(x). (On Mac OS X 10.9 (Mavericks), I get a
warning that gamma() is obsolete and should be replaced by tgamma().)


 Doc/Zsh/mod_mathfunc.yo | 2 +-
 Src/Modules/mathfunc.c  | 9 +++++++++
 Test/V03mathfunc.ztst   | 6 ++++++
 configure.ac            | 2 +-
 4 files changed, 17 insertions(+), 2 deletions(-)

diff --git a/Doc/Zsh/mod_mathfunc.yo b/Doc/Zsh/mod_mathfunc.yo
index dda4f36..5239da5 100644
--- a/Doc/Zsh/mod_mathfunc.yo
+++ b/Doc/Zsh/mod_mathfunc.yo
@@ -33,7 +33,7 @@ returns an integer.
 The function tt(signgam) takes no arguments, and returns an integer, which
 is the C variable of the same name, as described in manref(gamma)(3).  Note
 that it is therefore only useful immediately after a call to tt(gamma) or
-tt(lgamma).  Note also that `tt(signgam+LPAR()RPAR)' and `tt(signgam)' are
+tt(lgamma).  Note also that `tt(signgam+LPAR()RPAR())' and `tt(signgam)' are
 distinct expressions.
 
 The following functions take two floating point arguments: tt(copysign),
diff --git a/Src/Modules/mathfunc.c b/Src/Modules/mathfunc.c
index 04483b5..efadd86 100644
--- a/Src/Modules/mathfunc.c
+++ b/Src/Modules/mathfunc.c
@@ -340,7 +340,16 @@ math_func(char *name, int argc, mnumber *argv, int id)
       break;
 
   case MF_GAMMA:
+#ifdef HAVE_TGAMMA
+      retd = tgamma(argd);
+#else
+#ifdef HAVE_SIGNGAM
+      retd = lgamma(argd);
+      retd = signgam*exp(retd);
+#else /*XXX assume gamma(x) returns Gamma(x), not log(|Gamma(x)|) */
       retd = gamma(argd);
+#endif
+#endif
       break;
 
   case MF_HYPOT:
diff --git a/Test/V03mathfunc.ztst b/Test/V03mathfunc.ztst
index 069059d..ab383db 100644
--- a/Test/V03mathfunc.ztst
+++ b/Test/V03mathfunc.ztst
@@ -136,3 +136,9 @@ F:This test fails if your math library doesn't have erand48().
   done
   (( ok ))
 0:Test random number generator distributions are not grossly broken
+
+  float -F 5 g l
+  (( g = gamma(2), l = lgamma(2) ))
+  print $g, $l
+0:Test Gamma function gamma and lgamma
+>1.00000, 0.00000
diff --git a/configure.ac b/configure.ac
index c3debd8..c3093f2 100644
--- a/configure.ac
+++ b/configure.ac
@@ -1222,7 +1222,7 @@ AC_CHECK_FUNCS(strftime strptime mktime timelocal \
 	       getrlimit getrusage \
 	       setlocale \
 	       uname \
-	       signgam \
+	       signgam tgamma \
 	       putenv getenv setenv unsetenv xw\
 	       brk sbrk \
 	       pathconf sysconf \
-- 
1.8.3.4 (Apple Git-47)




Messages sorted by: Reverse Date, Date, Thread, Author