Make log gamma symbolic #10075

kcrisman · 2010-10-05T13:59:32Z

Currently, there is no way to send log_gamma to Maxima, for instance. This can be fixed by following the models in the functions/ directory; it should be possible to make it a GinacFunction. Before doing so, though, one will have to resolve #10072, since the evaluation will be wrong (?) otherwise.

Apply only attachment: trac_10075.patch.

Depends on #12507
Depends on #9130

CC: @sagetrac-ktkohl @benjaminfjones

Component: symbolics

Keywords: sd35.5

Author: Karen Kohl, Karl-Dieter Crisman

Reviewer: Karl-Dieter Crisman, Benjamin Jones

Merged: sage-5.0.beta7

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

The text was updated successfully, but these errors were encountered:

sagetrac-ktkohl · 2012-01-10T14:59:12Z

Changed keywords from none to sd35.5

sagetrac-ktkohl · 2012-01-11T04:27:14Z

comment:3

Sage gives this error message on startup:

ValueError: cannot find GiNaC function with name lgamma and 1 arguments

with this change in functions/other.py:

class Function_log_gamma(GinacFunction):

  def __init__(self):
    GinacFunction.__init__(self, "log_gamma", latex_name=r'\log\Gamma',
      ginac_name='lgamma', conversions={'mathematica':'LogGamma','maxima':'log_gamma'})

kcrisman · 2012-01-11T15:48:59Z

comment:4

More precisely,


class Function_log_gamma(GinacFunction):
    def __init__(self):
        GinacFunction.__init__(self, "log_gamma", latex_name=r'\log\Gamma',
          ginac_name='lgamma', conversions={'mathematica':'LogGamma','maxima':'log_gamma'})

log_gamma = Function_log_gamma()

causes this failure.

burcin · 2012-01-11T18:13:34Z

comment:5

It works if you drop the ginac_name argument. The function is named log_gamma in pynac.

kcrisman · 2012-01-11T19:01:20Z

comment:6

Weird. So what about things like

unsigned lgamma_serial "GiNaC::lgamma_SERIAL::serial" # logarithm of gamma function

?

sagetrac-ktkohl · 2012-01-11T20:52:14Z

Attachment: trac_10075_log_gamma.patch.gz

symbolic log_gamma (with modification of functions.rst in case merged before #9130)

sagetrac-ktkohl · 2012-01-11T20:52:36Z

symbolic log_gamma (with modification of functions.rst in case merged after #9130)

sagetrac-ktkohl · 2012-01-11T20:56:25Z

Author: Karen T. Kohl

sagetrac-ktkohl · 2012-01-11T20:56:25Z

comment:7

Attachment: trac_10075_log_gamma_without_functions.rst.patch.gz

Load one of the above two patches depending on whether the functions.rst documentation file has been modified already (as in the combined patch for #9130).

The second patch file above (without functions.rst) was edited by hand from the first.

kcrisman · 2012-01-12T16:30:28Z

comment:8

Burcin points out that

sage: log_gamma(-2.1) 
NaN

is not good. Sage itself does

sage: log(gamma(-2.1))
1.53171380819509 + 3.14159265358979*I

but Wolfram Alpha says

1.53171... - 9.42478... i

so the branches seem to differ even there.

kcrisman · 2012-02-16T15:54:21Z

comment:9

Since #9130 has positive review: Apply only attachment: trac_10075_log_gamma_without_functions.rst.patch.

One might think that Burcin's comment about -2.1 makes this 'needs work', but that is actually the current Sage behavior as well, so in principle that would be a different ticket, since making log_gamma symbolic would not introduce a regression...

In fact,

sage: log_gamma(-2.1)
NaN
sage: log_gamma(-3.1)
0.400311696703985
sage: log_gamma(-4.1)
NaN
sage: log_gamma(-5.1)
-2.63991581673655
sage: get_systems('log_gamma(2.1)')
['MPFR']

Apparently this is how MPFR deals with this function. So maybe all is well?

kcrisman · 2012-02-16T16:13:35Z

comment:10

Apparently this is how MPFR deals with this function. So maybe all is well?

I mean, for this ticket. Though we should not claim that it is evaluated by Ginac, because it isn't (all the above is in Sage with or without this patch).

Believe it or not:

Not any negative value, but in lngamma.c:

  /* if x < 0 and -2k-1 <= x <= -2k, then lngamma(x) = NaN */

probably because the gamma value is negative. This is because MPFR defines
lngamma as log(gamma(x)) while the C standard defines it as log|gamma(x)|. I
wonder if this should be regarded as a bug or if a new function (say,
mpfr_lgamma) should be defined in MPFR (in which case, not before 2.3.0). Do
other standards (other languages) define such a function, either as
log(gamma(x)) or as log|gamma(x)|?

I'm cc:ing Paul Z. just to confirm that this is intended MPFR behavior. We should then open another ticket to make sure to use mpmath or ginac or something to get complex answers. We currently somehow use PARI to get the complex versions.

sage: log_gamma(CC(-2.1))
1.53171380819509 + 3.14159265358979*I
sage: from sage.misc.citation import get_systems
sage: get_systems('log_gamma(CC(-2.1))')
['PARI', 'MPFR']

kcrisman · 2012-02-16T16:19:43Z

Reviewer: Karl-Dieter Crisman

kcrisman · 2012-02-16T16:19:43Z

comment:11

Ok, here we go.

sage: log_gamma(-21/10).n()
NaN
sage: get_systems('log_gamma(-21/10).n()')
['ginac']

So both give NaN, but we end up using RR.log_gamma() as in the GinacFunction code.

sage: log_gamma(-31/10).n()
0.400311696703985
sage: log_gamma(-3.1)
0.400311696703985
sage: a = RR(5)
sage: a.log_gamma()
3.17805383034795

I don't see anything holding this up except cosmetics. I'll try to make a refreshed patch momentarily.

kcrisman · 2012-02-16T17:10:51Z

comment:12

Okay, I've been messing with this for too long today.

sage: get_systems('log_gamma(SR(6))')
['ginac', 'GMP']
sage: get_systems('log_gamma(RR(6))')
[]
sage: get_systems('log_gamma(CC(6))')
['PARI', 'MPFR']
sage: get_systems('log_gamma(6.)')
['MPFR']

See also #10072, where a lot of the numerical evaluation was fixed. Anyway, updated patch with more explanation and other information coming up. It needs light review; no code was changed, only doctests.

I'm not sure I like the last doctest either

           sage: conjugate(log_gamma(-2))
            conjugate(+Infinity)

What is the conjugate of plus infinity? But I'll leave it for now, just to document it, unless someone has an objection, since we have in vanilla Sage

sage: conjugate(+Infinity)
conjugate(+Infinity)

I've opened #12521 for the evaluation at negative input with even ceiling function issue (i.e., log_gamma(-2.1)).

kcrisman · 2012-02-16T17:35:00Z

Changed author from Karen T. Kohl to Karen T. Kohl, Karl-Dieter Crisman

kcrisman · 2012-02-16T17:35:00Z

Dependencies: #12507, #9130

kcrisman · 2012-02-16T17:35:00Z

comment:13

I'm marking this as 'needs review', because I did change a fair number of tests. This definitely depends on #9130 because of some doc fixes. Also, I am marking this as depending on #12507, because I don't want to bother fixing that doctest if no one else is either. However, I don't really care either way.

benjaminfjones · 2012-02-17T02:06:56Z

comment:14

The latest patch attachment: trac_10075.patch failed to apply on top of 5.0.beta4 with this patch queue:

trac_12507_v2.patch
trac_9130-beta_function.2.patch
trac_9130-py_float_segfault.take2.patch
trac_9130-reviewer.patch

the failure seems to be in all.py

~/sage/latest/devel/sage> hg qpush -v
applying trac_10075.patch
patching file sage/functions/all.py
Hunk #1 FAILED at 15
1 out of 2 hunks FAILED -- saving rejects to file sage/functions/all.py.rej

Here is sage/functions/all.py.rej

--- all.py
+++ all.py
@@ -16,7 +16,7 @@
 
 
 from other import ( ceil, floor, gamma, psi, factorial,
-                    abs_symbolic, erf, sqrt,
+                    abs_symbolic, erf, sqrt, log_gamma,
                     gamma_inc, incomplete_gamma,
                     arg, real_part, real,
                     imag_part, imag, imaginary, conjugate)

benjaminfjones · 2012-02-17T02:06:56Z

Changed reviewer from Karl-Dieter Crisman to Karl-Dieter Crisman, Benjamin Jones

kcrisman · 2012-02-17T03:31:20Z

comment:15

Makes sense, since beta is in that list now. I was not careful enough about the dependencies, I guess. Coming up.

kcrisman · 2012-02-17T03:31:40Z

Attachment: trac_10075.patch.gz

Apply only this patch

kcrisman · 2012-02-17T03:33:59Z

comment:16

Okay, all should be well now? Sorry about that.

benjaminfjones · 2012-02-18T01:30:45Z

comment:17

OK! All looks good now.

The patch now applies cleanly to 5.0.beta4 on top of the patch queue in my last comment. I've tested everything in sage/functions and sage/symbolic and running all tests now (I don't expect any problems). The docs look good too. Positive review.

benjaminfjones · 2012-02-18T03:32:01Z

comment:18

Actually, one test did fail, but I don't think it's due to this patch (right?)

File "/home/jonesbe/sage/sage-5.0.beta4/devel/sage/sage/misc/trace.py", line 61:
    sage: print s.before[s.before.find('-'):]
Expected:
    ---...
    ipdb> c
    2 * 5
Got:
    ----------------------------------------------------------------------
    | Sage Version 5.0.beta4, Release Date: 2012-02-14                   |
    | Type notebook() for the GUI, and license() for information.        |
    ----------------------------------------------------------------------
    **********************************************************************
    *                                                                    *
    * Warning: this is a prerelease version, and it may be unstable.     *
    *                                                                    *
    **********************************************************************
    trace('print factor(10)'); print 3+97
    s
    c
    Loading Sage library. Current Mercurial branch is: 
**********************************************************************
1 items had failures:
   1 of  11 in __main__.example_1
***Test Failed*** 1 failures.
For whitespace errors, see the file /home/jonesbe/.sage//tmp/trace_30044.py
	 [2.2 s]

I haven't seen that failure before.

jdemeyer · 2012-03-04T21:19:46Z

Merged: sage-5.0.beta7

fchapoton · 2016-04-15T20:08:09Z

Changed author from Karen T. Kohl, Karl-Dieter Crisman to Karen Kohl, Karl-Dieter Crisman

kcrisman added this to the sage-4.8 milestone Oct 5, 2010

kcrisman added c: symbolics labels Oct 5, 2010

kcrisman assigned burcin Oct 5, 2010

burcin modified the milestones: sage-4.8, sage-5.0 Jan 11, 2012

sagetrac-ktkohl mannequin added the s: needs review label Jan 11, 2012

This comment has been minimized.

Sign in to view

benjaminfjones added s: positive review and removed s: needs review labels Feb 18, 2012

jdemeyer removed the s: positive review label Mar 4, 2012

jdemeyer closed this as completed Mar 4, 2012

kcrisman mentioned this issue Jun 24, 2016

evaluate log gamma for complex input #12521

Closed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Make log gamma symbolic #10075

Make log gamma symbolic #10075

kcrisman commented Oct 5, 2010

sagetrac-ktkohl mannequin commented Jan 10, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

kcrisman commented Jan 11, 2012

burcin commented Jan 11, 2012

kcrisman commented Jan 11, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

kcrisman commented Jan 12, 2012

kcrisman commented Feb 16, 2012

This comment has been minimized.

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

benjaminfjones commented Feb 17, 2012

benjaminfjones commented Feb 17, 2012

kcrisman commented Feb 17, 2012

kcrisman commented Feb 17, 2012

This comment has been minimized.

kcrisman commented Feb 17, 2012

benjaminfjones commented Feb 18, 2012

benjaminfjones commented Feb 18, 2012

jdemeyer commented Mar 4, 2012

fchapoton commented Apr 15, 2016

Make log gamma symbolic #10075

Make log gamma symbolic #10075

Comments

kcrisman commented Oct 5, 2010

sagetrac-ktkohl mannequin commented Jan 10, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

kcrisman commented Jan 11, 2012

burcin commented Jan 11, 2012

kcrisman commented Jan 11, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

sagetrac-ktkohl mannequin commented Jan 11, 2012

kcrisman commented Jan 12, 2012

kcrisman commented Feb 16, 2012

This comment has been minimized.

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

kcrisman commented Feb 16, 2012

benjaminfjones commented Feb 17, 2012

benjaminfjones commented Feb 17, 2012

kcrisman commented Feb 17, 2012

kcrisman commented Feb 17, 2012

This comment has been minimized.

kcrisman commented Feb 17, 2012

benjaminfjones commented Feb 18, 2012

benjaminfjones commented Feb 18, 2012

jdemeyer commented Mar 4, 2012

fchapoton commented Apr 15, 2016