projects / thirdparty / glibc.git / blame
| Commit | Line | Data |
|---|---|---|
| 0ecb606c JJ |
1 | /* s_scalbnl.c -- long double version of s_scalbn.c. |
| 2 | * Conversion to IEEE quad long double by Jakub Jelinek, jj@ultra.linux.cz. | |
| 3 | */ | |
| 4 | ||
| 5 | /* @(#)s_scalbn.c 5.1 93/09/24 */ | |
| 6 | /* | |
| 7 | * ==================================================== | |
| 8 | * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. | |
| 9 | * | |
| 10 | * Developed at SunPro, a Sun Microsystems, Inc. business. | |
| 11 | * Permission to use, copy, modify, and distribute this | |
| 12 | * software is freely granted, provided that this notice | |
| 13 | * is preserved. | |
| 14 | * ==================================================== | |
| 15 | */ | |
| 16 | ||
| 17 | #if defined(LIBM_SCCS) && !defined(lint) | |
| 18 | static char rcsid[] = "$NetBSD: $"; | |
| 19 | #endif | |
| 20 | ||
| 21 | /* | |
| 22 | * scalbnl (long double x, int n) | |
| 23 | * scalbnl(x,n) returns x* 2**n computed by exponent | |
| 24 | * manipulation rather than by actually performing an | |
| 25 | * exponentiation or a multiplication. | |
| 26 | */ | |
| 27 | ||
| 28 | #include "math.h" | |
| 29 | #include "math_private.h" | |
| 30 | #include <math_ldbl_opt.h> | |
| 31 | ||
| 32 | #ifdef __STDC__ | |
| 33 | static const long double | |
| 34 | #else | |
| 35 | static long double | |
| 36 | #endif | |
| 37 | twolm54 = 5.55111512312578270212e-17, /* 0x3C90000000000000, 0 */ | |
| 38 | huge = 1.0E+300L, | |
| 39 | tiny = 1.0E-300L; | |
| 40 | #ifdef __STDC__ | |
| 41 | static const double | |
| 42 | #else | |
| 43 | static double | |
| 44 | #endif | |
| 45 | two54 = 1.80143985094819840000e+16, /* 0x4350000000000000 */ | |
| 46 | twom54 = 5.55111512312578270212e-17; /* 0x3C90000000000000 */ | |
| 47 | ||
| 48 | #ifdef __STDC__ | |
| 49 | long double __scalbnl (long double x, int n) | |
| 50 | #else | |
| 51 | long double __scalbnl (x,n) | |
| 52 | long double x; int n; | |
| 53 | #endif | |
| 54 | { | |
| 55 | int64_t k,l,hx,lx; | |
| 56 | union { int64_t i; double d; } u; | |
| 57 | GET_LDOUBLE_WORDS64(hx,lx,x); | |
| 58 | k = (hx>>52)&0x7ff; /* extract exponent */ | |
| 59 | l = (lx>>52)&0x7ff; | |
| 60 | if (k==0) { /* 0 or subnormal x */ | |
| 61 | if (((hx|lx)&0x7fffffffffffffffULL)==0) return x; /* +-0 */ | |
| 62 | u.i = hx; | |
| 63 | u.d *= two54; | |
| 64 | hx = u.i; | |
| 65 | k = ((hx>>52)&0x7ff) - 54; | |
| 66 | } | |
| 67 | else if (k==0x7ff) return x+x; /* NaN or Inf */ | |
| 68 | k = k+n; | |
| 69 | if (n> 50000 || k > 0x7fe) | |
| 70 | return huge*__copysignl(huge,x); /* overflow */ | |
| 71 | if (n< -50000) return tiny*__copysignl(tiny,x); /*underflow */ | |
| 72 | if (k > 0) { /* normal result */ | |
| 73 | hx = (hx&0x800fffffffffffffULL)|(k<<52); | |
| 74 | if ((lx & 0x7fffffffffffffffULL) == 0) { /* low part +-0 */ | |
| 75 | SET_LDOUBLE_WORDS64(x,hx,lx); | |
| 76 | return x; | |
| 77 | } | |
| 78 | if (l == 0) { /* low part subnormal */ | |
| 79 | u.i = lx; | |
| 80 | u.d *= two54; | |
| 81 | lx = u.i; | |
| 82 | l = ((lx>>52)&0x7ff) - 54; | |
| 83 | } | |
| 84 | l = l + n; | |
| 85 | if (l > 0) | |
| 86 | lx = (lx&0x800fffffffffffffULL)|(l<<52); | |
| 87 | else if (l <= -54) | |
| 88 | lx = (lx&0x8000000000000000ULL); | |
| 89 | else { | |
| 90 | l += 54; | |
| 91 | u.i = (lx&0x800fffffffffffffULL)|(l<<52); | |
| 92 | u.d *= twom54; | |
| 93 | lx = u.i; | |
| 94 | } | |
| 95 | SET_LDOUBLE_WORDS64(x,hx,lx); | |
| 96 | return x; | |
| 97 | } | |
| 98 | if (k <= -54) | |
| 99 | return tiny*__copysignl(tiny,x); /*underflow*/ | |
| 100 | k += 54; /* subnormal result */ | |
| 101 | lx &= 0x8000000000000000ULL; | |
| 102 | SET_LDOUBLE_WORDS64(x,(hx&0x800fffffffffffffULL)|(k<<52),lx); | |
| 103 | return x*twolm54; | |
| 104 | } | |
| 105 | #ifdef IS_IN_libm | |
| 106 | long_double_symbol (libm, __scalbnl, scalbnl); | |
| 107 | #else | |
| 108 | long_double_symbol (libc, __scalbnl, scalbnl); | |
| 109 | #endif |