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