Handling of $gammalim and $factlim in simpgamma



The value of $gammalim is by default set to 1,000,000

I have not tried to finish the calculation gamma((gammalim-1)+1/2) with this
quit big value of $gammalim. But for $gammalim set to 100,000 we get a result in
about 60 seconds with GCL and a 2,2 GHZ CPU.

Second, $gammalim and $factlim are not independent because of the algorithm in
simpgamma. Thus e.g. for $gammalim > $factlim we get gamma (gammlim) =
(gammlim-1)! and not the noun form for gamma itself.

Another point is that the handling of negative integers is a bit confusing. We
get an error for negative integers, but only for absolute values which are
smaller than $gammalim. So gamma(-gammalim) gives an error, but
gamma(-gammalim-1) gives a noun form. Thus $gammalim not only controls the
expansion for rational numbers but the handling of negative integers.

I would like to suggest the following changes:

1. Move definition of $gammalim from csimp.lisp to csimp2.lisp.
   That is the only file we need the definition.
2. Set $gammalim to a more moderate value 100,000.
3. Remove an early test in simpgamma against $gammalim.
   With this change we can use $factlim and $gammalim independently.
4. Extending the test for integer values: Add a test against $factlim
   and return a noun form for positve integers greater than $factlim.
   Now we get an error for all negative integers, independent from the 
   value of $gammalim.
5. Extending the test for rational numbers: Check for $gammalim.
   Now $gammalim only controls the expansion for rational numbers.

These are the changes:

Index: csimp2.lisp
===================================================================
RCS file: /cvsroot/maxima/maxima/src/csimp2.lisp,v
retrieving revision 1.30
diff -u -r1.30 csimp2.lisp
--- csimp2.lisp   21 Oct 2008 12:09:25 -0000   1.30
+++ csimp2.lisp   21 Oct 2008 16:49:21 -0000
@@ -16,6 +16,9 @@
 
 (declare-top (special var %p%i varlist plogabs half%pi nn* dn*))
 
+(defmvar $gammalim 10000
+  "Controls simplification of gamma for rational number arguments.")
+
 (defvar $gamma_expand nil
   "Expand gamma(z+n) for n an integer when T.") 
 
@@ -225,17 +228,27 @@
                      ;; and not (1-z)*(2-z)*... 
                      ($factor
                        (simplify (list '($pochhammer) (sub 1 z) n))))))))
-     ((or (not (mnump j))
-          (ratgreaterp (simplify (list '(mabs) j)) $gammalim))
-      (eqtest (list '(%gamma) j) x))
      ((integerp j)
-      (cond ((> j 0) (simplify (list '(mfactorial) (1- j))))
+      (cond ((> j 0)
+             (cond ((<= j $factlim)
+                    ;; Positive integer less than $factlim. Evaluate.
+                    (simplify (list '(mfactorial) (1- j))))
+                   ;; Positive integer greater $factlim. Noun form.
+                   (t eqtest (list '(%gamma) j) x)))
+            ;; Negative integer. Throw a Maxima error.
             (errorsw (throw 'errorsw t))
             (t (merror "gamma(~:M) is undefined" j))))
      ($numer (gammafloat (fpcofrat j)))
      ((alike1 j '((rat) 1 2))
       (list '(mexpt simp) '$%pi j))
-     ((or (ratgreaterp j 1) (ratgreaterp 0 j)) (gammared j))
+     ((and (mnump j)
+           (ratgreaterp $gammalim (simplify (list '(mabs) j)))
+           (or (ratgreaterp j 1) (ratgreaterp 0 j)))
+      ;; Expand for rational number less than $gammalim.
+      (gammared j))
      (t (eqtest (list '(%gamma) j) x)))))

Success, CVS operation completed

So what do you think? Should we commit these changes?

Dieter Kaiser