0
|
1 /* Software floating-point emulation. Common operations.
|
|
2 Copyright (C) 1997,1998,1999,2006,2007 Free Software Foundation, Inc.
|
|
3 This file is part of the GNU C Library.
|
|
4 Contributed by Richard Henderson (rth@cygnus.com),
|
|
5 Jakub Jelinek (jj@ultra.linux.cz),
|
|
6 David S. Miller (davem@redhat.com) and
|
|
7 Peter Maydell (pmaydell@chiark.greenend.org.uk).
|
|
8
|
|
9 The GNU C Library is free software; you can redistribute it and/or
|
|
10 modify it under the terms of the GNU Lesser General Public
|
|
11 License as published by the Free Software Foundation; either
|
|
12 version 2.1 of the License, or (at your option) any later version.
|
|
13
|
|
14 In addition to the permissions in the GNU Lesser General Public
|
|
15 License, the Free Software Foundation gives you unlimited
|
|
16 permission to link the compiled version of this file into
|
|
17 combinations with other programs, and to distribute those
|
|
18 combinations without any restriction coming from the use of this
|
|
19 file. (The Lesser General Public License restrictions do apply in
|
|
20 other respects; for example, they cover modification of the file,
|
|
21 and distribution when not linked into a combine executable.)
|
|
22
|
|
23 The GNU C Library is distributed in the hope that it will be useful,
|
|
24 but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
25 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
|
|
26 Lesser General Public License for more details.
|
|
27
|
|
28 You should have received a copy of the GNU Lesser General Public
|
|
29 License along with the GNU C Library; if not, write to the Free
|
|
30 Software Foundation, 51 Franklin Street, Fifth Floor, Boston,
|
|
31 MA 02110-1301, USA. */
|
|
32
|
|
33 #define _FP_DECL(wc, X) \
|
|
34 _FP_I_TYPE X##_c __attribute__((unused)), X##_s, X##_e; \
|
|
35 _FP_FRAC_DECL_##wc(X)
|
|
36
|
|
37 /*
|
|
38 * Finish truely unpacking a native fp value by classifying the kind
|
|
39 * of fp value and normalizing both the exponent and the fraction.
|
|
40 */
|
|
41
|
|
42 #define _FP_UNPACK_CANONICAL(fs, wc, X) \
|
|
43 do { \
|
|
44 switch (X##_e) \
|
|
45 { \
|
|
46 default: \
|
|
47 _FP_FRAC_HIGH_RAW_##fs(X) |= _FP_IMPLBIT_##fs; \
|
|
48 _FP_FRAC_SLL_##wc(X, _FP_WORKBITS); \
|
|
49 X##_e -= _FP_EXPBIAS_##fs; \
|
|
50 X##_c = FP_CLS_NORMAL; \
|
|
51 break; \
|
|
52 \
|
|
53 case 0: \
|
|
54 if (_FP_FRAC_ZEROP_##wc(X)) \
|
|
55 X##_c = FP_CLS_ZERO; \
|
|
56 else \
|
|
57 { \
|
|
58 /* a denormalized number */ \
|
|
59 _FP_I_TYPE _shift; \
|
|
60 _FP_FRAC_CLZ_##wc(_shift, X); \
|
|
61 _shift -= _FP_FRACXBITS_##fs; \
|
|
62 _FP_FRAC_SLL_##wc(X, (_shift+_FP_WORKBITS)); \
|
|
63 X##_e -= _FP_EXPBIAS_##fs - 1 + _shift; \
|
|
64 X##_c = FP_CLS_NORMAL; \
|
|
65 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
66 } \
|
|
67 break; \
|
|
68 \
|
|
69 case _FP_EXPMAX_##fs: \
|
|
70 if (_FP_FRAC_ZEROP_##wc(X)) \
|
|
71 X##_c = FP_CLS_INF; \
|
|
72 else \
|
|
73 { \
|
|
74 X##_c = FP_CLS_NAN; \
|
|
75 /* Check for signaling NaN */ \
|
|
76 if (!(_FP_FRAC_HIGH_RAW_##fs(X) & _FP_QNANBIT_##fs)) \
|
|
77 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
78 } \
|
|
79 break; \
|
|
80 } \
|
|
81 } while (0)
|
|
82
|
|
83 /* Finish unpacking an fp value in semi-raw mode: the mantissa is
|
|
84 shifted by _FP_WORKBITS but the implicit MSB is not inserted and
|
|
85 other classification is not done. */
|
|
86 #define _FP_UNPACK_SEMIRAW(fs, wc, X) _FP_FRAC_SLL_##wc(X, _FP_WORKBITS)
|
|
87
|
|
88 /* A semi-raw value has overflowed to infinity. Adjust the mantissa
|
|
89 and exponent appropriately. */
|
|
90 #define _FP_OVERFLOW_SEMIRAW(fs, wc, X) \
|
|
91 do { \
|
|
92 if (FP_ROUNDMODE == FP_RND_NEAREST \
|
|
93 || (FP_ROUNDMODE == FP_RND_PINF && !X##_s) \
|
|
94 || (FP_ROUNDMODE == FP_RND_MINF && X##_s)) \
|
|
95 { \
|
|
96 X##_e = _FP_EXPMAX_##fs; \
|
|
97 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
|
|
98 } \
|
|
99 else \
|
|
100 { \
|
|
101 X##_e = _FP_EXPMAX_##fs - 1; \
|
|
102 _FP_FRAC_SET_##wc(X, _FP_MAXFRAC_##wc); \
|
|
103 } \
|
|
104 FP_SET_EXCEPTION(FP_EX_INEXACT); \
|
|
105 FP_SET_EXCEPTION(FP_EX_OVERFLOW); \
|
|
106 } while (0)
|
|
107
|
|
108 /* Check for a semi-raw value being a signaling NaN and raise the
|
|
109 invalid exception if so. */
|
|
110 #define _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X) \
|
|
111 do { \
|
|
112 if (X##_e == _FP_EXPMAX_##fs \
|
|
113 && !_FP_FRAC_ZEROP_##wc(X) \
|
|
114 && !(_FP_FRAC_HIGH_##fs(X) & _FP_QNANBIT_SH_##fs)) \
|
|
115 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
116 } while (0)
|
|
117
|
|
118 /* Choose a NaN result from an operation on two semi-raw NaN
|
|
119 values. */
|
|
120 #define _FP_CHOOSENAN_SEMIRAW(fs, wc, R, X, Y, OP) \
|
|
121 do { \
|
|
122 /* _FP_CHOOSENAN expects raw values, so shift as required. */ \
|
|
123 _FP_FRAC_SRL_##wc(X, _FP_WORKBITS); \
|
|
124 _FP_FRAC_SRL_##wc(Y, _FP_WORKBITS); \
|
|
125 _FP_CHOOSENAN(fs, wc, R, X, Y, OP); \
|
|
126 _FP_FRAC_SLL_##wc(R, _FP_WORKBITS); \
|
|
127 } while (0)
|
|
128
|
|
129 /* Test whether a biased exponent is normal (not zero or maximum). */
|
|
130 #define _FP_EXP_NORMAL(fs, wc, X) (((X##_e + 1) & _FP_EXPMAX_##fs) > 1)
|
|
131
|
|
132 /* Prepare to pack an fp value in semi-raw mode: the mantissa is
|
|
133 rounded and shifted right, with the rounding possibly increasing
|
|
134 the exponent (including changing a finite value to infinity). */
|
|
135 #define _FP_PACK_SEMIRAW(fs, wc, X) \
|
|
136 do { \
|
|
137 _FP_ROUND(wc, X); \
|
|
138 if (_FP_FRAC_HIGH_##fs(X) \
|
|
139 & (_FP_OVERFLOW_##fs >> 1)) \
|
|
140 { \
|
|
141 _FP_FRAC_HIGH_##fs(X) &= ~(_FP_OVERFLOW_##fs >> 1); \
|
|
142 X##_e++; \
|
|
143 if (X##_e == _FP_EXPMAX_##fs) \
|
|
144 _FP_OVERFLOW_SEMIRAW(fs, wc, X); \
|
|
145 } \
|
|
146 _FP_FRAC_SRL_##wc(X, _FP_WORKBITS); \
|
|
147 if (!_FP_EXP_NORMAL(fs, wc, X) && !_FP_FRAC_ZEROP_##wc(X)) \
|
|
148 { \
|
|
149 if (X##_e == 0) \
|
|
150 FP_SET_EXCEPTION(FP_EX_UNDERFLOW); \
|
|
151 else \
|
|
152 { \
|
|
153 if (!_FP_KEEPNANFRACP) \
|
|
154 { \
|
|
155 _FP_FRAC_SET_##wc(X, _FP_NANFRAC_##fs); \
|
|
156 X##_s = _FP_NANSIGN_##fs; \
|
|
157 } \
|
|
158 else \
|
|
159 _FP_FRAC_HIGH_RAW_##fs(X) |= _FP_QNANBIT_##fs; \
|
|
160 } \
|
|
161 } \
|
|
162 } while (0)
|
|
163
|
|
164 /*
|
|
165 * Before packing the bits back into the native fp result, take care
|
|
166 * of such mundane things as rounding and overflow. Also, for some
|
|
167 * kinds of fp values, the original parts may not have been fully
|
|
168 * extracted -- but that is ok, we can regenerate them now.
|
|
169 */
|
|
170
|
|
171 #define _FP_PACK_CANONICAL(fs, wc, X) \
|
|
172 do { \
|
|
173 switch (X##_c) \
|
|
174 { \
|
|
175 case FP_CLS_NORMAL: \
|
|
176 X##_e += _FP_EXPBIAS_##fs; \
|
|
177 if (X##_e > 0) \
|
|
178 { \
|
|
179 _FP_ROUND(wc, X); \
|
|
180 if (_FP_FRAC_OVERP_##wc(fs, X)) \
|
|
181 { \
|
|
182 _FP_FRAC_CLEAR_OVERP_##wc(fs, X); \
|
|
183 X##_e++; \
|
|
184 } \
|
|
185 _FP_FRAC_SRL_##wc(X, _FP_WORKBITS); \
|
|
186 if (X##_e >= _FP_EXPMAX_##fs) \
|
|
187 { \
|
|
188 /* overflow */ \
|
|
189 switch (FP_ROUNDMODE) \
|
|
190 { \
|
|
191 case FP_RND_NEAREST: \
|
|
192 X##_c = FP_CLS_INF; \
|
|
193 break; \
|
|
194 case FP_RND_PINF: \
|
|
195 if (!X##_s) X##_c = FP_CLS_INF; \
|
|
196 break; \
|
|
197 case FP_RND_MINF: \
|
|
198 if (X##_s) X##_c = FP_CLS_INF; \
|
|
199 break; \
|
|
200 } \
|
|
201 if (X##_c == FP_CLS_INF) \
|
|
202 { \
|
|
203 /* Overflow to infinity */ \
|
|
204 X##_e = _FP_EXPMAX_##fs; \
|
|
205 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
|
|
206 } \
|
|
207 else \
|
|
208 { \
|
|
209 /* Overflow to maximum normal */ \
|
|
210 X##_e = _FP_EXPMAX_##fs - 1; \
|
|
211 _FP_FRAC_SET_##wc(X, _FP_MAXFRAC_##wc); \
|
|
212 } \
|
|
213 FP_SET_EXCEPTION(FP_EX_OVERFLOW); \
|
|
214 FP_SET_EXCEPTION(FP_EX_INEXACT); \
|
|
215 } \
|
|
216 } \
|
|
217 else \
|
|
218 { \
|
|
219 /* we've got a denormalized number */ \
|
|
220 X##_e = -X##_e + 1; \
|
|
221 if (X##_e <= _FP_WFRACBITS_##fs) \
|
|
222 { \
|
|
223 _FP_FRAC_SRS_##wc(X, X##_e, _FP_WFRACBITS_##fs); \
|
|
224 _FP_ROUND(wc, X); \
|
|
225 if (_FP_FRAC_HIGH_##fs(X) \
|
|
226 & (_FP_OVERFLOW_##fs >> 1)) \
|
|
227 { \
|
|
228 X##_e = 1; \
|
|
229 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
|
|
230 } \
|
|
231 else \
|
|
232 { \
|
|
233 X##_e = 0; \
|
|
234 _FP_FRAC_SRL_##wc(X, _FP_WORKBITS); \
|
|
235 FP_SET_EXCEPTION(FP_EX_UNDERFLOW); \
|
|
236 } \
|
|
237 } \
|
|
238 else \
|
|
239 { \
|
|
240 /* underflow to zero */ \
|
|
241 X##_e = 0; \
|
|
242 if (!_FP_FRAC_ZEROP_##wc(X)) \
|
|
243 { \
|
|
244 _FP_FRAC_SET_##wc(X, _FP_MINFRAC_##wc); \
|
|
245 _FP_ROUND(wc, X); \
|
|
246 _FP_FRAC_LOW_##wc(X) >>= (_FP_WORKBITS); \
|
|
247 } \
|
|
248 FP_SET_EXCEPTION(FP_EX_UNDERFLOW); \
|
|
249 } \
|
|
250 } \
|
|
251 break; \
|
|
252 \
|
|
253 case FP_CLS_ZERO: \
|
|
254 X##_e = 0; \
|
|
255 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
|
|
256 break; \
|
|
257 \
|
|
258 case FP_CLS_INF: \
|
|
259 X##_e = _FP_EXPMAX_##fs; \
|
|
260 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
|
|
261 break; \
|
|
262 \
|
|
263 case FP_CLS_NAN: \
|
|
264 X##_e = _FP_EXPMAX_##fs; \
|
|
265 if (!_FP_KEEPNANFRACP) \
|
|
266 { \
|
|
267 _FP_FRAC_SET_##wc(X, _FP_NANFRAC_##fs); \
|
|
268 X##_s = _FP_NANSIGN_##fs; \
|
|
269 } \
|
|
270 else \
|
|
271 _FP_FRAC_HIGH_RAW_##fs(X) |= _FP_QNANBIT_##fs; \
|
|
272 break; \
|
|
273 } \
|
|
274 } while (0)
|
|
275
|
|
276 /* This one accepts raw argument and not cooked, returns
|
|
277 * 1 if X is a signaling NaN.
|
|
278 */
|
|
279 #define _FP_ISSIGNAN(fs, wc, X) \
|
|
280 ({ \
|
|
281 int __ret = 0; \
|
|
282 if (X##_e == _FP_EXPMAX_##fs) \
|
|
283 { \
|
|
284 if (!_FP_FRAC_ZEROP_##wc(X) \
|
|
285 && !(_FP_FRAC_HIGH_RAW_##fs(X) & _FP_QNANBIT_##fs)) \
|
|
286 __ret = 1; \
|
|
287 } \
|
|
288 __ret; \
|
|
289 })
|
|
290
|
|
291
|
|
292
|
|
293
|
|
294
|
|
295 /* Addition on semi-raw values. */
|
|
296 #define _FP_ADD_INTERNAL(fs, wc, R, X, Y, OP) \
|
|
297 do { \
|
|
298 if (X##_s == Y##_s) \
|
|
299 { \
|
|
300 /* Addition. */ \
|
|
301 R##_s = X##_s; \
|
|
302 int ediff = X##_e - Y##_e; \
|
|
303 if (ediff > 0) \
|
|
304 { \
|
|
305 R##_e = X##_e; \
|
|
306 if (Y##_e == 0) \
|
|
307 { \
|
|
308 /* Y is zero or denormalized. */ \
|
|
309 if (_FP_FRAC_ZEROP_##wc(Y)) \
|
|
310 { \
|
|
311 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
|
|
312 _FP_FRAC_COPY_##wc(R, X); \
|
|
313 goto add_done; \
|
|
314 } \
|
|
315 else \
|
|
316 { \
|
|
317 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
318 ediff--; \
|
|
319 if (ediff == 0) \
|
|
320 { \
|
|
321 _FP_FRAC_ADD_##wc(R, X, Y); \
|
|
322 goto add3; \
|
|
323 } \
|
|
324 if (X##_e == _FP_EXPMAX_##fs) \
|
|
325 { \
|
|
326 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
|
|
327 _FP_FRAC_COPY_##wc(R, X); \
|
|
328 goto add_done; \
|
|
329 } \
|
|
330 goto add1; \
|
|
331 } \
|
|
332 } \
|
|
333 else if (X##_e == _FP_EXPMAX_##fs) \
|
|
334 { \
|
|
335 /* X is NaN or Inf, Y is normal. */ \
|
|
336 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
|
|
337 _FP_FRAC_COPY_##wc(R, X); \
|
|
338 goto add_done; \
|
|
339 } \
|
|
340 \
|
|
341 /* Insert implicit MSB of Y. */ \
|
|
342 _FP_FRAC_HIGH_##fs(Y) |= _FP_IMPLBIT_SH_##fs; \
|
|
343 \
|
|
344 add1: \
|
|
345 /* Shift the mantissa of Y to the right EDIFF steps; \
|
|
346 remember to account later for the implicit MSB of X. */ \
|
|
347 if (ediff <= _FP_WFRACBITS_##fs) \
|
|
348 _FP_FRAC_SRS_##wc(Y, ediff, _FP_WFRACBITS_##fs); \
|
|
349 else if (!_FP_FRAC_ZEROP_##wc(Y)) \
|
|
350 _FP_FRAC_SET_##wc(Y, _FP_MINFRAC_##wc); \
|
|
351 _FP_FRAC_ADD_##wc(R, X, Y); \
|
|
352 } \
|
|
353 else if (ediff < 0) \
|
|
354 { \
|
|
355 ediff = -ediff; \
|
|
356 R##_e = Y##_e; \
|
|
357 if (X##_e == 0) \
|
|
358 { \
|
|
359 /* X is zero or denormalized. */ \
|
|
360 if (_FP_FRAC_ZEROP_##wc(X)) \
|
|
361 { \
|
|
362 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
|
|
363 _FP_FRAC_COPY_##wc(R, Y); \
|
|
364 goto add_done; \
|
|
365 } \
|
|
366 else \
|
|
367 { \
|
|
368 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
369 ediff--; \
|
|
370 if (ediff == 0) \
|
|
371 { \
|
|
372 _FP_FRAC_ADD_##wc(R, Y, X); \
|
|
373 goto add3; \
|
|
374 } \
|
|
375 if (Y##_e == _FP_EXPMAX_##fs) \
|
|
376 { \
|
|
377 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
|
|
378 _FP_FRAC_COPY_##wc(R, Y); \
|
|
379 goto add_done; \
|
|
380 } \
|
|
381 goto add2; \
|
|
382 } \
|
|
383 } \
|
|
384 else if (Y##_e == _FP_EXPMAX_##fs) \
|
|
385 { \
|
|
386 /* Y is NaN or Inf, X is normal. */ \
|
|
387 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
|
|
388 _FP_FRAC_COPY_##wc(R, Y); \
|
|
389 goto add_done; \
|
|
390 } \
|
|
391 \
|
|
392 /* Insert implicit MSB of X. */ \
|
|
393 _FP_FRAC_HIGH_##fs(X) |= _FP_IMPLBIT_SH_##fs; \
|
|
394 \
|
|
395 add2: \
|
|
396 /* Shift the mantissa of X to the right EDIFF steps; \
|
|
397 remember to account later for the implicit MSB of Y. */ \
|
|
398 if (ediff <= _FP_WFRACBITS_##fs) \
|
|
399 _FP_FRAC_SRS_##wc(X, ediff, _FP_WFRACBITS_##fs); \
|
|
400 else if (!_FP_FRAC_ZEROP_##wc(X)) \
|
|
401 _FP_FRAC_SET_##wc(X, _FP_MINFRAC_##wc); \
|
|
402 _FP_FRAC_ADD_##wc(R, Y, X); \
|
|
403 } \
|
|
404 else \
|
|
405 { \
|
|
406 /* ediff == 0. */ \
|
|
407 if (!_FP_EXP_NORMAL(fs, wc, X)) \
|
|
408 { \
|
|
409 if (X##_e == 0) \
|
|
410 { \
|
|
411 /* X and Y are zero or denormalized. */ \
|
|
412 R##_e = 0; \
|
|
413 if (_FP_FRAC_ZEROP_##wc(X)) \
|
|
414 { \
|
|
415 if (!_FP_FRAC_ZEROP_##wc(Y)) \
|
|
416 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
417 _FP_FRAC_COPY_##wc(R, Y); \
|
|
418 goto add_done; \
|
|
419 } \
|
|
420 else if (_FP_FRAC_ZEROP_##wc(Y)) \
|
|
421 { \
|
|
422 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
423 _FP_FRAC_COPY_##wc(R, X); \
|
|
424 goto add_done; \
|
|
425 } \
|
|
426 else \
|
|
427 { \
|
|
428 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
429 _FP_FRAC_ADD_##wc(R, X, Y); \
|
|
430 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
|
|
431 { \
|
|
432 /* Normalized result. */ \
|
|
433 _FP_FRAC_HIGH_##fs(R) \
|
|
434 &= ~(_FP_W_TYPE)_FP_IMPLBIT_SH_##fs; \
|
|
435 R##_e = 1; \
|
|
436 } \
|
|
437 goto add_done; \
|
|
438 } \
|
|
439 } \
|
|
440 else \
|
|
441 { \
|
|
442 /* X and Y are NaN or Inf. */ \
|
|
443 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
|
|
444 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
|
|
445 R##_e = _FP_EXPMAX_##fs; \
|
|
446 if (_FP_FRAC_ZEROP_##wc(X)) \
|
|
447 _FP_FRAC_COPY_##wc(R, Y); \
|
|
448 else if (_FP_FRAC_ZEROP_##wc(Y)) \
|
|
449 _FP_FRAC_COPY_##wc(R, X); \
|
|
450 else \
|
|
451 _FP_CHOOSENAN_SEMIRAW(fs, wc, R, X, Y, OP); \
|
|
452 goto add_done; \
|
|
453 } \
|
|
454 } \
|
|
455 /* The exponents of X and Y, both normal, are equal. The \
|
|
456 implicit MSBs will always add to increase the \
|
|
457 exponent. */ \
|
|
458 _FP_FRAC_ADD_##wc(R, X, Y); \
|
|
459 R##_e = X##_e + 1; \
|
|
460 _FP_FRAC_SRS_##wc(R, 1, _FP_WFRACBITS_##fs); \
|
|
461 if (R##_e == _FP_EXPMAX_##fs) \
|
|
462 /* Overflow to infinity (depending on rounding mode). */ \
|
|
463 _FP_OVERFLOW_SEMIRAW(fs, wc, R); \
|
|
464 goto add_done; \
|
|
465 } \
|
|
466 add3: \
|
|
467 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
|
|
468 { \
|
|
469 /* Overflow. */ \
|
|
470 _FP_FRAC_HIGH_##fs(R) &= ~(_FP_W_TYPE)_FP_IMPLBIT_SH_##fs; \
|
|
471 R##_e++; \
|
|
472 _FP_FRAC_SRS_##wc(R, 1, _FP_WFRACBITS_##fs); \
|
|
473 if (R##_e == _FP_EXPMAX_##fs) \
|
|
474 /* Overflow to infinity (depending on rounding mode). */ \
|
|
475 _FP_OVERFLOW_SEMIRAW(fs, wc, R); \
|
|
476 } \
|
|
477 add_done: ; \
|
|
478 } \
|
|
479 else \
|
|
480 { \
|
|
481 /* Subtraction. */ \
|
|
482 int ediff = X##_e - Y##_e; \
|
|
483 if (ediff > 0) \
|
|
484 { \
|
|
485 R##_e = X##_e; \
|
|
486 R##_s = X##_s; \
|
|
487 if (Y##_e == 0) \
|
|
488 { \
|
|
489 /* Y is zero or denormalized. */ \
|
|
490 if (_FP_FRAC_ZEROP_##wc(Y)) \
|
|
491 { \
|
|
492 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
|
|
493 _FP_FRAC_COPY_##wc(R, X); \
|
|
494 goto sub_done; \
|
|
495 } \
|
|
496 else \
|
|
497 { \
|
|
498 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
499 ediff--; \
|
|
500 if (ediff == 0) \
|
|
501 { \
|
|
502 _FP_FRAC_SUB_##wc(R, X, Y); \
|
|
503 goto sub3; \
|
|
504 } \
|
|
505 if (X##_e == _FP_EXPMAX_##fs) \
|
|
506 { \
|
|
507 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
|
|
508 _FP_FRAC_COPY_##wc(R, X); \
|
|
509 goto sub_done; \
|
|
510 } \
|
|
511 goto sub1; \
|
|
512 } \
|
|
513 } \
|
|
514 else if (X##_e == _FP_EXPMAX_##fs) \
|
|
515 { \
|
|
516 /* X is NaN or Inf, Y is normal. */ \
|
|
517 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
|
|
518 _FP_FRAC_COPY_##wc(R, X); \
|
|
519 goto sub_done; \
|
|
520 } \
|
|
521 \
|
|
522 /* Insert implicit MSB of Y. */ \
|
|
523 _FP_FRAC_HIGH_##fs(Y) |= _FP_IMPLBIT_SH_##fs; \
|
|
524 \
|
|
525 sub1: \
|
|
526 /* Shift the mantissa of Y to the right EDIFF steps; \
|
|
527 remember to account later for the implicit MSB of X. */ \
|
|
528 if (ediff <= _FP_WFRACBITS_##fs) \
|
|
529 _FP_FRAC_SRS_##wc(Y, ediff, _FP_WFRACBITS_##fs); \
|
|
530 else if (!_FP_FRAC_ZEROP_##wc(Y)) \
|
|
531 _FP_FRAC_SET_##wc(Y, _FP_MINFRAC_##wc); \
|
|
532 _FP_FRAC_SUB_##wc(R, X, Y); \
|
|
533 } \
|
|
534 else if (ediff < 0) \
|
|
535 { \
|
|
536 ediff = -ediff; \
|
|
537 R##_e = Y##_e; \
|
|
538 R##_s = Y##_s; \
|
|
539 if (X##_e == 0) \
|
|
540 { \
|
|
541 /* X is zero or denormalized. */ \
|
|
542 if (_FP_FRAC_ZEROP_##wc(X)) \
|
|
543 { \
|
|
544 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
|
|
545 _FP_FRAC_COPY_##wc(R, Y); \
|
|
546 goto sub_done; \
|
|
547 } \
|
|
548 else \
|
|
549 { \
|
|
550 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
551 ediff--; \
|
|
552 if (ediff == 0) \
|
|
553 { \
|
|
554 _FP_FRAC_SUB_##wc(R, Y, X); \
|
|
555 goto sub3; \
|
|
556 } \
|
|
557 if (Y##_e == _FP_EXPMAX_##fs) \
|
|
558 { \
|
|
559 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
|
|
560 _FP_FRAC_COPY_##wc(R, Y); \
|
|
561 goto sub_done; \
|
|
562 } \
|
|
563 goto sub2; \
|
|
564 } \
|
|
565 } \
|
|
566 else if (Y##_e == _FP_EXPMAX_##fs) \
|
|
567 { \
|
|
568 /* Y is NaN or Inf, X is normal. */ \
|
|
569 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
|
|
570 _FP_FRAC_COPY_##wc(R, Y); \
|
|
571 goto sub_done; \
|
|
572 } \
|
|
573 \
|
|
574 /* Insert implicit MSB of X. */ \
|
|
575 _FP_FRAC_HIGH_##fs(X) |= _FP_IMPLBIT_SH_##fs; \
|
|
576 \
|
|
577 sub2: \
|
|
578 /* Shift the mantissa of X to the right EDIFF steps; \
|
|
579 remember to account later for the implicit MSB of Y. */ \
|
|
580 if (ediff <= _FP_WFRACBITS_##fs) \
|
|
581 _FP_FRAC_SRS_##wc(X, ediff, _FP_WFRACBITS_##fs); \
|
|
582 else if (!_FP_FRAC_ZEROP_##wc(X)) \
|
|
583 _FP_FRAC_SET_##wc(X, _FP_MINFRAC_##wc); \
|
|
584 _FP_FRAC_SUB_##wc(R, Y, X); \
|
|
585 } \
|
|
586 else \
|
|
587 { \
|
|
588 /* ediff == 0. */ \
|
|
589 if (!_FP_EXP_NORMAL(fs, wc, X)) \
|
|
590 { \
|
|
591 if (X##_e == 0) \
|
|
592 { \
|
|
593 /* X and Y are zero or denormalized. */ \
|
|
594 R##_e = 0; \
|
|
595 if (_FP_FRAC_ZEROP_##wc(X)) \
|
|
596 { \
|
|
597 _FP_FRAC_COPY_##wc(R, Y); \
|
|
598 if (_FP_FRAC_ZEROP_##wc(Y)) \
|
|
599 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
|
|
600 else \
|
|
601 { \
|
|
602 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
603 R##_s = Y##_s; \
|
|
604 } \
|
|
605 goto sub_done; \
|
|
606 } \
|
|
607 else if (_FP_FRAC_ZEROP_##wc(Y)) \
|
|
608 { \
|
|
609 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
610 _FP_FRAC_COPY_##wc(R, X); \
|
|
611 R##_s = X##_s; \
|
|
612 goto sub_done; \
|
|
613 } \
|
|
614 else \
|
|
615 { \
|
|
616 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
617 _FP_FRAC_SUB_##wc(R, X, Y); \
|
|
618 R##_s = X##_s; \
|
|
619 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
|
|
620 { \
|
|
621 /* |X| < |Y|, negate result. */ \
|
|
622 _FP_FRAC_SUB_##wc(R, Y, X); \
|
|
623 R##_s = Y##_s; \
|
|
624 } \
|
|
625 else if (_FP_FRAC_ZEROP_##wc(R)) \
|
|
626 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
|
|
627 goto sub_done; \
|
|
628 } \
|
|
629 } \
|
|
630 else \
|
|
631 { \
|
|
632 /* X and Y are NaN or Inf, of opposite signs. */ \
|
|
633 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
|
|
634 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
|
|
635 R##_e = _FP_EXPMAX_##fs; \
|
|
636 if (_FP_FRAC_ZEROP_##wc(X)) \
|
|
637 { \
|
|
638 if (_FP_FRAC_ZEROP_##wc(Y)) \
|
|
639 { \
|
|
640 /* Inf - Inf. */ \
|
|
641 R##_s = _FP_NANSIGN_##fs; \
|
|
642 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
|
|
643 _FP_FRAC_SLL_##wc(R, _FP_WORKBITS); \
|
|
644 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
645 } \
|
|
646 else \
|
|
647 { \
|
|
648 /* Inf - NaN. */ \
|
|
649 R##_s = Y##_s; \
|
|
650 _FP_FRAC_COPY_##wc(R, Y); \
|
|
651 } \
|
|
652 } \
|
|
653 else \
|
|
654 { \
|
|
655 if (_FP_FRAC_ZEROP_##wc(Y)) \
|
|
656 { \
|
|
657 /* NaN - Inf. */ \
|
|
658 R##_s = X##_s; \
|
|
659 _FP_FRAC_COPY_##wc(R, X); \
|
|
660 } \
|
|
661 else \
|
|
662 { \
|
|
663 /* NaN - NaN. */ \
|
|
664 _FP_CHOOSENAN_SEMIRAW(fs, wc, R, X, Y, OP); \
|
|
665 } \
|
|
666 } \
|
|
667 goto sub_done; \
|
|
668 } \
|
|
669 } \
|
|
670 /* The exponents of X and Y, both normal, are equal. The \
|
|
671 implicit MSBs cancel. */ \
|
|
672 R##_e = X##_e; \
|
|
673 _FP_FRAC_SUB_##wc(R, X, Y); \
|
|
674 R##_s = X##_s; \
|
|
675 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
|
|
676 { \
|
|
677 /* |X| < |Y|, negate result. */ \
|
|
678 _FP_FRAC_SUB_##wc(R, Y, X); \
|
|
679 R##_s = Y##_s; \
|
|
680 } \
|
|
681 else if (_FP_FRAC_ZEROP_##wc(R)) \
|
|
682 { \
|
|
683 R##_e = 0; \
|
|
684 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
|
|
685 goto sub_done; \
|
|
686 } \
|
|
687 goto norm; \
|
|
688 } \
|
|
689 sub3: \
|
|
690 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
|
|
691 { \
|
|
692 int diff; \
|
|
693 /* Carry into most significant bit of larger one of X and Y, \
|
|
694 canceling it; renormalize. */ \
|
|
695 _FP_FRAC_HIGH_##fs(R) &= _FP_IMPLBIT_SH_##fs - 1; \
|
|
696 norm: \
|
|
697 _FP_FRAC_CLZ_##wc(diff, R); \
|
|
698 diff -= _FP_WFRACXBITS_##fs; \
|
|
699 _FP_FRAC_SLL_##wc(R, diff); \
|
|
700 if (R##_e <= diff) \
|
|
701 { \
|
|
702 /* R is denormalized. */ \
|
|
703 diff = diff - R##_e + 1; \
|
|
704 _FP_FRAC_SRS_##wc(R, diff, _FP_WFRACBITS_##fs); \
|
|
705 R##_e = 0; \
|
|
706 } \
|
|
707 else \
|
|
708 { \
|
|
709 R##_e -= diff; \
|
|
710 _FP_FRAC_HIGH_##fs(R) &= ~(_FP_W_TYPE)_FP_IMPLBIT_SH_##fs; \
|
|
711 } \
|
|
712 } \
|
|
713 sub_done: ; \
|
|
714 } \
|
|
715 } while (0)
|
|
716
|
|
717 #define _FP_ADD(fs, wc, R, X, Y) _FP_ADD_INTERNAL(fs, wc, R, X, Y, '+')
|
|
718 #define _FP_SUB(fs, wc, R, X, Y) \
|
|
719 do { \
|
|
720 if (!(Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(Y))) Y##_s ^= 1; \
|
|
721 _FP_ADD_INTERNAL(fs, wc, R, X, Y, '-'); \
|
|
722 } while (0)
|
|
723
|
|
724
|
|
725 /*
|
|
726 * Main negation routine. FIXME -- when we care about setting exception
|
|
727 * bits reliably, this will not do. We should examine all of the fp classes.
|
|
728 */
|
|
729
|
|
730 #define _FP_NEG(fs, wc, R, X) \
|
|
731 do { \
|
|
732 _FP_FRAC_COPY_##wc(R, X); \
|
|
733 R##_c = X##_c; \
|
|
734 R##_e = X##_e; \
|
|
735 R##_s = 1 ^ X##_s; \
|
|
736 } while (0)
|
|
737
|
|
738
|
|
739 /*
|
|
740 * Main multiplication routine. The input values should be cooked.
|
|
741 */
|
|
742
|
|
743 #define _FP_MUL(fs, wc, R, X, Y) \
|
|
744 do { \
|
|
745 R##_s = X##_s ^ Y##_s; \
|
|
746 switch (_FP_CLS_COMBINE(X##_c, Y##_c)) \
|
|
747 { \
|
|
748 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_NORMAL): \
|
|
749 R##_c = FP_CLS_NORMAL; \
|
|
750 R##_e = X##_e + Y##_e + 1; \
|
|
751 \
|
|
752 _FP_MUL_MEAT_##fs(R,X,Y); \
|
|
753 \
|
|
754 if (_FP_FRAC_OVERP_##wc(fs, R)) \
|
|
755 _FP_FRAC_SRS_##wc(R, 1, _FP_WFRACBITS_##fs); \
|
|
756 else \
|
|
757 R##_e--; \
|
|
758 break; \
|
|
759 \
|
|
760 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_NAN): \
|
|
761 _FP_CHOOSENAN(fs, wc, R, X, Y, '*'); \
|
|
762 break; \
|
|
763 \
|
|
764 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_NORMAL): \
|
|
765 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_INF): \
|
|
766 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_ZERO): \
|
|
767 R##_s = X##_s; \
|
|
768 \
|
|
769 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_INF): \
|
|
770 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_NORMAL): \
|
|
771 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_NORMAL): \
|
|
772 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_ZERO): \
|
|
773 _FP_FRAC_COPY_##wc(R, X); \
|
|
774 R##_c = X##_c; \
|
|
775 break; \
|
|
776 \
|
|
777 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_NAN): \
|
|
778 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_NAN): \
|
|
779 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_NAN): \
|
|
780 R##_s = Y##_s; \
|
|
781 \
|
|
782 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_INF): \
|
|
783 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_ZERO): \
|
|
784 _FP_FRAC_COPY_##wc(R, Y); \
|
|
785 R##_c = Y##_c; \
|
|
786 break; \
|
|
787 \
|
|
788 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_ZERO): \
|
|
789 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_INF): \
|
|
790 R##_s = _FP_NANSIGN_##fs; \
|
|
791 R##_c = FP_CLS_NAN; \
|
|
792 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
|
|
793 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
794 break; \
|
|
795 \
|
|
796 default: \
|
|
797 abort(); \
|
|
798 } \
|
|
799 } while (0)
|
|
800
|
|
801
|
|
802 /*
|
|
803 * Main division routine. The input values should be cooked.
|
|
804 */
|
|
805
|
|
806 #define _FP_DIV(fs, wc, R, X, Y) \
|
|
807 do { \
|
|
808 R##_s = X##_s ^ Y##_s; \
|
|
809 switch (_FP_CLS_COMBINE(X##_c, Y##_c)) \
|
|
810 { \
|
|
811 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_NORMAL): \
|
|
812 R##_c = FP_CLS_NORMAL; \
|
|
813 R##_e = X##_e - Y##_e; \
|
|
814 \
|
|
815 _FP_DIV_MEAT_##fs(R,X,Y); \
|
|
816 break; \
|
|
817 \
|
|
818 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_NAN): \
|
|
819 _FP_CHOOSENAN(fs, wc, R, X, Y, '/'); \
|
|
820 break; \
|
|
821 \
|
|
822 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_NORMAL): \
|
|
823 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_INF): \
|
|
824 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_ZERO): \
|
|
825 R##_s = X##_s; \
|
|
826 _FP_FRAC_COPY_##wc(R, X); \
|
|
827 R##_c = X##_c; \
|
|
828 break; \
|
|
829 \
|
|
830 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_NAN): \
|
|
831 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_NAN): \
|
|
832 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_NAN): \
|
|
833 R##_s = Y##_s; \
|
|
834 _FP_FRAC_COPY_##wc(R, Y); \
|
|
835 R##_c = Y##_c; \
|
|
836 break; \
|
|
837 \
|
|
838 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_INF): \
|
|
839 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_INF): \
|
|
840 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_NORMAL): \
|
|
841 R##_c = FP_CLS_ZERO; \
|
|
842 break; \
|
|
843 \
|
|
844 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_ZERO): \
|
|
845 FP_SET_EXCEPTION(FP_EX_DIVZERO); \
|
|
846 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_ZERO): \
|
|
847 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_NORMAL): \
|
|
848 R##_c = FP_CLS_INF; \
|
|
849 break; \
|
|
850 \
|
|
851 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_INF): \
|
|
852 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_ZERO): \
|
|
853 R##_s = _FP_NANSIGN_##fs; \
|
|
854 R##_c = FP_CLS_NAN; \
|
|
855 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
|
|
856 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
857 break; \
|
|
858 \
|
|
859 default: \
|
|
860 abort(); \
|
|
861 } \
|
|
862 } while (0)
|
|
863
|
|
864
|
|
865 /*
|
|
866 * Main differential comparison routine. The inputs should be raw not
|
|
867 * cooked. The return is -1,0,1 for normal values, 2 otherwise.
|
|
868 */
|
|
869
|
|
870 #define _FP_CMP(fs, wc, ret, X, Y, un) \
|
|
871 do { \
|
|
872 /* NANs are unordered */ \
|
|
873 if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(X)) \
|
|
874 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(Y))) \
|
|
875 { \
|
|
876 ret = un; \
|
|
877 } \
|
|
878 else \
|
|
879 { \
|
|
880 int __is_zero_x; \
|
|
881 int __is_zero_y; \
|
|
882 \
|
|
883 __is_zero_x = (!X##_e && _FP_FRAC_ZEROP_##wc(X)) ? 1 : 0; \
|
|
884 __is_zero_y = (!Y##_e && _FP_FRAC_ZEROP_##wc(Y)) ? 1 : 0; \
|
|
885 \
|
|
886 if (__is_zero_x && __is_zero_y) \
|
|
887 ret = 0; \
|
|
888 else if (__is_zero_x) \
|
|
889 ret = Y##_s ? 1 : -1; \
|
|
890 else if (__is_zero_y) \
|
|
891 ret = X##_s ? -1 : 1; \
|
|
892 else if (X##_s != Y##_s) \
|
|
893 ret = X##_s ? -1 : 1; \
|
|
894 else if (X##_e > Y##_e) \
|
|
895 ret = X##_s ? -1 : 1; \
|
|
896 else if (X##_e < Y##_e) \
|
|
897 ret = X##_s ? 1 : -1; \
|
|
898 else if (_FP_FRAC_GT_##wc(X, Y)) \
|
|
899 ret = X##_s ? -1 : 1; \
|
|
900 else if (_FP_FRAC_GT_##wc(Y, X)) \
|
|
901 ret = X##_s ? 1 : -1; \
|
|
902 else \
|
|
903 ret = 0; \
|
|
904 } \
|
|
905 } while (0)
|
|
906
|
|
907
|
|
908 /* Simplification for strict equality. */
|
|
909
|
|
910 #define _FP_CMP_EQ(fs, wc, ret, X, Y) \
|
|
911 do { \
|
|
912 /* NANs are unordered */ \
|
|
913 if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(X)) \
|
|
914 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(Y))) \
|
|
915 { \
|
|
916 ret = 1; \
|
|
917 } \
|
|
918 else \
|
|
919 { \
|
|
920 ret = !(X##_e == Y##_e \
|
|
921 && _FP_FRAC_EQ_##wc(X, Y) \
|
|
922 && (X##_s == Y##_s || (!X##_e && _FP_FRAC_ZEROP_##wc(X)))); \
|
|
923 } \
|
|
924 } while (0)
|
|
925
|
|
926 /* Version to test unordered. */
|
|
927
|
|
928 #define _FP_CMP_UNORD(fs, wc, ret, X, Y) \
|
|
929 do { \
|
|
930 ret = ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(X)) \
|
|
931 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(Y))); \
|
|
932 } while (0)
|
|
933
|
|
934 /*
|
|
935 * Main square root routine. The input value should be cooked.
|
|
936 */
|
|
937
|
|
938 #define _FP_SQRT(fs, wc, R, X) \
|
|
939 do { \
|
|
940 _FP_FRAC_DECL_##wc(T); _FP_FRAC_DECL_##wc(S); \
|
|
941 _FP_W_TYPE q; \
|
|
942 switch (X##_c) \
|
|
943 { \
|
|
944 case FP_CLS_NAN: \
|
|
945 _FP_FRAC_COPY_##wc(R, X); \
|
|
946 R##_s = X##_s; \
|
|
947 R##_c = FP_CLS_NAN; \
|
|
948 break; \
|
|
949 case FP_CLS_INF: \
|
|
950 if (X##_s) \
|
|
951 { \
|
|
952 R##_s = _FP_NANSIGN_##fs; \
|
|
953 R##_c = FP_CLS_NAN; /* NAN */ \
|
|
954 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
|
|
955 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
956 } \
|
|
957 else \
|
|
958 { \
|
|
959 R##_s = 0; \
|
|
960 R##_c = FP_CLS_INF; /* sqrt(+inf) = +inf */ \
|
|
961 } \
|
|
962 break; \
|
|
963 case FP_CLS_ZERO: \
|
|
964 R##_s = X##_s; \
|
|
965 R##_c = FP_CLS_ZERO; /* sqrt(+-0) = +-0 */ \
|
|
966 break; \
|
|
967 case FP_CLS_NORMAL: \
|
|
968 R##_s = 0; \
|
|
969 if (X##_s) \
|
|
970 { \
|
|
971 R##_c = FP_CLS_NAN; /* sNAN */ \
|
|
972 R##_s = _FP_NANSIGN_##fs; \
|
|
973 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
|
|
974 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
975 break; \
|
|
976 } \
|
|
977 R##_c = FP_CLS_NORMAL; \
|
|
978 if (X##_e & 1) \
|
|
979 _FP_FRAC_SLL_##wc(X, 1); \
|
|
980 R##_e = X##_e >> 1; \
|
|
981 _FP_FRAC_SET_##wc(S, _FP_ZEROFRAC_##wc); \
|
|
982 _FP_FRAC_SET_##wc(R, _FP_ZEROFRAC_##wc); \
|
|
983 q = _FP_OVERFLOW_##fs >> 1; \
|
|
984 _FP_SQRT_MEAT_##wc(R, S, T, X, q); \
|
|
985 } \
|
|
986 } while (0)
|
|
987
|
|
988 /*
|
|
989 * Convert from FP to integer. Input is raw.
|
|
990 */
|
|
991
|
|
992 /* RSIGNED can have following values:
|
|
993 * 0: the number is required to be 0..(2^rsize)-1, if not, NV is set plus
|
|
994 * the result is either 0 or (2^rsize)-1 depending on the sign in such
|
|
995 * case.
|
|
996 * 1: the number is required to be -(2^(rsize-1))..(2^(rsize-1))-1, if not,
|
|
997 * NV is set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
|
|
998 * depending on the sign in such case.
|
|
999 * -1: the number is required to be -(2^(rsize-1))..(2^rsize)-1, if not, NV is
|
|
1000 * set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
|
|
1001 * depending on the sign in such case.
|
|
1002 */
|
|
1003 #define _FP_TO_INT(fs, wc, r, X, rsize, rsigned) \
|
|
1004 do { \
|
|
1005 if (X##_e < _FP_EXPBIAS_##fs) \
|
|
1006 { \
|
|
1007 r = 0; \
|
|
1008 if (X##_e == 0) \
|
|
1009 { \
|
|
1010 if (!_FP_FRAC_ZEROP_##wc(X)) \
|
|
1011 { \
|
|
1012 FP_SET_EXCEPTION(FP_EX_INEXACT); \
|
|
1013 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
1014 } \
|
|
1015 } \
|
|
1016 else \
|
|
1017 FP_SET_EXCEPTION(FP_EX_INEXACT); \
|
|
1018 } \
|
|
1019 else if (X##_e >= _FP_EXPBIAS_##fs + rsize - (rsigned > 0 || X##_s) \
|
|
1020 || (!rsigned && X##_s)) \
|
|
1021 { \
|
|
1022 /* Overflow or converting to the most negative integer. */ \
|
|
1023 if (rsigned) \
|
|
1024 { \
|
|
1025 r = 1; \
|
|
1026 r <<= rsize - 1; \
|
|
1027 r -= 1 - X##_s; \
|
|
1028 } else { \
|
|
1029 r = 0; \
|
|
1030 if (X##_s) \
|
|
1031 r = ~r; \
|
|
1032 } \
|
|
1033 \
|
|
1034 if (rsigned && X##_s && X##_e == _FP_EXPBIAS_##fs + rsize - 1) \
|
|
1035 { \
|
|
1036 /* Possibly converting to most negative integer; check the \
|
|
1037 mantissa. */ \
|
|
1038 int inexact = 0; \
|
|
1039 (void)((_FP_FRACBITS_##fs > rsize) \
|
|
1040 ? ({ _FP_FRAC_SRST_##wc(X, inexact, \
|
|
1041 _FP_FRACBITS_##fs - rsize, \
|
|
1042 _FP_FRACBITS_##fs); 0; }) \
|
|
1043 : 0); \
|
|
1044 if (!_FP_FRAC_ZEROP_##wc(X)) \
|
|
1045 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
1046 else if (inexact) \
|
|
1047 FP_SET_EXCEPTION(FP_EX_INEXACT); \
|
|
1048 } \
|
|
1049 else \
|
|
1050 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
1051 } \
|
|
1052 else \
|
|
1053 { \
|
|
1054 _FP_FRAC_HIGH_RAW_##fs(X) |= _FP_IMPLBIT_##fs; \
|
|
1055 if (X##_e >= _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs - 1) \
|
|
1056 { \
|
|
1057 _FP_FRAC_ASSEMBLE_##wc(r, X, rsize); \
|
|
1058 r <<= X##_e - _FP_EXPBIAS_##fs - _FP_FRACBITS_##fs + 1; \
|
|
1059 } \
|
|
1060 else \
|
|
1061 { \
|
|
1062 int inexact; \
|
|
1063 _FP_FRAC_SRST_##wc(X, inexact, \
|
|
1064 (_FP_FRACBITS_##fs + _FP_EXPBIAS_##fs - 1 \
|
|
1065 - X##_e), \
|
|
1066 _FP_FRACBITS_##fs); \
|
|
1067 if (inexact) \
|
|
1068 FP_SET_EXCEPTION(FP_EX_INEXACT); \
|
|
1069 _FP_FRAC_ASSEMBLE_##wc(r, X, rsize); \
|
|
1070 } \
|
|
1071 if (rsigned && X##_s) \
|
|
1072 r = -r; \
|
|
1073 } \
|
|
1074 } while (0)
|
|
1075
|
|
1076 /* Convert integer to fp. Output is raw. RTYPE is unsigned even if
|
|
1077 input is signed. */
|
|
1078 #define _FP_FROM_INT(fs, wc, X, r, rsize, rtype) \
|
|
1079 do { \
|
|
1080 if (r) \
|
|
1081 { \
|
|
1082 rtype ur_; \
|
|
1083 \
|
|
1084 if ((X##_s = (r < 0))) \
|
|
1085 r = -(rtype)r; \
|
|
1086 \
|
|
1087 ur_ = (rtype) r; \
|
|
1088 (void)((rsize <= _FP_W_TYPE_SIZE) \
|
|
1089 ? ({ \
|
|
1090 int lz_; \
|
|
1091 __FP_CLZ(lz_, (_FP_W_TYPE)ur_); \
|
|
1092 X##_e = _FP_EXPBIAS_##fs + _FP_W_TYPE_SIZE - 1 - lz_; \
|
|
1093 }) \
|
|
1094 : ((rsize <= 2 * _FP_W_TYPE_SIZE) \
|
|
1095 ? ({ \
|
|
1096 int lz_; \
|
|
1097 __FP_CLZ_2(lz_, (_FP_W_TYPE)(ur_ >> _FP_W_TYPE_SIZE), \
|
|
1098 (_FP_W_TYPE)ur_); \
|
|
1099 X##_e = (_FP_EXPBIAS_##fs + 2 * _FP_W_TYPE_SIZE - 1 \
|
|
1100 - lz_); \
|
|
1101 }) \
|
|
1102 : (abort(), 0))); \
|
|
1103 \
|
|
1104 if (rsize - 1 + _FP_EXPBIAS_##fs >= _FP_EXPMAX_##fs \
|
|
1105 && X##_e >= _FP_EXPMAX_##fs) \
|
|
1106 { \
|
|
1107 /* Exponent too big; overflow to infinity. (May also \
|
|
1108 happen after rounding below.) */ \
|
|
1109 _FP_OVERFLOW_SEMIRAW(fs, wc, X); \
|
|
1110 goto pack_semiraw; \
|
|
1111 } \
|
|
1112 \
|
|
1113 if (rsize <= _FP_FRACBITS_##fs \
|
|
1114 || X##_e < _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs) \
|
|
1115 { \
|
|
1116 /* Exactly representable; shift left. */ \
|
|
1117 _FP_FRAC_DISASSEMBLE_##wc(X, ur_, rsize); \
|
|
1118 _FP_FRAC_SLL_##wc(X, (_FP_EXPBIAS_##fs \
|
|
1119 + _FP_FRACBITS_##fs - 1 - X##_e)); \
|
|
1120 } \
|
|
1121 else \
|
|
1122 { \
|
|
1123 /* More bits in integer than in floating type; need to \
|
|
1124 round. */ \
|
|
1125 if (_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 < X##_e) \
|
|
1126 ur_ = ((ur_ >> (X##_e - _FP_EXPBIAS_##fs \
|
|
1127 - _FP_WFRACBITS_##fs + 1)) \
|
|
1128 | ((ur_ << (rsize - (X##_e - _FP_EXPBIAS_##fs \
|
|
1129 - _FP_WFRACBITS_##fs + 1))) \
|
|
1130 != 0)); \
|
|
1131 _FP_FRAC_DISASSEMBLE_##wc(X, ur_, rsize); \
|
|
1132 if ((_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 - X##_e) > 0) \
|
|
1133 _FP_FRAC_SLL_##wc(X, (_FP_EXPBIAS_##fs \
|
|
1134 + _FP_WFRACBITS_##fs - 1 - X##_e)); \
|
|
1135 _FP_FRAC_HIGH_##fs(X) &= ~(_FP_W_TYPE)_FP_IMPLBIT_SH_##fs; \
|
|
1136 pack_semiraw: \
|
|
1137 _FP_PACK_SEMIRAW(fs, wc, X); \
|
|
1138 } \
|
|
1139 } \
|
|
1140 else \
|
|
1141 { \
|
|
1142 X##_s = 0; \
|
|
1143 X##_e = 0; \
|
|
1144 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
|
|
1145 } \
|
|
1146 } while (0)
|
|
1147
|
|
1148
|
|
1149 /* Extend from a narrower floating-point format to a wider one. Input
|
|
1150 and output are raw. */
|
|
1151 #define FP_EXTEND(dfs,sfs,dwc,swc,D,S) \
|
|
1152 do { \
|
|
1153 if (_FP_FRACBITS_##dfs < _FP_FRACBITS_##sfs \
|
|
1154 || (_FP_EXPMAX_##dfs - _FP_EXPBIAS_##dfs \
|
|
1155 < _FP_EXPMAX_##sfs - _FP_EXPBIAS_##sfs) \
|
|
1156 || (_FP_EXPBIAS_##dfs < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1 \
|
|
1157 && _FP_EXPBIAS_##dfs != _FP_EXPBIAS_##sfs)) \
|
|
1158 abort(); \
|
|
1159 D##_s = S##_s; \
|
|
1160 _FP_FRAC_COPY_##dwc##_##swc(D, S); \
|
|
1161 if (_FP_EXP_NORMAL(sfs, swc, S)) \
|
|
1162 { \
|
|
1163 D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs; \
|
|
1164 _FP_FRAC_SLL_##dwc(D, (_FP_FRACBITS_##dfs - _FP_FRACBITS_##sfs)); \
|
|
1165 } \
|
|
1166 else \
|
|
1167 { \
|
|
1168 if (S##_e == 0) \
|
|
1169 { \
|
|
1170 if (_FP_FRAC_ZEROP_##swc(S)) \
|
|
1171 D##_e = 0; \
|
|
1172 else if (_FP_EXPBIAS_##dfs \
|
|
1173 < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1) \
|
|
1174 { \
|
|
1175 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
1176 _FP_FRAC_SLL_##dwc(D, (_FP_FRACBITS_##dfs \
|
|
1177 - _FP_FRACBITS_##sfs)); \
|
|
1178 D##_e = 0; \
|
|
1179 } \
|
|
1180 else \
|
|
1181 { \
|
|
1182 int _lz; \
|
|
1183 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
1184 _FP_FRAC_CLZ_##swc(_lz, S); \
|
|
1185 _FP_FRAC_SLL_##dwc(D, \
|
|
1186 _lz + _FP_FRACBITS_##dfs \
|
|
1187 - _FP_FRACTBITS_##sfs); \
|
|
1188 D##_e = (_FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs + 1 \
|
|
1189 + _FP_FRACXBITS_##sfs - _lz); \
|
|
1190 } \
|
|
1191 } \
|
|
1192 else \
|
|
1193 { \
|
|
1194 D##_e = _FP_EXPMAX_##dfs; \
|
|
1195 if (!_FP_FRAC_ZEROP_##swc(S)) \
|
|
1196 { \
|
|
1197 if (!(_FP_FRAC_HIGH_RAW_##sfs(S) & _FP_QNANBIT_##sfs)) \
|
|
1198 FP_SET_EXCEPTION(FP_EX_INVALID); \
|
|
1199 _FP_FRAC_SLL_##dwc(D, (_FP_FRACBITS_##dfs \
|
|
1200 - _FP_FRACBITS_##sfs)); \
|
|
1201 } \
|
|
1202 } \
|
|
1203 } \
|
|
1204 } while (0)
|
|
1205
|
|
1206 /* Truncate from a wider floating-point format to a narrower one.
|
|
1207 Input and output are semi-raw. */
|
|
1208 #define FP_TRUNC(dfs,sfs,dwc,swc,D,S) \
|
|
1209 do { \
|
|
1210 if (_FP_FRACBITS_##sfs < _FP_FRACBITS_##dfs \
|
|
1211 || (_FP_EXPBIAS_##sfs < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1 \
|
|
1212 && _FP_EXPBIAS_##sfs != _FP_EXPBIAS_##dfs)) \
|
|
1213 abort(); \
|
|
1214 D##_s = S##_s; \
|
|
1215 if (_FP_EXP_NORMAL(sfs, swc, S)) \
|
|
1216 { \
|
|
1217 D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs; \
|
|
1218 if (D##_e >= _FP_EXPMAX_##dfs) \
|
|
1219 _FP_OVERFLOW_SEMIRAW(dfs, dwc, D); \
|
|
1220 else \
|
|
1221 { \
|
|
1222 if (D##_e <= 0) \
|
|
1223 { \
|
|
1224 if (D##_e < 1 - _FP_FRACBITS_##dfs) \
|
|
1225 { \
|
|
1226 _FP_FRAC_SET_##swc(S, _FP_ZEROFRAC_##swc); \
|
|
1227 _FP_FRAC_LOW_##swc(S) |= 1; \
|
|
1228 } \
|
|
1229 else \
|
|
1230 { \
|
|
1231 _FP_FRAC_HIGH_##sfs(S) |= _FP_IMPLBIT_SH_##sfs; \
|
|
1232 _FP_FRAC_SRS_##swc(S, (_FP_WFRACBITS_##sfs \
|
|
1233 - _FP_WFRACBITS_##dfs + 1 - D##_e), \
|
|
1234 _FP_WFRACBITS_##sfs); \
|
|
1235 } \
|
|
1236 D##_e = 0; \
|
|
1237 } \
|
|
1238 else \
|
|
1239 _FP_FRAC_SRS_##swc(S, (_FP_WFRACBITS_##sfs \
|
|
1240 - _FP_WFRACBITS_##dfs), \
|
|
1241 _FP_WFRACBITS_##sfs); \
|
|
1242 _FP_FRAC_COPY_##dwc##_##swc(D, S); \
|
|
1243 } \
|
|
1244 } \
|
|
1245 else \
|
|
1246 { \
|
|
1247 if (S##_e == 0) \
|
|
1248 { \
|
|
1249 D##_e = 0; \
|
|
1250 if (_FP_FRAC_ZEROP_##swc(S)) \
|
|
1251 _FP_FRAC_SET_##dwc(D, _FP_ZEROFRAC_##dwc); \
|
|
1252 else \
|
|
1253 { \
|
|
1254 FP_SET_EXCEPTION(FP_EX_DENORM); \
|
|
1255 if (_FP_EXPBIAS_##sfs \
|
|
1256 < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1) \
|
|
1257 { \
|
|
1258 _FP_FRAC_SRS_##swc(S, (_FP_WFRACBITS_##sfs \
|
|
1259 - _FP_WFRACBITS_##dfs), \
|
|
1260 _FP_WFRACBITS_##sfs); \
|
|
1261 _FP_FRAC_COPY_##dwc##_##swc(D, S); \
|
|
1262 } \
|
|
1263 else \
|
|
1264 { \
|
|
1265 _FP_FRAC_SET_##dwc(D, _FP_ZEROFRAC_##dwc); \
|
|
1266 _FP_FRAC_LOW_##dwc(D) |= 1; \
|
|
1267 } \
|
|
1268 } \
|
|
1269 } \
|
|
1270 else \
|
|
1271 { \
|
|
1272 D##_e = _FP_EXPMAX_##dfs; \
|
|
1273 if (_FP_FRAC_ZEROP_##swc(S)) \
|
|
1274 _FP_FRAC_SET_##dwc(D, _FP_ZEROFRAC_##dwc); \
|
|
1275 else \
|
|
1276 { \
|
|
1277 _FP_CHECK_SIGNAN_SEMIRAW(sfs, swc, S); \
|
|
1278 _FP_FRAC_SRL_##swc(S, (_FP_WFRACBITS_##sfs \
|
|
1279 - _FP_WFRACBITS_##dfs)); \
|
|
1280 _FP_FRAC_COPY_##dwc##_##swc(D, S); \
|
|
1281 /* Semi-raw NaN must have all workbits cleared. */ \
|
|
1282 _FP_FRAC_LOW_##dwc(D) \
|
|
1283 &= ~(_FP_W_TYPE) ((1 << _FP_WORKBITS) - 1); \
|
|
1284 _FP_FRAC_HIGH_##dfs(D) |= _FP_QNANBIT_SH_##dfs; \
|
|
1285 } \
|
|
1286 } \
|
|
1287 } \
|
|
1288 } while (0)
|
|
1289
|
|
1290 /*
|
|
1291 * Helper primitives.
|
|
1292 */
|
|
1293
|
|
1294 /* Count leading zeros in a word. */
|
|
1295
|
|
1296 #ifndef __FP_CLZ
|
|
1297 /* GCC 3.4 and later provide the builtins for us. */
|
|
1298 #define __FP_CLZ(r, x) \
|
|
1299 do { \
|
|
1300 if (sizeof (_FP_W_TYPE) == sizeof (unsigned int)) \
|
|
1301 r = __builtin_clz (x); \
|
|
1302 else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long)) \
|
|
1303 r = __builtin_clzl (x); \
|
|
1304 else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long long)) \
|
|
1305 r = __builtin_clzll (x); \
|
|
1306 else \
|
|
1307 abort (); \
|
|
1308 } while (0)
|
|
1309 #endif /* ndef __FP_CLZ */
|
|
1310
|
|
1311 #define _FP_DIV_HELP_imm(q, r, n, d) \
|
|
1312 do { \
|
|
1313 q = n / d, r = n % d; \
|
|
1314 } while (0)
|
|
1315
|
|
1316
|
|
1317 /* A restoring bit-by-bit division primitive. */
|
|
1318
|
|
1319 #define _FP_DIV_MEAT_N_loop(fs, wc, R, X, Y) \
|
|
1320 do { \
|
|
1321 int count = _FP_WFRACBITS_##fs; \
|
|
1322 _FP_FRAC_DECL_##wc (u); \
|
|
1323 _FP_FRAC_DECL_##wc (v); \
|
|
1324 _FP_FRAC_COPY_##wc (u, X); \
|
|
1325 _FP_FRAC_COPY_##wc (v, Y); \
|
|
1326 _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc); \
|
|
1327 /* Normalize U and V. */ \
|
|
1328 _FP_FRAC_SLL_##wc (u, _FP_WFRACXBITS_##fs); \
|
|
1329 _FP_FRAC_SLL_##wc (v, _FP_WFRACXBITS_##fs); \
|
|
1330 /* First round. Since the operands are normalized, either the \
|
|
1331 first or second bit will be set in the fraction. Produce a \
|
|
1332 normalized result by checking which and adjusting the loop \
|
|
1333 count and exponent accordingly. */ \
|
|
1334 if (_FP_FRAC_GE_1 (u, v)) \
|
|
1335 { \
|
|
1336 _FP_FRAC_SUB_##wc (u, u, v); \
|
|
1337 _FP_FRAC_LOW_##wc (R) |= 1; \
|
|
1338 count--; \
|
|
1339 } \
|
|
1340 else \
|
|
1341 R##_e--; \
|
|
1342 /* Subsequent rounds. */ \
|
|
1343 do { \
|
|
1344 int msb = (_FP_WS_TYPE) _FP_FRAC_HIGH_##wc (u) < 0; \
|
|
1345 _FP_FRAC_SLL_##wc (u, 1); \
|
|
1346 _FP_FRAC_SLL_##wc (R, 1); \
|
|
1347 if (msb || _FP_FRAC_GE_1 (u, v)) \
|
|
1348 { \
|
|
1349 _FP_FRAC_SUB_##wc (u, u, v); \
|
|
1350 _FP_FRAC_LOW_##wc (R) |= 1; \
|
|
1351 } \
|
|
1352 } while (--count > 0); \
|
|
1353 /* If there's anything left in U, the result is inexact. */ \
|
|
1354 _FP_FRAC_LOW_##wc (R) |= !_FP_FRAC_ZEROP_##wc (u); \
|
|
1355 } while (0)
|
|
1356
|
|
1357 #define _FP_DIV_MEAT_1_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 1, R, X, Y)
|
|
1358 #define _FP_DIV_MEAT_2_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 2, R, X, Y)
|
|
1359 #define _FP_DIV_MEAT_4_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 4, R, X, Y)
|