applying full_simplify() to gamma functions causes an error

[sagetrac-tomc]

Applying full_simplify() to the gamma function sometimes causes an error. This example works:

sage: gamma(4/3).full_simplify()
1/3*gamma(1/3)

but this example does not:

sage: gamma(1/3).full_simplify()
ERROR: An unexpected error occurred while tokenizing input
The following traceback may be corrupted or invalid
The error message is: ('EOF in multi-line statement', (1254, 0))

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)

/Users/tomc/sage-4.4.1/<ipython console> in <module>()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/symbolic/expression.so in sage.symbolic.expression.Expression.simplify_full (sage/symbolic/expression.cpp:21549)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/symbolic/expression.so in sage.symbolic.expression.Expression.simplify_factorial (sage/symbolic/expression.cpp:22240)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/structure/parent.so in sage.structure.parent.Parent.__call__ (sage/structure/parent.c:6332)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/structure/coerce_maps.so in sage.structure.coerce_maps.NamedConvertMap._call_ (sage/structure/coerce_maps.c:4053)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/interfaces/maxima.pyc in _symbolic_(self, R)
   1810             sqrt(2)
   1811         """
-> 1812         return R(self._sage_())
   1813 
   1814     def __complex__(self):

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/interfaces/maxima.pyc in _sage_(self)
   1791         import sage.calculus.calculus as calculus
   1792         return calculus.symbolic_expression_from_maxima_string(self.name(),
-> 1793                 maxima=self.parent())
   1794 
   1795     def _symbolic_(self, R):

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/calculus/calculus.pyc in symbolic_expression_from_maxima_string(x, equals_sub, maxima)
   1524         # evaluation of maxima code are assumed pre-simplified
   1525         is_simplified = True
-> 1526         return symbolic_expression_from_string(s, syms, accept_sequence=True)
   1527     except SyntaxError:
   1528         raise TypeError, "unable to make sense of Maxima expression '%s' in Sage"%s

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/calculus/calculus.pyc in symbolic_expression_from_string(s, syms, accept_sequence)
   1692             global _augmented_syms
   1693             _augmented_syms = syms
-> 1694             return parse_func(s)
   1695         finally:
   1696             _augmented_syms = {}

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.parse_sequence (sage/misc/parser.c:3855)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.parse_sequence (sage/misc/parser.c:3747)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.p_sequence (sage/misc/parser.c:4376)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.p_tuple (sage/misc/parser.c:5032)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.p_eqn (sage/misc/parser.c:5145)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.p_expr (sage/misc/parser.c:5465)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.p_term (sage/misc/parser.c:5690)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.p_factor (sage/misc/parser.c:6053)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/misc/parser.so in sage.misc.parser.Parser.p_power (sage/misc/parser.c:6264)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/symbolic/function.so in sage.symbolic.function.GinacFunction.__call__ (sage/symbolic/function.cpp:6321)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/symbolic/expression.so in sage.symbolic.expression.Expression.factorial (sage/symbolic/expression.cpp:20595)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/symbolic/pynac.so in sage.symbolic.pynac.py_factorial (sage/symbolic/pynac.cpp:9156)()

/Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/rings/arith.pyc in factorial(n, algorithm)
    403     """
    404     if n < 0:
--> 405         raise ValueError, "factorial -- must be nonnegative"
    406     if algorithm == 'gmp':
    407         return ZZ(n).factorial()

ValueError: factorial -- must be nonnegative

I am running Sage 4.4.1 on Mac OS X version 10.6 (Snow Leopard), built from source. But the second example also fails on Sage 4.3.5 on 64-bit Linux, built from source. Looking at the source code suggests that the second example will fail on all platforms.

The problem occurs because full_simplify() here runs the following commands in Maxima:

(%i1) minfactorial(factcomb(makefact(gamma(1/3))));
                                       2
(%o1)                               (- -)!
                                       3

and then the Maxima interface converts this to Sage as factorial(-2/3). This causes an error. For Sage, factorial(x) is only defined if x is a non-negative integer, whereas for Maxima factorial(x) is equivalent to gamma(1+x) and so makes sense whenever x is not in {-1, -2, -3, ...}

Apply: trac_9240_full_simplify_gamma.patch, trac_9240-factorial_evaluation.patch

Depends on #11415

CC: @kcrisman

Component: symbolics

Keywords: gamma function, full_simplify, factorial

Author: Tom Coates, Burcin Erocal

Reviewer: Dan Drake, Karl-Dieter Crisman, François Bissey

Merged: sage-4.7.1.alpha4

Issue created by migration from https://trac.sagemath.org/ticket/9240

[sagetrac-tomc]

Attachment: trac_9240_full_simplify_gamma.patch.gz

based on Sage 4.4.1

[dandrake]

comment:1

The patch applies cleanly to 4.4.4.alpha0, fixes the problem and includes some nice tests. I would give this a positive review, except that I really don't know the Pynac integration code.

Oh, and you say "For Sage, factorial(x) is only defined if x is a non-negative integer", but (with a clean 4.4.4.alpha0), it's defined the same way as Maxima:

sage: factorial(5)
120
sage: factorial(1/2)
1/2*sqrt(pi)
sage: factorial(1/21)
gamma(22/21)

[sagetrac-tomc]

comment:2

OK: I had misunderstood the docstring, which says:

sage: ? factorial

String Form:    factorial
Namespace:      Interactive
File:           /Users/tomc/sage-4.4.1/local/lib/python2.6/site-packages/sage/functions/other.py
Definition:     factorial(self, *args, coerce=True, hold=False, dont_call_method_on_arg=False)
Docstring:
       Returns the factorial of n.
    
       INPUT:
    
       * ``n`` - an integer, or symbolic expression
...

I suppose that non-integer numerical values (rational numbers, real numbers, etc) count as symbolic expressions here, as they can be canonically coerced into the symbolic ring. But then there is definitely a bug in factorial(), because [in an unpatched version of Sage]:

sage: factorial(-1/2)
ERROR: An unexpected error occurred while tokenizing input
...
ValueError: factorial -- must be nonnegative

The patch fixes this.

[dandrake]

comment:3

Replying to @sagetrac-tomc:

OK: I had misunderstood the docstring,

Well, of course you misunderstood the docstring, since it's wrong, or at least misleading. I've opened #9248 to fix the docstring for factorial().

[kcrisman]

comment:4

See also this thread, where it came up again.

[burcin]

comment:5

I suppose this needs review.

[dandrake]

comment:6

For some reason, the patchbot isn't applying this patch, but it does apply cleanly and pass doctests in 4.7.rc1. (On 64-bit Ubuntu 10.10)

Looking at it again, I might recommend that in the TESTS, the first few examples be moved, since they already work and aren't testing to see if the problem in this ticket has been fixed. Only the full_simplify() tests are, strictly speaking, relevant to this ticket.