1 /*
2 * ====================================================
3 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
4 *
5 * Developed at SunPro, a Sun Microsystems, Inc. business.
6 * Permission to use, copy, modify, and distribute this
7 * software is freely granted, provided that this notice
8 * is preserved.
9 * ====================================================
10 */
11
12 /*
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.
17 */
18
19 #include <math.h>
20 #include <math_private.h>
21
22 static const double
23 two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
24 twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
25 huge = 1.0e+300,
26 tiny = 1.0e-300;
27
28 double
__scalbln(double x,long int n)29 __scalbln (double x, long int n)
30 {
31 int64_t ix;
32 int64_t k;
33 EXTRACT_WORDS64(ix,x);
34 k = (ix >> 52) & 0x7ff; /* extract exponent */
35 if (__builtin_expect(k==0, 0)) { /* 0 or subnormal x */
36 if ((ix & UINT64_C(0xfffffffffffff))==0) return x; /* +-0 */
37 x *= two54;
38 EXTRACT_WORDS64(ix,x);
39 k = ((ix >> 52) & 0x7ff) - 54;
40 }
41 if (__builtin_expect(k==0x7ff, 0)) return x+x; /* NaN or Inf */
42 if (__builtin_expect(n< -50000, 0))
43 return tiny*copysign(tiny,x); /*underflow*/
44 if (__builtin_expect(n> 50000 || k+n > 0x7fe, 0))
45 return huge*copysign(huge,x); /* overflow */
46 /* Now k and n are bounded we know that k = k+n does not
47 overflow. */
48 k = k+n;
49 if (__builtin_expect(k > 0, 1)) /* normal result */
50 {INSERT_WORDS64(x,(ix&UINT64_C(0x800fffffffffffff))|(k<<52));
51 return x;}
52 if (k <= -54)
53 return tiny*copysign(tiny,x); /*underflow*/
54 k += 54; /* subnormal result */
55 INSERT_WORDS64(x,(ix&INT64_C(0x800fffffffffffff))|(k<<52));
56 return x*twom54;
57 }
58 #ifdef NO_LONG_DOUBLE
59 strong_alias (__scalbln, __scalblnl)
60 #endif
61