math: Gamma(x) special cases are wrong

[robpike]

from -1 through -33 the answer is always +Inf; after that it oscillates between -Inf and
+Inf. i'm not sure it matters but it doesn't seem right.

this is the Go 1 behavior on Feb 10, 2012.

package main

import "fmt"
import "math"

func main() {
    for i := -1; i > -100; i-- { fmt.Println(i, ":",math.Gamma(float64(i))) }
}

-1 : +Inf
-2 : +Inf
-3 : +Inf
-4 : +Inf
-5 : +Inf
-6 : +Inf
-7 : +Inf
-8 : +Inf
-9 : +Inf
-10 : +Inf
-11 : +Inf
-12 : +Inf
-13 : +Inf
-14 : +Inf
-15 : +Inf
-16 : +Inf
-17 : +Inf
-18 : +Inf
-19 : +Inf
-20 : +Inf
-21 : +Inf
-22 : +Inf
-23 : +Inf
-24 : +Inf
-25 : +Inf
-26 : +Inf
-27 : +Inf
-28 : +Inf
-29 : +Inf
-30 : +Inf
-31 : +Inf
-32 : +Inf
-33 : +Inf
-34 : -Inf
-35 : +Inf
-36 : -Inf
-37 : +Inf
-38 : -Inf
-39 : +Inf
-40 : -Inf

[dsymonds]

Comment 1:

My opinion is that it should return NaN for all non-positive integers, which makes
mathematical sense (Γ is not analytic at those points) as well as matching the return
value of C's gamma function*.
* except that C's gamma(0) returns -Inf.

[rsc]

Comment 2:

From C99 rationale: 7.12.8.4 The tgamma functions.  New feature for C99.  In many other
standards, the meaning of gamma has changed over the  years.  Originally, it computed
the logarithm of the absolute value of the mathematical gamma function, with an external
int, signgam, being set to the sign of the gamma function.  Then gamma was replaced with
lgamma, and gamma was slated to be withdrawn.  About that time, NCEG changed gamma to
compute the mathematical gamma function, and that is what was adopted into C99 CD1; but
it appears that the old meaning of gamma has not yet been withdrawn, so there would have
been a conflict between C99 and current industry practice.  C99 therefore changed the
name in the FCD to tgamma, meaning “true gamma,” to avoid this conflict.
C99 originally had tgamma(0) = NaN, but by TC2 (probably also TC1) they had changed it
to tgamma(+0) = +Inf, tgamma(-0) = -Inf.  Other non-positive integers do not have
"sides", so they are NaN instead.
http://std.dkuug.dk/jtc1/sc22/wg14/www/docs/dr_244.htm looks to be the defect report
that led to the change.
I believe we should match (the new) C99 in the absence of expert arguments to the
contrary, meaning:
    Special cases are:
    Gamma(+Inf) = +Inf
    Gamma(+0) = +Inf
    Gamma(-0) = -Inf
    Gamma(x) = NaN for integer x < 0
    Gamma(-Inf) = NaN
    Gamma(NaN) = NaN
The attached program prints:
    netlib  gamma   lgamma  tgamma
inf inf inf inf inf
-inf    nan nan inf nan
0   nan inf inf inf
-0  nan -inf    inf -inf
-1  nan nan inf nan
-2  nan nan inf nan
-3  nan nan inf nan
-4  nan nan inf nan
-5  nan nan inf nan
-6  nan nan inf nan
-7  nan nan inf nan
-8  nan nan inf nan
-9  nan nan inf nan
The netlib column is the netlib source that the Go file claims to have started with,
although it in a different place on netlib now and I cut out all the #ifdefs and ran it
through cb.

Attachments:

1. gamma.c (3434 bytes)
2. Screen shot 2012-02-10 at 12.25.47 AM.png (100773 bytes)
3. Screen shot 2012-02-10 at 12.26.33 AM.png (104552 bytes)

[rsc]

Comment 3:

This issue was closed by revision 5496e94.

Status changed to Fixed.