Print this page
11210 libm should be cstyle(1ONBLD) clean
| Split |
Close |
| Expand all |
| Collapse all |
--- old/usr/src/lib/libm/common/Q/fmodl.c
+++ new/usr/src/lib/libm/common/Q/fmodl.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 __fmodl = fmodl
31 32
32 33 #include "libm.h"
33 34
34 -static const int
35 - is = -0x7fffffff - 1,
36 - im = 0x0000ffff,
37 - iu = 0x00010000;
38 -
39 -static const long double
40 - zero = 0.0L,
41 - one = 1.0L;
35 +static const int is = -0x7fffffff - 1, im = 0x0000ffff, iu = 0x00010000;
36 +static const long double zero = 0.0L, one = 1.0L;
42 37
43 38 #ifdef __LITTLE_ENDIAN
44 -#define __H0(x) *(3 + (int *) &x)
45 -#define __H1(x) *(2 + (int *) &x)
46 -#define __H2(x) *(1 + (int *) &x)
47 -#define __H3(x) *(0 + (int *) &x)
39 +#define __H0(x) *(3 + (int *)&x)
40 +#define __H1(x) *(2 + (int *)&x)
41 +#define __H2(x) *(1 + (int *)&x)
42 +#define __H3(x) *(0 + (int *)&x)
48 43 #else
49 -#define __H0(x) *(0 + (int *) &x)
50 -#define __H1(x) *(1 + (int *) &x)
51 -#define __H2(x) *(2 + (int *) &x)
52 -#define __H3(x) *(3 + (int *) &x)
44 +#define __H0(x) *(0 + (int *)&x)
45 +#define __H1(x) *(1 + (int *)&x)
46 +#define __H2(x) *(2 + (int *)&x)
47 +#define __H3(x) *(3 + (int *)&x)
53 48 #endif
54 49
55 50 long double
56 -fmodl(long double x, long double y) {
51 +fmodl(long double x, long double y)
52 +{
57 53 long double a, b;
58 54 int n, ix, iy, k, sx;
59 55 int hx;
60 56 int x0, y0, z0, carry;
61 57 unsigned x1, x2, x3, y1, y2, y3, z1, z2, z3;
62 58
63 59 hx = __H0(x);
64 60 x1 = __H1(x);
65 61 x2 = __H2(x);
66 62 x3 = __H3(x);
67 63 y0 = __H0(y);
68 64 y1 = __H1(y);
69 65 y2 = __H2(y);
70 66 y3 = __H3(y);
71 67
72 68 sx = hx & 0x80000000;
73 69 x0 = hx ^ sx;
74 70 y0 &= 0x7fffffff;
75 71
76 72 /* purge off exception values */
77 - if (x0 >= 0x7fff0000 || /* !finitel(x) */
73 + if (x0 >= 0x7fff0000 || /* !finitel(x) */
78 74 (y0 > 0x7fff0000) || (y0 == 0x7fff0000 && ((y1 | y2 | y3) != 0)) ||
79 75 (y0 | y1 | y2 | y3) == 0) /* isnanl(y) || y = 0 */
80 76 return ((x * y) / (x * y));
77 +
81 78 a = fabsl(x);
82 79 b = fabsl(y);
80 +
83 81 if (a <= b) {
84 82 if (a < b)
85 83 return (x);
86 84 else
87 85 return (zero * x);
88 86 }
87 +
89 88 /* determine ix = ilogbl(x) */
90 - if (x0 < iu) { /* subnormal x */
89 + if (x0 < iu) { /* subnormal x */
91 90 ix = -16382;
91 +
92 92 while (x0 == 0) {
93 93 ix -= 16;
94 94 x0 = x1 >> 16;
95 95 x1 = (x1 << 16) | (x2 >> 16);
96 96 x2 = (x2 << 16) | (x3 >> 16);
97 97 x3 = (x3 << 16);
98 98 }
99 +
99 100 while (x0 < iu) {
100 101 ix -= 1;
101 102 x0 = (x0 << 1) | (x1 >> 31);
102 103 x1 = (x1 << 1) | (x2 >> 31);
103 104 x2 = (x2 << 1) | (x3 >> 31);
104 105 x3 <<= 1;
105 106 }
106 107 } else {
107 108 ix = (x0 >> 16) - 16383;
108 109 x0 = iu | (x0 & im);
109 110 }
110 111
111 112 /* determine iy = ilogbl(y) */
112 - if (y0 < iu) { /* subnormal y */
113 + if (y0 < iu) { /* subnormal y */
113 114 iy = -16382;
115 +
114 116 while (y0 == 0) {
115 117 iy -= 16;
116 118 y0 = y1 >> 16;
117 119 y1 = (y1 << 16) | (y2 >> 16);
118 120 y2 = (y2 << 16) | (y3 >> 16);
119 121 y3 = (y3 << 16);
120 122 }
123 +
121 124 while (y0 < iu) {
122 125 iy -= 1;
123 126 y0 = (y0 << 1) | (y1 >> 31);
124 127 y1 = (y1 << 1) | (y2 >> 31);
125 128 y2 = (y2 << 1) | (y3 >> 31);
126 129 y3 <<= 1;
127 130 }
128 131 } else {
129 132 iy = (y0 >> 16) - 16383;
130 133 y0 = iu | (y0 & im);
131 134 }
132 135
133 136 /* fix point fmod */
134 137 n = ix - iy;
138 +
135 139 while (n--) {
136 140 while (x0 == 0 && n >= 16) {
137 141 n -= 16;
138 142 x0 = x1 >> 16;
139 143 x1 = (x1 << 16) | (x2 >> 16);
140 144 x2 = (x2 << 16) | (x3 >> 16);
141 145 x3 = (x3 << 16);
142 146 }
147 +
143 148 while (x0 < iu && n >= 1) {
144 149 n -= 1;
145 150 x0 = (x0 << 1) | (x1 >> 31);
146 151 x1 = (x1 << 1) | (x2 >> 31);
147 152 x2 = (x2 << 1) | (x3 >> 31);
148 153 x3 = (x3 << 1);
149 154 }
155 +
150 156 carry = 0;
151 157 z3 = x3 - y3;
152 158 carry = (z3 > x3);
159 +
153 160 if (carry == 0) {
154 161 z2 = x2 - y2;
155 162 carry = (z2 > x2);
156 163 } else {
157 164 z2 = x2 - y2 - 1;
158 165 carry = (z2 >= x2);
159 166 }
167 +
160 168 if (carry == 0) {
161 169 z1 = x1 - y1;
162 170 carry = (z1 > x1);
163 171 } else {
164 172 z1 = x1 - y1 - 1;
165 173 carry = (z1 >= x1);
166 174 }
175 +
167 176 z0 = x0 - y0 - carry;
168 - if (z0 < 0) { /* double x */
177 +
178 + if (z0 < 0) { /* double x */
169 179 x0 = x0 + x0 + ((x1 & is) != 0);
170 180 x1 = x1 + x1 + ((x2 & is) != 0);
171 181 x2 = x2 + x2 + ((x3 & is) != 0);
172 182 x3 = x3 + x3;
173 183 } else {
174 184 if (z0 == 0) {
175 185 if ((z1 | z2 | z3) == 0) { /* 0: done */
176 186 __H0(a) = hx & is;
177 187 __H1(a) = __H2(a) = __H3(a) = 0;
178 188 return (a);
179 189 }
180 190 }
191 +
181 192 /* x = z << 1 */
182 193 z0 = z0 + z0 + ((z1 & is) != 0);
183 194 z1 = z1 + z1 + ((z2 & is) != 0);
184 195 z2 = z2 + z2 + ((z3 & is) != 0);
185 196 z3 = z3 + z3;
186 197 x0 = z0;
187 198 x1 = z1;
188 199 x2 = z2;
189 200 x3 = z3;
190 201 }
191 202 }
192 203
193 204 carry = 0;
194 205 z3 = x3 - y3;
195 206 carry = (z3 > x3);
207 +
196 208 if (carry == 0) {
197 209 z2 = x2 - y2;
198 210 carry = (z2 > x2);
199 211 } else {
200 212 z2 = x2 - y2 - 1;
201 213 carry = (z2 >= x2);
202 214 }
215 +
203 216 if (carry == 0) {
204 217 z1 = x1 - y1;
205 218 carry = (z1 > x1);
206 219 } else {
207 220 z1 = x1 - y1 - 1;
208 221 carry = (z1 >= x1);
209 222 }
223 +
210 224 z0 = x0 - y0 - carry;
225 +
211 226 if (z0 >= 0) {
212 227 x0 = z0;
213 228 x1 = z1;
214 229 x2 = z2;
215 230 x3 = z3;
216 231 }
232 +
217 233 /* convert back to floating value and restore the sign */
218 234 if ((x0 | x1 | x2 | x3) == 0) {
219 235 __H0(a) = hx & is;
220 236 __H1(a) = __H2(a) = __H3(a) = 0;
221 237 return (a);
222 238 }
239 +
223 240 while (x0 < iu) {
224 241 if (x0 == 0) {
225 242 iy -= 16;
226 243 x0 = x1 >> 16;
227 244 x1 = (x1 << 16) | (x2 >> 16);
228 245 x2 = (x2 << 16) | (x3 >> 16);
229 246 x3 = (x3 << 16);
230 247 } else {
231 248 x0 = x0 + x0 + ((x1 & is) != 0);
232 249 x1 = x1 + x1 + ((x2 & is) != 0);
233 250 x2 = x2 + x2 + ((x3 & is) != 0);
234 251 x3 = x3 + x3;
|
↓ open down ↓ |
2 lines elided |
↑ open up ↑ |
235 252 iy -= 1;
236 253 }
237 254 }
238 255
239 256 /* normalize output */
240 257 if (iy >= -16382) {
241 258 __H0(a) = sx | (x0 - iu) | ((iy + 16383) << 16);
242 259 __H1(a) = x1;
243 260 __H2(a) = x2;
244 261 __H3(a) = x3;
245 - } else { /* subnormal output */
262 + } else { /* subnormal output */
246 263 n = -16382 - iy;
247 264 k = n & 31;
265 +
248 266 if (k != 0) {
249 267 if (k <= 16) {
250 268 x3 = (x2 << (32 - k)) | (x3 >> k);
251 269 x2 = (x1 << (32 - k)) | (x2 >> k);
252 270 x1 = (x0 << (32 - k)) | (x1 >> k);
253 271 x0 >>= k;
254 272 } else {
255 273 x3 = (x2 << (32 - k)) | (x3 >> k);
256 274 x2 = (x1 << (32 - k)) | (x2 >> k);
257 275 x1 = (x0 << (32 - k)) | (x1 >> k);
258 276 x0 = 0;
259 277 }
260 278 }
279 +
261 280 while (n >= 32) {
262 281 n -= 32;
263 282 x3 = x2;
264 283 x2 = x1;
265 284 x1 = x0;
266 285 x0 = 0;
267 286 }
287 +
268 288 __H0(a) = x0 | sx;
269 289 __H1(a) = x1;
270 290 __H2(a) = x2;
271 291 __H3(a) = x3;
272 292 a *= one;
273 293 }
294 +
274 295 return (a);
275 296 }
XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX