This is the mail archive of the libc-alpha@sourceware.org mailing list for the glibc project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

Re: [PATCH] Fix for logb/logbf/logbl (bz 13954/13955/13956)


On 04/30/2012 09:43 PM, Adhemerval Zanella wrote:
POSIX 2008 states that if the input for 'logb[f|l]' is a subnormal number, it should
be treated as it were normalized. It means the implementation should calculate the
log2 of mantissa and add it subnormal exponent number (-126 for float and -1022 for
double and IBM long double). This patch takes care of that.

On which platforms was this tested?


---

2012-04-30 Adhemerval Zanella<azanella@linux.vnet.ibm.com>

[BZ #13954, BZ#13955, and BZ #13956]

We normally write: [BZ #13954] [BZ #13955] [BZ #13956]

	* sysdeps/ieee754/dbl-64/s_logb.c (__logb): Fix for subnormal number.
	* sysdeps/ieee754/flt-32/s_logbf.c (__logf): Likewise.
	* sysdeps/ieee754/ldbl-128ibm/s_logbl.c (__logbl): Likewise.
	* math/libm-test.inc (logb_test) : Additional logb tests.

What about the dbl-64/wordsize-64, ldbl-96 and ldbl-128 routines? Do they have the same problem and need fixing as well?


Your code change looks fine. Could you check the other implementations as well and resubmit then?

Thanks,
Andreas



diff --git a/math/libm-test.inc b/math/libm-test.inc
index 1b5d1c7..5ef088b 100644
--- a/math/libm-test.inc
+++ b/math/libm-test.inc
@@ -5098,6 +5098,18 @@ logb_test (void)
    TEST_f_f (logb, 1024, 10);
    TEST_f_f (logb, -2000, 10);

+  TEST_f_f (logb, 0x0.1p-127, -131);
+  TEST_f_f (logb, 0x0.01p-127, -135);
+  TEST_f_f (logb, 0x0.011p-127, -135);
+#ifndef TEST_FLOAT
+  TEST_f_f (logb, 0x0.8p-1022, -1023);
+  TEST_f_f (logb, 0x0.1p-1022, -1026);
+  TEST_f_f (logb, 0x0.00111p-1022, -1034);
+  TEST_f_f (logb, 0x0.00001p-1022, -1042);
+  TEST_f_f (logb, 0x0.000011p-1022, -1042);
+  TEST_f_f (logb, 0x0.0000000000001p-1022, -1074);
+#endif
+
    END (logb);
  }

diff --git a/sysdeps/ieee754/dbl-64/s_logb.c b/sysdeps/ieee754/dbl-64/s_logb.c
index 2382fbb..baa35e1 100644
--- a/sysdeps/ieee754/dbl-64/s_logb.c
+++ b/sysdeps/ieee754/dbl-64/s_logb.c
@@ -10,10 +10,6 @@
   * ====================================================
   */

-#if defined(LIBM_SCCS)&&  !defined(lint)
-static char rcsid[] = "$NetBSD: s_logb.c,v 1.8 1995/05/10 20:47:50 jtc Exp $";
-#endif
-
  /*
   * double logb(x)
   * IEEE 754 logb. Included to pass IEEE test suite. Not recommend.
@@ -23,20 +19,29 @@ static char rcsid[] = "$NetBSD: s_logb.c,v 1.8 1995/05/10 20:47:50 jtc Exp $";
  #include<math.h>
  #include<math_private.h>

-double __logb(double x)
+double
+__logb (double x)
  {
-	int32_t lx,ix;
-	EXTRACT_WORDS(ix,lx,x);
-	ix&= 0x7fffffff;			/* high |x| */
-	if((ix|lx)==0) return -1.0/fabs(x);
-	if(ix>=0x7ff00000) return x*x;
-	if((ix>>=20)==0) 			/* IEEE 754 logb */
-		return -1022.0;
-	else
-		return (double) (ix-1023);
+  int32_t lx, ix, rix;
+
+  EXTRACT_WORDS (ix, lx, x);
+  ix&= 0x7fffffff;		/* high |x| */
+  if ((ix | lx) == 0)
+    return -1.0 / fabs (x);
+  if (ix>= 0x7ff00000)
+    return x * x;
+  if (__builtin_expect ((rix = ix>>  20) == 0, 0))
+    {
+      /* POSIX specifies that denormal number is treated as
+         though it were normalized.  */
+      int m1 = (ix == 0) ? 0 : __builtin_clz (ix);
+      int m2 = (lx == 0) ? 0 : __builtin_clz (lx);
+      int ma = (m1 == 0) ? m2 + 32 : m1;
+      return -1022.0 + (double)(11 - ma);
+    }
+  return (double) (rix - 1023);
  }
  weak_alias (__logb, logb)
  #ifdef NO_LONG_DOUBLE
-strong_alias (__logb, __logbl)
-weak_alias (__logb, logbl)
+strong_alias (__logb, __logbl) weak_alias (__logb, logbl)
  #endif
diff --git a/sysdeps/ieee754/flt-32/s_logbf.c b/sysdeps/ieee754/flt-32/s_logbf.c
index b6aa0f0..166d1a4 100644
--- a/sysdeps/ieee754/flt-32/s_logbf.c
+++ b/sysdeps/ieee754/flt-32/s_logbf.c
@@ -13,23 +13,27 @@
   * ====================================================
   */

-#if defined(LIBM_SCCS)&&  !defined(lint)
-static char rcsid[] = "$NetBSD: s_logbf.c,v 1.4 1995/05/10 20:47:51 jtc Exp $";
-#endif
-
  #include<math.h>
  #include<math_private.h>

-float __logbf(float x)
+float
+__logbf (float x)
  {
-	int32_t ix;
-	GET_FLOAT_WORD(ix,x);
-	ix&= 0x7fffffff;			/* high |x| */
-	if(ix==0) return (float)-1.0/fabsf(x);
-	if(ix>=0x7f800000) return x*x;
-	if((ix>>=23)==0) 			/* IEEE 754 logb */
-		return -126.0;
-	else
-		return (float) (ix-127);
+  int32_t ix, rix;
+
+  GET_FLOAT_WORD (ix, x);
+  ix&= 0x7fffffff;		/* high |x| */
+  if (ix == 0)
+    return (float) -1.0 / fabsf (x);
+  if (ix>= 0x7f800000)
+    return x * x;
+  if (__builtin_expect ((rix = ix>>  23) == 0, 0))
+    {
+      /* POSIX specifies that denormal number is treated as
+         though it were normalized.  */
+      int m = (ix == 0) ? 0 : __builtin_clz (ix);
+      return -126.0 + (double)(8 -m);
+    }
+  return (float) (rix - 127);
  }
  weak_alias (__logbf, logbf)
diff --git a/sysdeps/ieee754/ldbl-128ibm/s_logbl.c b/sysdeps/ieee754/ldbl-128ibm/s_logbl.c
index f38b129..accd2c2 100644
--- a/sysdeps/ieee754/ldbl-128ibm/s_logbl.c
+++ b/sysdeps/ieee754/ldbl-128ibm/s_logbl.c
@@ -13,10 +13,6 @@
   * ====================================================
   */

-#if defined(LIBM_SCCS)&&  !defined(lint)
-static char rcsid[] = "$NetBSD: $";
-#endif
-
  /*
   * long double logbl(x)
   * IEEE 754 logb. Included to pass IEEE test suite. Not recommend.
@@ -27,16 +23,27 @@ static char rcsid[] = "$NetBSD: $";
  #include<math_private.h>
  #include<math_ldbl_opt.h>

-long double __logbl(long double x)
+long double
+__logbl (long double x)
  {
-	int64_t lx,hx;
-	GET_LDOUBLE_WORDS64(hx,lx,x);
-	hx&= 0x7fffffffffffffffLL;		/* high |x| */
-	if((hx|(lx&0x7fffffffffffffffLL))==0) return -1.0/fabs(x);
-	if(hx>=0x7ff0000000000000LL) return x*x;
-	if((hx>>=52)==0) 			/* IEEE 754 logb */
-		return -1022.0;
-	else
-		return (long double) (hx-0x3ff);
+  int64_t lx, hx, rhx;
+
+  GET_LDOUBLE_WORDS64 (hx, lx, x);
+  hx&= 0x7fffffffffffffffLL;	/* high |x| */
+  if ((hx | (lx&  0x7fffffffffffffffLL)) == 0)
+    return -1.0 / fabs (x);
+  if (hx>= 0x7ff0000000000000LL)
+    return x * x;
+  if (__builtin_expect ((rhx = hx>>  52) == 0, 0))
+    {
+      /* POSIX specifies that denormal number is treated as
+         though it were normalized.  */
+      int m1 = (hx == 0) ? 0 : __builtin_clzl (hx);
+      int m2 = (lx == 0) ? 0 : __builtin_clzl (lx);
+      int ma = (m1 == 0) ? m2 + 64 : m1;
+      return -1022.0 + (double)(11 - ma);
+    }
+  return (long double) (rhx - 1023);
  }
+
  long_double_symbol (libm, __logbl, logbl);


--
 Andreas Jaeger aj@{suse.com,opensuse.org} Twitter/Identica: jaegerandi
  SUSE LINUX Products GmbH, Maxfeldstr. 5, 90409 Nürnberg, Germany
   GF: Jeff Hawn,Jennifer Guild,Felix Imendörffer,HRB16746 (AG Nürnberg)
    GPG fingerprint = 93A3 365E CE47 B889 DF7F  FED1 389A 563C C272 A126


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]