[BZ #16823] Fix log1pl returning wrong infinity sign

diff --git a/ChangeLog b/ChangeLog
index 4007295..803d114 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,11 @@
+2014-04-29  Stefan Liebler  <stli@linux.vnet.ibm.com>
+
+       [BZ #16823]
+       * sysdeps/ieee754/ldbl-128/s_log1pl.c (__log1pl):
+       Always divide by positive zero when computing -Inf result.
+       * sysdeps/ieee754/dbl-64/s_log1p.c (__log1p): Likewise.
+       * sysdeps/ieee754/flt-32/s_log1pf.c (__log1pf): Likewise.
+
 2014-04-28  Adhemerval Zanella  <azanella@linux.vnet.ibm.com>
 
        * sysdeps/powerpc/fpu/fclrexcpt.c (__feclearexcept): Do not update

diff --git a/sysdeps/ieee754/dbl-64/s_log1p.c b/sysdeps/ieee754/dbl-64/s_log1p.c
index fd4dce5..c922148 100644 (file)
--- a/sysdeps/ieee754/dbl-64/s_log1p.c
+++ b/sysdeps/ieee754/dbl-64/s_log1p.c
@@ -110,7 +110,7 @@ __log1p (double x)
       if (__glibc_unlikely (ax >= 0x3ff00000))           /* x <= -1.0 */
        {
          if (x == -1.0)
-           return -two54 / (x - x);            /* log1p(-1)=+inf */
+           return -two54 / zero;               /* log1p(-1)=-inf */
          else
            return (x - x) / (x - x);           /* log1p(x<-1)=NaN */
        }

diff --git a/sysdeps/ieee754/flt-32/s_log1pf.c b/sysdeps/ieee754/flt-32/s_log1pf.c
index 0307277..5f00feb 100644 (file)
--- a/sysdeps/ieee754/flt-32/s_log1pf.c
+++ b/sysdeps/ieee754/flt-32/s_log1pf.c
@@ -42,7 +42,7 @@ __log1pf(float x)
        k = 1;
        if (hx < 0x3ed413d7) {                  /* x < 0.41422  */
            if(ax>=0x3f800000) {                /* x <= -1.0 */
-               if(x==(float)-1.0) return -two25/(x-x); /* log1p(-1)=+inf */
+               if(x==(float)-1.0) return -two25/zero; /* log1p(-1)=-inf */
                else return (x-x)/(x-x);        /* log1p(x<-1)=NaN */
            }
            if(ax<0x31000000) {                 /* |x| < 2**-29 */

diff --git a/sysdeps/ieee754/ldbl-128/s_log1pl.c b/sysdeps/ieee754/ldbl-128/s_log1pl.c
index d991e8a..d8d89f0 100644 (file)
--- a/sysdeps/ieee754/ldbl-128/s_log1pl.c
+++ b/sysdeps/ieee754/ldbl-128/s_log1pl.c
@@ -150,7 +150,7 @@ __log1pl (long double xm1)
   if (x <= 0.0L)
     {
       if (x == 0.0L)
-       return (-1.0L / (x - x));
+       return (-1.0L / zero);  /* log1p(-1) = -inf */
       else
        return (zero / (x - x));
     }