1 /* 2 * CDDL HEADER START 3 * 4 * The contents of this file are subject to the terms of the 5 * Common Development and Distribution License (the "License"). 6 * You may not use this file except in compliance with the License. 7 * 8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE 9 * or http://www.opensolaris.org/os/licensing. 10 * See the License for the specific language governing permissions 11 * and limitations under the License. 12 * 13 * When distributing Covered Code, include this CDDL HEADER in each 14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE. 15 * If applicable, add the following below this CDDL HEADER, with the 16 * fields enclosed by brackets "[]" replaced with your own identifying 17 * information: Portions Copyright [yyyy] [name of copyright owner] 18 * 19 * CDDL HEADER END 20 */ 21 /* 22 * Copyright 2011 Nexenta Systems, Inc. All rights reserved. 23 */ 24 /* 25 * Copyright 2006 Sun Microsystems, Inc. All rights reserved. 26 * Use is subject to license terms. 27 */ 28 29 #pragma weak scalbn = __scalbn 30 31 #include "libm.h" 32 33 static const double 34 one = 1.0, 35 huge = 1.0e300, 36 tiny = 1.0e-300, 37 twom54 = 5.5511151231257827021181583404541015625e-17; 38 39 #if defined(USE_FPSCALE) || defined(__x86) 40 static const double two52 = 4503599627370496.0; 41 #else 42 /* 43 * Normalize non-zero subnormal x and return biased exponent of x in [-51,0] 44 */ 45 static int 46 ilogb_biased(unsigned *px) { 47 int s = 52; 48 unsigned v = px[HIWORD] & ~0x80000000, w = px[LOWORD], t = v; 49 50 if (t) 51 s -= 32; 52 else 53 t = w; 54 if (t & 0xffff0000) 55 s -= 16, t >>= 16; 56 if (t & 0xff00) 57 s -= 8, t >>= 8; 58 if (t & 0xf0) 59 s -= 4, t >>= 4; 60 t <<= 1; 61 s -= (0xffffaa50 >> t) & 0x3; 62 if (s < 32) { 63 v = (v << s) | w >> (32 - s); 64 w <<= s; 65 } else { 66 v = w << (s - 32); 67 w = 0; 68 } 69 px[HIWORD] = (px[HIWORD] & 0x80000000) | v; 70 px[LOWORD] = w; 71 return (1 - s); 72 } 73 #endif /* defined(USE_FPSCALE) */ 74 75 double 76 scalbn(double x, int n) { 77 int *px, ix, hx, k; 78 79 px = (int *)&x; 80 ix = px[HIWORD]; 81 hx = ix & ~0x80000000; 82 k = hx >> 20; 83 84 if (k == 0x7ff) /* x is inf or NaN */ 85 return (x * one); 86 87 if (k == 0) { 88 if ((hx | px[LOWORD]) == 0 || n == 0) 89 return (x); 90 #if defined(USE_FPSCALE) || defined(__x86) 91 x *= two52; 92 ix = px[HIWORD]; 93 k = ((ix & ~0x80000000) >> 20) - 52; 94 #else 95 k = ilogb_biased((unsigned *)px); 96 ix = px[HIWORD]; 97 #endif 98 /* now k is in the range -51..0 */ 99 k += n; 100 if (k > n) /* integer overflow occurred */ 101 k = -100; 102 } else { 103 /* k is in the range 1..1023 */ 104 k += n; 105 if (k < n) /* integer overflow occurred */ 106 k = 0x7ff; 107 } 108 109 if (k > 0x7fe) 110 return (huge * ((ix < 0)? -huge : huge)); 111 if (k < 1) { 112 if (k <= -54) 113 return (tiny * ((ix < 0)? -tiny : tiny)); 114 k += 54; 115 px[HIWORD] = (ix & ~0x7ff00000) | (k << 20); 116 return (x * twom54); 117 } 118 px[HIWORD] = (ix & ~0x7ff00000) | (k << 20); 119 return (x); 120 } 121