Print this page
5261 libm should stop using synonyms.h
5298 fabs is 0-sized, confuses dis(1) and others
Reviewed by: Josef 'Jeff' Sipek <jeffpc@josefsipek.net>
Approved by: Gordon Ross <gwr@nexenta.com>
| Split |
Close |
| Expand all |
| Collapse all |
--- old/usr/src/lib/libm/common/m9x/remquof.c
+++ new/usr/src/lib/libm/common/m9x/remquof.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
|
↓ open down ↓ |
19 lines elided |
↑ open up ↑ |
20 20 */
21 21
22 22 /*
23 23 * Copyright 2011 Nexenta Systems, Inc. All rights reserved.
24 24 */
25 25 /*
26 26 * Copyright 2006 Sun Microsystems, Inc. All rights reserved.
27 27 * Use is subject to license terms.
28 28 */
29 29
30 -#pragma weak remquof = __remquof
30 +#pragma weak __remquof = remquof
31 31
32 32 /* INDENT OFF */
33 33 /*
34 34 * float remquof(float x, float y, int *quo) return remainderf(x,y) and an
35 35 * integer pointer quo such that *quo = N mod (2**31), where N is the
36 36 * exact integeral part of x/y rounded to nearest even.
37 37 *
38 38 * remquof call internal fmodquof
39 39 */
40 40
41 41 #include "libm.h"
42 -#include "libm_synonyms.h"
43 42 #include "libm_protos.h"
44 43 #include <math.h>
45 44 extern float fabsf(float);
46 45
47 46 static const int
48 47 is = (int) 0x80000000,
49 48 im = 0x007fffff,
50 49 ii = 0x7f800000,
51 50 iu = 0x00800000;
52 51
53 52 static const float zero = 0.0F, half = 0.5F;
54 53 /* INDENT ON */
55 54
56 55 static float
57 56 fmodquof(float x, float y, int *quo) {
58 57 float w;
59 58 int hx, ix, iy, iz, k, ny, nd, m, sq;
60 59
61 60 hx = *(int *) &x;
62 61 ix = hx & 0x7fffffff;
63 62 iy = *(int *) &y;
64 63 sq = (iy ^ hx) & is; /* sign of x/y */
65 64 iy &= 0x7fffffff;
66 65
67 66 /* purge off exception values */
68 67 *quo = 0;
69 68 if (ix >= ii || iy > ii || iy == 0) {
70 69 w = x * y;
71 70 w = w / w;
72 71 } else if (ix <= iy) {
73 72 if (ix < iy)
74 73 w = x; /* return x if |x|<|y| */
75 74 else {
76 75 *quo = 1 + (sq >> 30);
77 76 w = zero * x; /* return sign(x)*0.0 */
78 77 }
79 78 } else {
80 79 /* INDENT OFF */
81 80 /*
82 81 * scale x,y to "normal" with
83 82 * ny = exponent of y
84 83 * nd = exponent of x minus exponent of y
85 84 */
86 85 /* INDENT ON */
87 86 ny = iy >> 23;
88 87 k = ix >> 23;
89 88
90 89 /* special case for subnormal y or x */
91 90 if (ny == 0) {
92 91 ny = 1;
93 92 while (iy < iu) {
94 93 ny -= 1;
95 94 iy += iy;
96 95 }
97 96 nd = k - ny;
98 97 if (k == 0) {
99 98 nd += 1;
100 99 while (ix < iu) {
101 100 nd -= 1;
102 101 ix += ix;
103 102 }
104 103 } else
105 104 ix = iu | (ix & im);
106 105 } else {
107 106 nd = k - ny;
108 107 ix = iu | (ix & im);
109 108 iy = iu | (iy & im);
110 109 }
111 110 /* INDENT OFF */
112 111 /* fix point fmod for normalized ix and iy */
113 112 /*
114 113 * while (nd--) {
115 114 * iz = ix - iy;
116 115 * if (iz < 0)
117 116 * ix = ix + ix;
118 117 * else if (iz == 0) {
119 118 * *(int *) &w = is & hx;
120 119 * return w;
121 120 * } else
122 121 * ix = iz + iz;
123 122 * }
124 123 */
125 124 /* INDENT ON */
126 125 /* unroll the above loop 4 times to gain performance */
127 126 m = 0;
128 127 k = nd >> 2;
129 128 nd -= (k << 2);
130 129 while (k--) {
131 130 iz = ix - iy;
132 131 if (iz >= 0) {
133 132 m += 1;
134 133 ix = iz + iz;
135 134 } else
136 135 ix += ix;
137 136 m += m;
138 137 iz = ix - iy;
139 138 if (iz >= 0) {
140 139 m += 1;
141 140 ix = iz + iz;
142 141 } else
143 142 ix += ix;
144 143 m += m;
145 144 iz = ix - iy;
146 145 if (iz >= 0) {
147 146 m += 1;
148 147 ix = iz + iz;
149 148 } else
150 149 ix += ix;
151 150 m += m;
152 151 iz = ix - iy;
153 152 if (iz >= 0) {
154 153 m += 1;
155 154 ix = iz + iz;
156 155 } else
157 156 ix += ix;
158 157 m += m;
159 158 if (iz == 0) {
160 159 iz = (k << 2) + nd;
161 160 if (iz < 32)
162 161 m <<= iz;
163 162 else
164 163 m = 0;
165 164 m &= 0x7fffffff;
166 165 *quo = sq >= 0 ? m : -m;
167 166 *(int *) &w = is & hx;
168 167 return (w);
169 168 }
170 169 }
171 170 while (nd--) {
172 171 iz = ix - iy;
173 172 if (iz >= 0) {
174 173 m += 1;
175 174 ix = iz + iz;
176 175 } else
177 176 ix += ix;
178 177 m += m;
179 178 }
180 179 /* end of unrolling */
181 180
182 181 iz = ix - iy;
183 182 if (iz >= 0) {
184 183 m += 1;
185 184 ix = iz;
186 185 }
187 186 m &= 0x7fffffff;
188 187 *quo = sq >= 0 ? m : -m;
189 188
190 189 /* convert back to floating value and restore the sign */
191 190 if (ix == 0) {
192 191 *(int *) &w = is & hx;
193 192 return (w);
194 193 }
195 194 while (ix < iu) {
196 195 ix += ix;
197 196 ny -= 1;
198 197 }
199 198 while (ix > (iu + iu)) {
200 199 ny += 1;
201 200 ix >>= 1;
202 201 }
203 202 if (ny > 0)
204 203 *(int *) &w = (is & hx) | (ix & im) | (ny << 23);
205 204 else { /* subnormal output */
206 205 k = -ny + 1;
207 206 ix >>= k;
208 207 *(int *) &w = (is & hx) | ix;
209 208 }
210 209 }
211 210 return (w);
212 211 }
213 212
214 213 float
215 214 remquof(float x, float y, int *quo) {
216 215 int hx, hy, sx, sq;
217 216 float v;
218 217
219 218 hx = *(int *) &x; /* high word of x */
220 219 hy = *(int *) &y; /* high word of y */
221 220 sx = hx & is; /* sign of x */
222 221 sq = (hx ^ hy) & is; /* sign of x/y */
223 222 hx ^= sx; /* |x| */
224 223 hy &= 0x7fffffff; /* |y| */
225 224
226 225 /* purge off exception values: y is 0 or NaN, x is Inf or NaN */
227 226 *quo = 0;
228 227 if (hx >= ii || hy > ii || hy == 0) {
229 228 v = x * y;
230 229 return (v / v);
231 230 }
232 231
233 232 y = fabsf(y);
234 233 x = fabsf(x);
235 234 if (hy <= 0x7f7fffff) {
236 235 x = fmodquof(x, y + y, quo);
237 236 *quo = ((*quo) & 0x3fffffff) << 1;
238 237 }
239 238 if (hy < 0x01000000) {
240 239 if (x + x > y) {
241 240 *quo += 1;
242 241 if (x == y)
243 242 x = zero;
244 243 else
245 244 x -= y;
246 245 if (x + x >= y) {
247 246 x -= y;
248 247 *quo += 1;
249 248 }
250 249 }
251 250 } else {
252 251 v = half * y;
253 252 if (x > v) {
254 253 *quo += 1;
255 254 if (x == y)
256 255 x = zero;
257 256 else
258 257 x -= y;
259 258 if (x >= v) {
260 259 x -= y;
261 260 *quo += 1;
262 261 }
263 262 }
264 263 }
265 264 if (sq != 0)
266 265 *quo = -(*quo);
267 266 return (sx == 0 ? x : -x);
268 267 }
|
↓ open down ↓ |
216 lines elided |
↑ open up ↑ |
XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX