[PATCH] Fix for logb/logbf/logbl (bz 13954/13955/13956)
Andreas Jaeger
aj@suse.com
Tue May 1 07:48:00 GMT 2012
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
More information about the Libc-alpha
mailing list