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 /* 23 * Copyright 2011 Nexenta Systems, Inc. All rights reserved. 24 */ 25 /* 26 * Copyright 2006 Sun Microsystems, Inc. All rights reserved. 27 * Use is subject to license terms. 28 */ 29 30 #if defined(ELFOBJ) 31 #pragma weak frexp = __frexp 32 #endif 33 34 /* 35 * frexp(x, exp) returns the normalized significand of x and sets 36 * *exp so that x = r*2^(*exp) where r is the return value. If x 37 * is finite and nonzero, 1/2 <= |r| < 1. 38 * 39 * If x is zero, infinite or NaN, frexp returns x and sets *exp = 0. 40 * (The relevant standards do not specify *exp when x is infinite or 41 * NaN, but this code sets it anyway.) 42 * 43 * If x is a signaling NaN, this code returns x without attempting 44 * to raise the invalid operation exception. If x is subnormal, 45 * this code treats it as nonzero regardless of nonstandard mode. 46 */ 47 48 #include "libm.h" 49 50 double 51 __frexp(double x, int *exp) { 52 union { 53 unsigned i[2]; 54 double d; 55 } xx, yy; 56 double t; 57 unsigned hx; 58 int e; 59 60 xx.d = x; 61 hx = xx.i[HIWORD] & ~0x80000000; 62 63 if (hx >= 0x7ff00000) { /* x is infinite or NaN */ 64 *exp = 0; 65 return (x); 66 } 67 68 e = 0; 69 if (hx < 0x00100000) { /* x is subnormal or zero */ 70 if ((hx | xx.i[LOWORD]) == 0) { 71 *exp = 0; 72 return (x); 73 } 74 75 /* 76 * normalize x by regarding it as an integer 77 * 78 * Here we use 32-bit integer arithmetic to avoid trapping 79 * or emulating 64-bit arithmetic. If 64-bit arithmetic is 80 * available (e.g., in SPARC V9), do this instead: 81 * 82 * long lx = ((long) hx << 32) | xx.i[LOWORD]; 83 * xx.d = (xx.i[HIWORD] < 0)? -lx : lx; 84 * 85 * If subnormal arithmetic doesn't trap, just multiply x by 86 * a power of two. 87 */ 88 yy.i[HIWORD] = 0x43300000 | hx; 89 yy.i[LOWORD] = xx.i[LOWORD]; 90 t = yy.d; 91 yy.i[HIWORD] = 0x43300000; 92 yy.i[LOWORD] = 0; 93 t -= yy.d; /* t = |x| scaled */ 94 xx.d = ((int)xx.i[HIWORD] < 0)? -t : t; 95 hx = xx.i[HIWORD] & ~0x80000000; 96 e = -1074; 97 } 98 99 /* now xx.d is normal */ 100 xx.i[HIWORD] = (xx.i[HIWORD] & ~0x7ff00000) | 0x3fe00000; 101 *exp = e + (hx >> 20) - 0x3fe; 102 return (xx.d); 103 } 104