1 /* origin: FreeBSD /usr/src/lib/msun/src/s_scalbn.c */
3 * ====================================================
4 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
6 * Developed at SunPro, a Sun Microsystems, Inc. business.
7 * Permission to use, copy, modify, and distribute this
8 * software is freely granted, provided that this notice
10 * ====================================================
13 * scalbn (double x, int n)
14 * scalbn(x,n) returns x* 2**n computed by exponent
15 * manipulation rather than by actually performing an
16 * exponentiation or a multiplication.
22 two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
23 twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
27 double scalbn(double x, int n)
29 // FIXME: k+n check depends on signed int overflow.. use unsigned hx
30 // TODO: when long != int:
31 // scalbln(x,long n) { if(n>9999)n=9999; else if(n<-9999)n=-9999; return scalbn(x,n); }
32 // TODO: n < -50000 ...
35 EXTRACT_WORDS(hx, lx, x);
36 k = (hx&0x7ff00000)>>20; /* extract exponent */
37 if (k == 0) { /* 0 or subnormal x */
38 if ((lx|(hx&0x7fffffff)) == 0) /* +-0 */
42 k = ((hx&0x7ff00000)>>20) - 54;
44 return tiny*x; /*underflow*/
46 if (k == 0x7ff) /* NaN or Inf */
50 return huge*copysign(huge, x); /* overflow */
51 if (k > 0) { /* normal result */
52 SET_HIGH_WORD(x, (hx&0x800fffff)|(k<<20));
56 if (n > 50000) /* in case integer overflow in n+k */
57 return huge*copysign(huge, x); /*overflow*/
58 return tiny*copysign(tiny, x); /*underflow*/
59 k += 54; /* subnormal result */
60 SET_HIGH_WORD(x, (hx&0x800fffff)|(k<<20));