On Thu, Jan 21, 2016 at 9:01 PM, William Stein <wst...@gmail.com> wrote: > Hi Buck, > > I asked Mike Rubinstein, who has implemented the incomplete Gamma > function in various situations (e.g., I think he did the first ever > open source implementation that was general enough for me to use for > my research), written papers mentioning issues, etc., what he thinks, > and here's his response: "Sounds like he wants to implement the > incomplete gamma function. He doesn't provide any detail about > implementation, just the naming conventions. The different ways of > normalizing, or whether one is looking at gamma(a,x) or what he calls > lower_gamma(a,x) makes little difference for the end result... one can > get one from the other. > > Depending on the relative size of the parameters, it makes more sense > to compute gamma(a,x) or lower_gamma(a,x). > > The fact that he is using the letter x suggests he is thinking of > implementing it for real values of x. Might be worthwhile to allow > both parameters to be complex. Certainly that would be better from the > point of view of anyone who uses functions of a complex variable. > > The incomplete gamma function changes its behaviour when a is close to > x. That's where things get a bit tricky (for truncation bounds), > especially in the complex plane. > > Section 3.4.2 of my paper: > http://arxiv.org/pdf/math/0412181v1.pdf (published in Recent > Perspectives in Random Matrix Theory and Number Theory, F. Mezzadri > and N.C. Snaith editors, 407-483, 200) > has a discussion about how I've computed the incomplete gamma > function. Some of my code for that is here: > > https://code.google.com/p/l-calc/source/browse/include/Lgamma.h > https://code.google.com/p/l-calc/source/browse/src/Lgamma.cc > (there's actual code in the .h file because of templates) > > Regarding his 'open questions', ex about notation as in `log_gamma() = > gamma(log=True)', I don't have any preferences. > > Best, > Mike"
Arb (which is now in Sage) permits computing incomplete gamma functions with rigorous error bounds over arbitrary-precision real/complex (interval) fields. Supporting the different regularized versions is just a matter of dividing by a complete gamma function. Looking at the attachment, I don't know if the logarithmic versions of incomplete gamma functions are really useful. The regularized incomplete gamma functions P and Q are quite standard and used in statistics, I think. Personally, I think the usual gamma() function is too important in its own right to complicate by adding extra options. I would rather leave gamma() a single-argument function and have separate functions for logarithmic and incomplete gamma functions of various kinds. In mpmath, there is a single function gammainc() that covers all cases of incomplete gamma functions. But it's not necessary to do it in that particular way. If you want symbolic manipulation to be aware of different variants, you probably need to do something in Pynac/SymPy/Maxima rather than Sage itself. Fredrik -- You received this message because you are subscribed to the Google Groups "sage-devel" group. To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+unsubscr...@googlegroups.com. To post to this group, send email to sage-devel@googlegroups.com. Visit this group at https://groups.google.com/group/sage-devel. For more options, visit https://groups.google.com/d/optout.