Subject: Handling of $gammalim and $factlim in simpgamma
From: Dieter Kaiser
Date: Tue, 21 Oct 2008 20:03:10 +0200
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