lgamma_r and lgammaf_r return the wrong sign for -0.f

enh enh at google.com
Fri Sep 12 22:17:24 UTC 2014


On Fri, Sep 12, 2014 at 3:05 PM, Steve Kargl
<sgk at troutmask.apl.washington.edu> wrote:
> On Fri, Sep 12, 2014 at 02:15:37PM -0700, enh via freebsd-numerics wrote:
>> if I pass -0.f to lgammaf_r, the sign returned in *signgamp is 1. this
>> is incorrect --- it should be -1.
>>
>> both lgamma_r and lgammaf_r are affected, but the other special cases
>> in those functions look fine to me.
>>
>> this is fixed in OpenBSD and glibc, but FreeBSD and NetBSD both have
>> the same bug.
>
> Are you sure FreeBSD has a bug?  From n1256.pdf (committee draft of C99)

note that i'm reporting bugs in lgamma_r and lgammaf_r. lgamma and
lgammaf don't have an out argument for the sign.

>   F.9.5.3 The lgamma functions
>   -- lgamma(1) returns +0.
>   -- lgamma(2) returns +0.
>   -- lgamma(x) returns +inf and raises the ``divide-by-zero'' floating-pointr
>      exception for x a negative integer or zero.
>   -- lgamma(-inf) returns +inf.
>   -- lgamma(+inf) returns +inf.
>
> See the 3rd bullet.  -0 is 0 and -0 is a negative integer.
>
> POSIX appears to defer to ISO C.  n1570.pdf (committe draft for C11)
> has (almost?) identical text.
>
>> patch below (whitespace mangled courtesy of gmail). i'd prefer to wait
>> for this to be fixed in FreeBSD and pull down the fix rather than just
>> fix it locally.
>
> I have a bigger patch coming with ld80 and ld128 version of lgammal
> and lgammal_r.

sorry, i haven't even looked at the *l variants.

>> btw, it looks like you're missing coshl/sinhl/tanhl for ld128 now?
>> (they've been removed from imprecise.c without having ld128
>> implementations added afaics.)
>
> The ld80 and ld128 versions are in the same file.
>
> troutmask:fvwm:kargl[226] grep LDBL /usr/src/lib/msun/src/e_coshl.c
> #if LDBL_MAX_EXP != 0x4000
> #define BIAS    (LDBL_MAX_EXP - 1)
> #if LDBL_MANT_DIG == 64
> #elif LDBL_MANT_DIG == 113
> #endif /* LDBL_MANT_DIG == 64 */
>         GET_LDBL_EXPSIGN(ix,x);
>             if (ix<BIAS-(LDBL_MANT_DIG+1)/2)    /* |x| < TINY */
> #if LDBL_MANT_DIG == 64
> #elif LDBL_MANT_DIG == 113
>
> see also e_sinhl.c and s_tanhl.c.

d'oh. yes, i need to import those files!

 --elliott


More information about the freebsd-numerics mailing list