BUG: gammainc C code does not have a special case for np.inf

[amyoshino]

Describe the issue:

gammainc C code does not have a special case for np.inf.

Instead of returning 0, 1 or nan it always return nan.

The solution is to add

Reproducable code example:

import numpy as np
import pytensor
import scipy.special as sp

from pytensor import tensor as pt
from pytensor.scalar.math import gammaincinv, gammainc

x1 = pt.dscalar()
x2 = pt.dscalar()
x2 = pt.dscalar()
y = gammainc(x1, x2)
test_func = pytensor.function([x1, x2], y)
# Compare output from pytensor with scipy - this should return 0.0 but instead returns nan
a, b = np.inf, 3.8
print(f"{a=}, {b=}, {test_func(a, b)}, {sp.gammainc(a, b)}")

# Compare output from pytensor with scipy - this should return 1.0 but instead returns nan
a, b = 2.4, np.inf
print(f"{a=}, {b=}, {test_func(a, b)}, {sp.gammainc(a, b)}")

Error message:

No error messages. The code fails silently returning the wrong value.

PyTensor version information:

pytensor: rel-2.18.6

Context for the issue:

This issue appeared while implementing pymc-devs/pymc#6845

Will submit a PR to address this bug.