math: Gamma(x) for -ve half integral values < -170.5 returns +Inf not 0

[ncw]

Se here for a demo: http://play.golang.org/p/4npToJjm79

    for i := -100; i >= -200; i-- {
        x := float64(i) - 0.5
        fmt.Printf("gamma(%g) = %g\n", x, math.Gamma(x))
    }

Returns (trimmed)

gamma(-166.5) = -2.7020668623696067e-299
gamma(-167.5) = 1.613174246190816e-301
gamma(-168.5) = -9.573734398758526e-304
gamma(-169.5) = 5.6482208842233304e-306
gamma(-170.5) = -3.3127395215386025e-308
gamma(-171.5) = +Inf
gamma(-172.5) = +Inf
gamma(-173.5) = +Inf
gamma(-174.5) = +Inf
gamma(-175.5) = +Inf

[andrewaustin]

This ended up being a small fix. See: https://go-review.googlesource.com/#/c/17236/

The Gamma function was checking for values < -170.5 and always returning +Inf instead of zero. Removing this check fixed the issue.

For sanity, I also compared the results to tgamma from cmath with the following code:

#include <iostream>
#include <cmath>

int main(int argc, char **argv)
{
    double a;
    for(a = 299.5; a > -300.5; a = a - 1.0)
    {
        std::cout << "gamma(" << a << ") = " << tgamma(a) << std::endl;
    }
    return 0;
}

You can see the results from my change here: https://docs.google.com/spreadsheets/d/14gBcKZxs0Oh4km40BatucqDrxqEgUREpFajCUJoLf1E/edit#gid=0&vpid=A1

We now match cmath's gamma for the 600 values tested.

[gopherbot]

CL https://golang.org/cl/17236 mentions this issue.

[andrewaustin]

On the CL it has been pointed out that there are some values that are still incorrect, namely we know of at least two values that are tested for in CPython's math test cases (https://github.com/python/cpython/blob/1fe0fd9feb6a4472a9a1b186502eb9c0b2366326/Lib/test/math_testcases.txt#L275) that fail:

gamma(-1000.5) should be -0.0, but we give NaN.
gamma(-1000000000.3) should be -0.0, but we give NaN

I looked deeper at the Cephes C library that the Go code was based and it too incorrectly returns NaN for these values:

gamma(-1000.5) = nan
gamma(-1000000000.3) = nan

The current Go 1.5.1 implementation always returns +Inf for values < -170.5.

[bradfitz]

Sorry, we lost this during the Go 1.7 cycle.

[gopherbot]

CL https://golang.org/cl/30146 mentions this issue.

[bradfitz]

This fix broke the 387 builder:

--- FAIL: TestGamma (0.00s)
    all_test.go:2228: Gamma(-171.5) = 1.9316265431712e-310, want 0
FAIL
FAIL    math    0.014s

[ncw]

That is actually the correct value according to python!

>>> math.gamma(-171.5)
1.9316265431712e-310

[rsc]

@ncw thanks for that clue.