Add flo:signed-lgamma.

diff --git a/doc/ref-manual/numbers.texi b/doc/ref-manual/numbers.texi
index fbd8b9f6608e1fa71231cbe96b748fc165ec8e7a..4975258bb29f1f0b816ded882c564eed3ce52e7b 100644 (file)
--- a/doc/ref-manual/numbers.texi
+++ b/doc/ref-manual/numbers.texi
@@ -2036,6 +2036,41 @@ This is the flonum version of @code{atan} with two arguments.  When
 compiled, it does not check the types of its arguments.
 @end deffn
 
+@deffn procedure flo:signed-lgamma x
+Returns two values,
+@iftex
+@tex
+$$m = \log \left|\Gamma(x)\right| \quad\hbox{and}\quad s = \mathop{\rm sign} \Gamma(x),$$
+@end tex
+@end iftex
+@ifnottex
+
+@example
+@group
+m = log(|Gamma(@var{x})|)
+
+and
+
+s = sign(Gamma(@var{x})),
+@end group
+@end example
+
+@end ifnottex
+respectively a flonum and an exact integer either @code{-1} or
+@code{1}, so that
+@iftex
+@tex
+$$\Gamma(x) = s \cdot e^m.$$
+@end tex
+@end iftex
+@ifnottex
+
+@example
+Gamma(x) = s * e^m.
+@end example
+
+@end ifnottex
+@end deffn
 
 @deffn procedure flo:min x1 x2
 @deffnx procedure flo:max x1 x2

diff --git a/src/microcode/flonum.c b/src/microcode/flonum.c
index 806c4f6978b77843a664b155ee35c72013143f67..b8c43fe0e36f859acfb85ab984078af05e8896e6 100644 (file)
--- a/src/microcode/flonum.c
+++ b/src/microcode/flonum.c
@@ -334,6 +334,26 @@ DEFINE_PRIMITIVE ("FLONUM-LGAMMA", Prim_flonum_lgamma, 1, 1, 0)
   FLONUM_RESULT (result);
 }
 
+DEFINE_PRIMITIVE ("FLONUM-SIGNED-LGAMMA", Prim_flonum_signed_lgamma, 1, 1, 0)
+{
+  double x;
+  double result;
+  int sign;
+  PRIMITIVE_HEADER (1);
+
+  x = (arg_flonum (1));
+#ifdef HAVE_LGAMMA_R
+  result = (lgamma_r (x, (&sign)));
+#else
+  result = (lgamma (x));
+  sign = signgam;
+#endif
+
+  assert (LONG_TO_FIXNUM_P (sign));
+  PRIMITIVE_RETURN
+    (cons ((double_to_flonum (result)), (LONG_TO_FIXNUM (sign))));
+}
+
 DEFINE_PRIMITIVE ("FLONUM-GAMMA", Prim_flonum_gamma, 1, 1, 0)
      SIMPLE_TRANSCENDENTAL_FUNCTION (tgamma)
 DEFINE_PRIMITIVE ("FLONUM-ERF", Prim_flonum_erf, 1, 1, 0)

diff --git a/src/runtime/primitive-arithmetic.scm b/src/runtime/primitive-arithmetic.scm
index 0c5a7c4fb3d22d1799c25cce08a2a04b3567d98e..b535bf1960e50d3dc90e3c29317c1ffca35af9af 100644 (file)
--- a/src/runtime/primitive-arithmetic.scm
+++ b/src/runtime/primitive-arithmetic.scm
@@ -274,7 +274,7 @@ USA.
 
 (define (flo:total-order-mag x y)
   (flo:total-order (flo:abs x) (flo:abs y)))
-
+\f
 (define (flo:total< x y)
   (if (or (flo:nan? x) (flo:nan? y))
       ;; Must handle NaNs first and carefully to avoid exception on
@@ -420,6 +420,10 @@ USA.
   (if (and (flo:finite? x) (not (flo:safe-zero? x)))
       (fix:- (cdr ((ucode-primitive flonum-normalize 1) x)) 1)
       (begin (flo:raise-exceptions! (flo:exception:invalid-operation)) #f)))
+
+(define (flo:signed-lgamma x)
+  (let ((p ((ucode-primitive flonum-signed-lgamma 1) x)))
+    (values (car p) (cdr p))))
 \f
 ;;;; Exact integers
 

diff --git a/src/runtime/runtime.pkg b/src/runtime/runtime.pkg
index 145ed4b582e1dcd1ee684a0ccbbc9a1003935b82..acb35ec2875236434fc211cd15d1f2bcae39f56a 100644 (file)
--- a/src/runtime/runtime.pkg
+++ b/src/runtime/runtime.pkg
@@ -359,6 +359,7 @@ USA.
          flo:safe>
          flo:safe>=
          flo:sign-negative?
+         flo:signed-lgamma
          flo:sin
          flo:sinh
          flo:snan

diff --git a/tests/runtime/test-arith.scm b/tests/runtime/test-arith.scm
index 712bc68fb5f72c1d48ba37ab630ba0c39714ae2e..5cc1f23707011532b7d97031319003148db4ef79 100644 (file)
--- a/tests/runtime/test-arith.scm
+++ b/tests/runtime/test-arith.scm
@@ -1097,4 +1097,12 @@ USA.
         (assert-eqv (exact->inexact x) y)
         (if (not (= x y))
             (assert-except/no-traps (flo:exception:inexact-result)
-                                    (lambda () (exact->inexact x))))))))
\ No newline at end of file
+                                    (lambda () (exact->inexact x))))))))
+
+(define-enumerated-test 'flo:lgamma
+  (list (list -0.123))
+  (lambda (x)
+    (receive (log-gamma sign) (flo:signed-lgamma x)
+      (assert-eqv (flo:lgamma x) log-gamma)
+      (let ((gamma (* sign (exp log-gamma))))
+       (assert-<= (relerr (flo:gamma x) gamma) 1e-15)))))
\ No newline at end of file