[SciPy-Dev] Use NPY_INFINITY in cephes code in scipy.special?

[Warren Weckesser]

scipy.special.zeta(1.0, 1) returns 1.7976931348623157e+308
to indicate infinity. Here's the relevant snippet from
zeta.c:

-----
if( x == 1.0 )
    goto retinf;

if( x < 1.0 )
    {
domerr:
    mtherr( "zeta", DOMAIN );
    return(NPY_NAN);
    }

if( q <= 0.0 )
    {
    if(q == floor(q))
        {
        mtherr( "zeta", SING );
retinf:
        return( MAXNUM );
        }
    if( x != floor(x) )
        goto domerr; /* because q^-x not defined */
    }
-----

Is there any problem with changing MAXNUM to NPY_INFINITY?
Note that the use of NPY_NAN was added in r5880, with the
check-in comment "Use npy_math NAN and INFINITY instead of
hacked/buggy version of cephes.".  David, if you're out
there, can you see any reason to not use NPY_INFINITY here
instead of MAXNUM?

Similar (ab)uses of MAXNUM occur in other special
functions.  For example, in ndtri() (the inverse of the
normal distribution's CDF), we have:

-----
if( y0 <= 0.0 )
    {
    mtherr( "ndtri", DOMAIN );
    return( -MAXNUM );
    }
if( y0 >= 1.0 )
    {
    mtherr( "ndtri", DOMAIN );
    return( MAXNUM );
    }
-----

but this makes more sense to me (and it matches the
behavior of scipy.stats.norm.ppf()):

-----
if (y0 < 0.0 || y > 1.0)
    {
    mtherr("ndtri", DOMAIN);
    return(NPY_NAN);
    }
if (y0 == 0.0)
    {
    return(-NPY_INFINITY);
    }
if (y0 == 1.0)
    {
    return(NPY_INFINITY);
    }
-----

Before I pursue this further, is there any reason
to *not* make changes like these?

Warren
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20101204/8602029a/attachment.html>