Print this page
5262 libm needs to be carefully unifdef'd
5268 libm doesn't need to hide symbols which are already local
| Split |
Close |
| Expand all |
| Collapse all |
--- old/usr/src/lib/libm/common/C/sqrt.c
+++ new/usr/src/lib/libm/common/C/sqrt.c
1 1 /*
2 2 * CDDL HEADER START
3 3 *
4 4 * The contents of this file are subject to the terms of the
5 5 * Common Development and Distribution License (the "License").
6 6 * You may not use this file except in compliance with the License.
7 7 *
8 8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
9 9 * or http://www.opensolaris.org/os/licensing.
10 10 * See the License for the specific language governing permissions
11 11 * and limitations under the License.
12 12 *
13 13 * When distributing Covered Code, include this CDDL HEADER in each
14 14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
15 15 * If applicable, add the following below this CDDL HEADER, with the
16 16 * fields enclosed by brackets "[]" replaced with your own identifying
17 17 * information: Portions Copyright [yyyy] [name of copyright owner]
18 18 *
19 19 * CDDL HEADER END
20 20 */
21 21 /*
22 22 * Copyright 2011 Nexenta Systems, Inc. All rights reserved.
|
↓ open down ↓ |
22 lines elided |
↑ open up ↑ |
23 23 */
24 24 /*
25 25 * Copyright 2005 Sun Microsystems, Inc. All rights reserved.
26 26 * Use is subject to license terms.
27 27 */
28 28
29 29 #pragma weak sqrt = __sqrt
30 30
31 31 #include "libm.h"
32 32
33 -#ifdef __INLINE
34 33
35 34 extern double __inline_sqrt(double);
36 35
37 36 double
38 37 sqrt(double x) {
39 38 double z = __inline_sqrt(x);
40 39
41 40 if (isnan(x))
42 41 return (z);
43 42 return ((x < 0.0)? _SVID_libm_err(x, x, 26) : z);
44 43 }
45 44
46 -#else /* defined(__INLINE) */
47 -
48 -/*
49 - * Warning: This correctly rounded sqrt is extremely slow because it computes
50 - * the sqrt bit by bit using integer arithmetic.
51 - */
52 -
53 -static const double big = 1.0e30, small = 1.0e-30;
54 -
55 -double
56 -sqrt(double x)
57 -{
58 - double z;
59 - unsigned r, t1, s1, ix1, q1;
60 - int ix0, s0, j, q, m, n, t;
61 - int *px = (int *)&x, *pz = (int *)&z;
62 -
63 - ix0 = px[HIWORD];
64 - ix1 = px[LOWORD];
65 - if ((ix0 & 0x7ff00000) == 0x7ff00000) { /* x is inf or NaN */
66 - if (ix0 == 0xfff00000 && ix1 == 0)
67 - return (_SVID_libm_err(x, x, 26));
68 - return (x + x);
69 - }
70 - if (((ix0 & 0x7fffffff) | ix1) == 0) /* x is zero */
71 - return (x);
72 -
73 - /* extract exponent and significand */
74 - m = ilogb(x);
75 - z = scalbn(x, -m);
76 - ix0 = (pz[HIWORD] & 0x000fffff) | 0x00100000;
77 - ix1 = pz[LOWORD];
78 - n = m >> 1;
79 - if (n + n != m) {
80 - ix0 = (ix0 << 1) | (ix1 >> 31);
81 - ix1 <<= 1;
82 - m -= 1;
83 - }
84 -
85 - /* generate sqrt(x) bit by bit */
86 - ix0 = (ix0 << 1) | (ix1 >> 31);
87 - ix1 <<= 1;
88 - q = q1 = s0 = s1 = 0;
89 - r = 0x00200000;
90 -
91 - for (j = 1; j <= 22; j++) {
92 - t = s0 + r;
93 - if (t <= ix0) {
94 - s0 = t + r;
95 - ix0 -= t;
96 - q += r;
97 - }
98 - ix0 = (ix0 << 1) | (ix1 >> 31);
99 - ix1 <<= 1;
100 - r >>= 1;
101 - }
102 -
103 - r = 0x80000000;
104 - for (j = 1; j <= 32; j++) {
105 - t1 = s1 + r;
106 - t = s0;
107 - if (t < ix0 || (t == ix0 && t1 <= ix1)) {
108 - s1 = t1 + r;
109 - if ((t1 & 0x80000000) == 0x80000000 &&
110 - (s1 & 0x80000000) == 0)
111 - s0 += 1;
112 - ix0 -= t;
113 - if (ix1 < t1)
114 - ix0 -= 1;
115 - ix1 -= t1;
116 - q1 += r;
117 - }
118 - ix0 = (ix0 << 1) | (ix1 >> 31);
119 - ix1 <<= 1;
120 - r >>= 1;
121 - }
122 -
123 - /* round */
124 - if ((ix0 | ix1) == 0)
125 - goto done;
126 - z = big - small; /* trigger inexact flag */
127 - if (z < big)
128 - goto done;
129 - if (q1 == 0xffffffff) {
130 - q1 = 0;
131 - q += 1;
132 - goto done;
133 - }
134 - z = big + small;
135 - if (z > big) {
136 - if (q1 == 0xfffffffe)
137 - q += 1;
138 - q1 += 2;
139 - goto done;
140 - }
141 - q1 += (q1 & 1);
142 -done:
143 - pz[HIWORD] = (q >> 1) + 0x3fe00000;
144 - pz[LOWORD] = q1 >> 1;
145 - if ((q & 1) == 1)
146 - pz[LOWORD] |= 0x80000000;
147 - return (scalbn(z, n));
148 -}
149 -
150 -#endif /* defined(__INLINE) */
XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX