diff options
Diffstat (limited to 'sysdeps/ieee754/ldbl-128')
| -rw-r--r-- | sysdeps/ieee754/ldbl-128/s_logbl.c | 31 |
1 files changed, 21 insertions, 10 deletions
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) |