diff options
| -rw-r--r-- | ChangeLog | 13 | ||||
| -rw-r--r-- | NEWS | 5 | ||||
| -rw-r--r-- | math/libm-test.inc | 16 | ||||
| -rw-r--r-- | sysdeps/ieee754/dbl-64/s_logb.c | 37 | ||||
| -rw-r--r-- | sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c | 13 | ||||
| -rw-r--r-- | sysdeps/ieee754/flt-32/s_logbf.c | 32 | ||||
| -rw-r--r-- | sysdeps/ieee754/ldbl-128/s_logbl.c | 31 | ||||
| -rw-r--r-- | sysdeps/ieee754/ldbl-128ibm/s_logbl.c | 35 | ||||
| -rw-r--r-- | sysdeps/ieee754/ldbl-96/s_logbl.c | 35 |
9 files changed, 143 insertions, 74 deletions
diff --git a/ChangeLog b/ChangeLog index 8bceaab8be..6848fb4838 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,16 @@ +2012-05-09 Adhemerval Zanella <azanella@linux.vnet.ibm.com> + + [BZ #13954] + [BZ #13955] + [BZ #13956] + * sysdeps/ieee754/dbl-64/s_logb.c (__logb): Fix for subnormal number. + * sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c (__logb): Likewise. + * sysdeps/ieee754/flt-32/s_logbf.c (__logf): Likewise. + * sysdeps/ieee754/ldbl-128/s_logbl.c (__logbl): Likewise. + * sysdeps/ieee754/ldbl-128ibm/s_logbl.c (__logbl): Likewise. + * sysdeps/ieee754/ldbl-96/s_logbl.c (__logbl): Likewise. + * math/libm-test.inc (logb_test) : Additional logb tests. + 2012-05-09 Andreas Schwab <schwab@linux-m68k.org> Andreas Jaeger <aj@suse.de> diff --git a/NEWS b/NEWS index 5c90f3da28..b3a1233482 100644 --- a/NEWS +++ b/NEWS @@ -23,8 +23,9 @@ Version 2.16 13852, 13854, 13871, 13872, 13873, 13879, 13883, 13884, 13885, 13886, 13892, 13895, 13908, 13910, 13911, 13912, 13913, 13914, 13915, 13916, 13917, 13918, 13919, 13920, 13921, 13922, 13923, 13924, 13926, 13927, - 13928, 13938, 13941, 13942, 13963, 13967, 13970, 13973, 13979, 13983, - 14027, 14033, 14034, 14040, 14049, 14053, 14055, 14064, 14080, 14083 + 13928, 13938, 13941, 13942, 13954, 13955, 13956, 13963, 13967, 13970, + 13973, 13979, 13983, 14027, 14033, 14034, 14040, 14049, 14053, 14055, + 14064, 14080, 14083 * ISO C11 support: diff --git a/math/libm-test.inc b/math/libm-test.inc index 542131dc7b..5a38dbf3a7 100644 --- a/math/libm-test.inc +++ b/math/libm-test.inc @@ -5376,6 +5376,22 @@ 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 +#if defined TEST_LDOUBLE && LDBL_MIN_EXP - LDBL_MANT_DIG <= -16400 + TEST_f_f (logb, 0x1p-16400L, -16400); + TEST_f_f (logb, 0x.00000000001p-16382L, -16426); +#endif + END (logb); } diff --git a/sysdeps/ieee754/dbl-64/s_logb.c b/sysdeps/ieee754/dbl-64/s_logb.c index 2382fbb414..baa35e14d8 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/dbl-64/wordsize-64/s_logb.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c index 2ad6c7ddbd..474eeef36b 100644 --- a/sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c +++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_logb.c @@ -1,5 +1,5 @@ /* Compute radix independent exponent. - Copyright (C) 2011 Free Software Foundation, Inc. + Copyright (C) 2011, 2012 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@gmail.com>, 2011. @@ -25,16 +25,21 @@ double __logb (double x) { - int64_t ix; + int64_t ix, ex; EXTRACT_WORDS64 (ix, x); ix &= UINT64_C(0x7fffffffffffffff); if (ix == 0) return -1.0 / fabs (x); - unsigned int ex = ix >> 52; + ex = ix >> 52; if (ex == 0x7ff) return x * x; - return ex == 0 ? -1022.0 : (double) (ex - 1023); + if (__builtin_expect (ex == 0, 0)) + { + int m = (ix == 0) ? 0 : __builtin_clzl (ix); + return -1022.0 + (double)(11 -m); + } + return (double) (ex - 1023); } weak_alias (__logb, logb) #ifdef NO_LONG_DOUBLE diff --git a/sysdeps/ieee754/flt-32/s_logbf.c b/sysdeps/ieee754/flt-32/s_logbf.c index b6aa0f057d..025c70de7e 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 + (float)(8 - m); + } + return (float) (rix - 127); } weak_alias (__logbf, logbf) diff --git a/sysdeps/ieee754/ldbl-128/s_logbl.c b/sysdeps/ieee754/ldbl-128/s_logbl.c index 0b09b289c2..cf6003e055 100644 --- a/sysdeps/ieee754/ldbl-128/s_logbl.c +++ b/sysdeps/ieee754/ldbl-128/s_logbl.c @@ -26,16 +26,27 @@ static char rcsid[] = "$NetBSD: $"; #include <math.h> #include <math_private.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)==0) return -1.0/fabs(x); - if(hx>=0x7fff000000000000LL) return x*x; - if((hx>>=48)==0) /* IEEE 754 logb */ - return -16382.0; - else - return (long double) (hx-0x3fff); + int64_t lx, hx, ex; + + GET_LDOUBLE_WORDS64 (hx, lx, x); + hx &= 0x7fffffffffffffffLL; /* high |x| */ + if ((hx | lx) == 0) + return -1.0 / fabs (x); + if (hx >= 0x7fff000000000000LL) + return x * x; + if ((ex = hx >> 48) == 0) /* IEEE 754 logb */ + { + /* POSIX specifies that denormal number is treated as + though it were normalized. */ + int m1 = (hx == 0) ? 0 : __builtin_clzll (hx); + int m2 = (lx == 0) ? 0 : __builtin_clzll (lx); + int ma = (m1 == 0) ? m2 + 64 : m1; + return -16382.0 + (long double)(15 - ma); + } + return (long double) (ex - 16383); } + weak_alias (__logbl, logbl) diff --git a/sysdeps/ieee754/ldbl-128ibm/s_logbl.c b/sysdeps/ieee754/ldbl-128ibm/s_logbl.c index f38b129971..678b6cad57 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_clzll (hx); + int m2 = (lx == 0) ? 0 : __builtin_clzll (lx); + int ma = (m1 == 0) ? m2 + 64 : m1; + return -1022.0 + (long double)(11 - ma); + } + return (long double) (rhx - 1023); } + long_double_symbol (libm, __logbl, logbl); diff --git a/sysdeps/ieee754/ldbl-96/s_logbl.c b/sysdeps/ieee754/ldbl-96/s_logbl.c index 95b644c030..d8ad4bcfcf 100644 --- a/sysdeps/ieee754/ldbl-96/s_logbl.c +++ b/sysdeps/ieee754/ldbl-96/s_logbl.c @@ -14,10 +14,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.h> #include <math_private.h> -long double __logbl(long double x) +long double +__logbl (long double x) { - int32_t es,lx,ix; - GET_LDOUBLE_WORDS(es,ix,lx,x); - es &= 0x7fff; /* exponent */ - if((es|ix|lx)==0) return -1.0/fabs(x); - if(es==0x7fff) return x*x; - if(es==0) /* IEEE 754 logb */ - return -16382.0; - else - return (long double) (es-0x3fff); + int32_t es, lx, ix; + + GET_LDOUBLE_WORDS (es, ix, lx, x); + es &= 0x7fff; /* exponent */ + if ((es | ix | lx) == 0) + return -1.0 / fabs (x); + if (es == 0x7fff) + return x * x; + if (es == 0) /* IEEE 754 logb */ + { + /* 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 -16382.0 - (long double)(ma); + } + return (long double) (es - 16383); } + weak_alias (__logbl, logbl) |