Properly handle unavailable elements in LC_MONETARY category
[glibc.git] / soft-fp / op-common.h
blob67cdc33b4cf9eef360780b8506a55af349d7d6b7
1 /* Software floating-point emulation. Common operations.
2 Copyright (C) 1997-2013 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).
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.
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.)
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.
28 You should have received a copy of the GNU Lesser General Public
29 License along with the GNU C Library; if not, see
30 <http://www.gnu.org/licenses/>. */
32 #define _FP_DECL(wc, X) \
33 _FP_I_TYPE X##_c __attribute__ ((unused)); \
34 _FP_I_TYPE X##_s __attribute__ ((unused)); \
35 _FP_I_TYPE X##_e; \
36 _FP_FRAC_DECL_##wc (X)
38 /* Test whether the qNaN bit denotes a signaling NaN. */
39 #define _FP_FRAC_SNANP(fs, X) \
40 ((_FP_QNANNEGATEDP) \
41 ? (_FP_FRAC_HIGH_RAW_##fs (X) & _FP_QNANBIT_##fs) \
42 : !(_FP_FRAC_HIGH_RAW_##fs (X) & _FP_QNANBIT_##fs))
43 #define _FP_FRAC_SNANP_SEMIRAW(fs, X) \
44 ((_FP_QNANNEGATEDP) \
45 ? (_FP_FRAC_HIGH_##fs (X) & _FP_QNANBIT_SH_##fs) \
46 : !(_FP_FRAC_HIGH_##fs (X) & _FP_QNANBIT_SH_##fs))
49 * Finish truly unpacking a native fp value by classifying the kind
50 * of fp value and normalizing both the exponent and the fraction.
53 #define _FP_UNPACK_CANONICAL(fs, wc, X) \
54 do \
55 { \
56 switch (X##_e) \
57 { \
58 default: \
59 _FP_FRAC_HIGH_RAW_##fs (X) |= _FP_IMPLBIT_##fs; \
60 _FP_FRAC_SLL_##wc (X, _FP_WORKBITS); \
61 X##_e -= _FP_EXPBIAS_##fs; \
62 X##_c = FP_CLS_NORMAL; \
63 break; \
65 case 0: \
66 if (_FP_FRAC_ZEROP_##wc (X)) \
67 X##_c = FP_CLS_ZERO; \
68 else \
69 { \
70 /* a denormalized number */ \
71 _FP_I_TYPE _shift; \
72 _FP_FRAC_CLZ_##wc (_shift, X); \
73 _shift -= _FP_FRACXBITS_##fs; \
74 _FP_FRAC_SLL_##wc (X, (_shift+_FP_WORKBITS)); \
75 X##_e -= _FP_EXPBIAS_##fs - 1 + _shift; \
76 X##_c = FP_CLS_NORMAL; \
77 FP_SET_EXCEPTION (FP_EX_DENORM); \
78 } \
79 break; \
81 case _FP_EXPMAX_##fs: \
82 if (_FP_FRAC_ZEROP_##wc (X)) \
83 X##_c = FP_CLS_INF; \
84 else \
85 { \
86 X##_c = FP_CLS_NAN; \
87 /* Check for signaling NaN */ \
88 if (_FP_FRAC_SNANP (fs, X)) \
89 FP_SET_EXCEPTION (FP_EX_INVALID); \
90 } \
91 break; \
92 } \
93 } \
94 while (0)
96 /* Finish unpacking an fp value in semi-raw mode: the mantissa is
97 shifted by _FP_WORKBITS but the implicit MSB is not inserted and
98 other classification is not done. */
99 #define _FP_UNPACK_SEMIRAW(fs, wc, X) _FP_FRAC_SLL_##wc (X, _FP_WORKBITS)
101 /* A semi-raw value has overflowed to infinity. Adjust the mantissa
102 and exponent appropriately. */
103 #define _FP_OVERFLOW_SEMIRAW(fs, wc, X) \
104 do \
106 if (FP_ROUNDMODE == FP_RND_NEAREST \
107 || (FP_ROUNDMODE == FP_RND_PINF && !X##_s) \
108 || (FP_ROUNDMODE == FP_RND_MINF && X##_s)) \
110 X##_e = _FP_EXPMAX_##fs; \
111 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
113 else \
115 X##_e = _FP_EXPMAX_##fs - 1; \
116 _FP_FRAC_SET_##wc (X, _FP_MAXFRAC_##wc); \
118 FP_SET_EXCEPTION (FP_EX_INEXACT); \
119 FP_SET_EXCEPTION (FP_EX_OVERFLOW); \
121 while (0)
123 /* Check for a semi-raw value being a signaling NaN and raise the
124 invalid exception if so. */
125 #define _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X) \
126 do \
128 if (X##_e == _FP_EXPMAX_##fs \
129 && !_FP_FRAC_ZEROP_##wc (X) \
130 && _FP_FRAC_SNANP_SEMIRAW (fs, X)) \
131 FP_SET_EXCEPTION (FP_EX_INVALID); \
133 while (0)
135 /* Choose a NaN result from an operation on two semi-raw NaN
136 values. */
137 #define _FP_CHOOSENAN_SEMIRAW(fs, wc, R, X, Y, OP) \
138 do \
140 /* _FP_CHOOSENAN expects raw values, so shift as required. */ \
141 _FP_FRAC_SRL_##wc (X, _FP_WORKBITS); \
142 _FP_FRAC_SRL_##wc (Y, _FP_WORKBITS); \
143 _FP_CHOOSENAN (fs, wc, R, X, Y, OP); \
144 _FP_FRAC_SLL_##wc (R, _FP_WORKBITS); \
146 while (0)
148 /* Make the fractional part a quiet NaN, preserving the payload
149 if possible, otherwise make it the canonical quiet NaN and set
150 the sign bit accordingly. */
151 #define _FP_SETQNAN(fs, wc, X) \
152 do \
154 if (_FP_QNANNEGATEDP) \
156 _FP_FRAC_HIGH_RAW_##fs (X) &= _FP_QNANBIT_##fs - 1; \
157 if (_FP_FRAC_ZEROP_##wc (X)) \
159 X##_s = _FP_NANSIGN_##fs; \
160 _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs); \
163 else \
164 _FP_FRAC_HIGH_RAW_##fs (X) |= _FP_QNANBIT_##fs; \
166 while (0)
167 #define _FP_SETQNAN_SEMIRAW(fs, wc, X) \
168 do \
170 if (_FP_QNANNEGATEDP) \
172 _FP_FRAC_HIGH_##fs (X) &= _FP_QNANBIT_SH_##fs - 1; \
173 if (_FP_FRAC_ZEROP_##wc (X)) \
175 X##_s = _FP_NANSIGN_##fs; \
176 _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs); \
177 _FP_FRAC_SLL_##wc (X, _FP_WORKBITS); \
180 else \
181 _FP_FRAC_HIGH_##fs (X) |= _FP_QNANBIT_SH_##fs; \
183 while (0)
185 /* Test whether a biased exponent is normal (not zero or maximum). */
186 #define _FP_EXP_NORMAL(fs, wc, X) (((X##_e + 1) & _FP_EXPMAX_##fs) > 1)
188 /* Prepare to pack an fp value in semi-raw mode: the mantissa is
189 rounded and shifted right, with the rounding possibly increasing
190 the exponent (including changing a finite value to infinity). */
191 #define _FP_PACK_SEMIRAW(fs, wc, X) \
192 do \
194 _FP_ROUND (wc, X); \
195 if (X##_e == 0 && !_FP_FRAC_ZEROP_##wc (X)) \
197 if ((FP_CUR_EXCEPTIONS & FP_EX_INEXACT) \
198 || (FP_TRAPPING_EXCEPTIONS & FP_EX_UNDERFLOW)) \
199 FP_SET_EXCEPTION (FP_EX_UNDERFLOW); \
201 if (_FP_FRAC_HIGH_##fs (X) \
202 & (_FP_OVERFLOW_##fs >> 1)) \
204 _FP_FRAC_HIGH_##fs (X) &= ~(_FP_OVERFLOW_##fs >> 1); \
205 X##_e++; \
206 if (X##_e == _FP_EXPMAX_##fs) \
207 _FP_OVERFLOW_SEMIRAW (fs, wc, X); \
209 _FP_FRAC_SRL_##wc (X, _FP_WORKBITS); \
210 if (X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
212 if (!_FP_KEEPNANFRACP) \
214 _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs); \
215 X##_s = _FP_NANSIGN_##fs; \
217 else \
218 _FP_SETQNAN (fs, wc, X); \
221 while (0)
224 * Before packing the bits back into the native fp result, take care
225 * of such mundane things as rounding and overflow. Also, for some
226 * kinds of fp values, the original parts may not have been fully
227 * extracted -- but that is ok, we can regenerate them now.
230 #define _FP_PACK_CANONICAL(fs, wc, X) \
231 do \
233 switch (X##_c) \
235 case FP_CLS_NORMAL: \
236 X##_e += _FP_EXPBIAS_##fs; \
237 if (X##_e > 0) \
239 _FP_ROUND (wc, X); \
240 if (_FP_FRAC_OVERP_##wc (fs, X)) \
242 _FP_FRAC_CLEAR_OVERP_##wc (fs, X); \
243 X##_e++; \
245 _FP_FRAC_SRL_##wc (X, _FP_WORKBITS); \
246 if (X##_e >= _FP_EXPMAX_##fs) \
248 /* overflow */ \
249 switch (FP_ROUNDMODE) \
251 case FP_RND_NEAREST: \
252 X##_c = FP_CLS_INF; \
253 break; \
254 case FP_RND_PINF: \
255 if (!X##_s) \
256 X##_c = FP_CLS_INF; \
257 break; \
258 case FP_RND_MINF: \
259 if (X##_s) \
260 X##_c = FP_CLS_INF; \
261 break; \
263 if (X##_c == FP_CLS_INF) \
265 /* Overflow to infinity */ \
266 X##_e = _FP_EXPMAX_##fs; \
267 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
269 else \
271 /* Overflow to maximum normal */ \
272 X##_e = _FP_EXPMAX_##fs - 1; \
273 _FP_FRAC_SET_##wc (X, _FP_MAXFRAC_##wc); \
275 FP_SET_EXCEPTION (FP_EX_OVERFLOW); \
276 FP_SET_EXCEPTION (FP_EX_INEXACT); \
279 else \
281 /* we've got a denormalized number */ \
282 X##_e = -X##_e + 1; \
283 if (X##_e <= _FP_WFRACBITS_##fs) \
285 _FP_FRAC_SRS_##wc (X, X##_e, _FP_WFRACBITS_##fs); \
286 _FP_ROUND (wc, X); \
287 if (_FP_FRAC_HIGH_##fs (X) \
288 & (_FP_OVERFLOW_##fs >> 1)) \
290 X##_e = 1; \
291 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
292 FP_SET_EXCEPTION (FP_EX_INEXACT); \
294 else \
296 X##_e = 0; \
297 _FP_FRAC_SRL_##wc (X, _FP_WORKBITS); \
299 if ((FP_CUR_EXCEPTIONS & FP_EX_INEXACT) \
300 || (FP_TRAPPING_EXCEPTIONS & FP_EX_UNDERFLOW)) \
301 FP_SET_EXCEPTION (FP_EX_UNDERFLOW); \
303 else \
305 /* underflow to zero */ \
306 X##_e = 0; \
307 if (!_FP_FRAC_ZEROP_##wc (X)) \
309 _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc); \
310 _FP_ROUND (wc, X); \
311 _FP_FRAC_LOW_##wc (X) >>= (_FP_WORKBITS); \
313 FP_SET_EXCEPTION (FP_EX_UNDERFLOW); \
316 break; \
318 case FP_CLS_ZERO: \
319 X##_e = 0; \
320 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
321 break; \
323 case FP_CLS_INF: \
324 X##_e = _FP_EXPMAX_##fs; \
325 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
326 break; \
328 case FP_CLS_NAN: \
329 X##_e = _FP_EXPMAX_##fs; \
330 if (!_FP_KEEPNANFRACP) \
332 _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs); \
333 X##_s = _FP_NANSIGN_##fs; \
335 else \
336 _FP_SETQNAN (fs, wc, X); \
337 break; \
340 while (0)
342 /* This one accepts raw argument and not cooked, returns
343 * 1 if X is a signaling NaN.
345 #define _FP_ISSIGNAN(fs, wc, X) \
346 ({ \
347 int __ret = 0; \
348 if (X##_e == _FP_EXPMAX_##fs) \
350 if (!_FP_FRAC_ZEROP_##wc (X) \
351 && _FP_FRAC_SNANP (fs, X)) \
352 __ret = 1; \
354 __ret; \
361 /* Addition on semi-raw values. */
362 #define _FP_ADD_INTERNAL(fs, wc, R, X, Y, OP) \
363 do \
365 if (X##_s == Y##_s) \
367 /* Addition. */ \
368 R##_s = X##_s; \
369 int ediff = X##_e - Y##_e; \
370 if (ediff > 0) \
372 R##_e = X##_e; \
373 if (Y##_e == 0) \
375 /* Y is zero or denormalized. */ \
376 if (_FP_FRAC_ZEROP_##wc (Y)) \
378 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
379 _FP_FRAC_COPY_##wc (R, X); \
380 goto add_done; \
382 else \
384 FP_SET_EXCEPTION (FP_EX_DENORM); \
385 ediff--; \
386 if (ediff == 0) \
388 _FP_FRAC_ADD_##wc (R, X, Y); \
389 goto add3; \
391 if (X##_e == _FP_EXPMAX_##fs) \
393 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
394 _FP_FRAC_COPY_##wc (R, X); \
395 goto add_done; \
397 goto add1; \
400 else if (X##_e == _FP_EXPMAX_##fs) \
402 /* X is NaN or Inf, Y is normal. */ \
403 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
404 _FP_FRAC_COPY_##wc (R, X); \
405 goto add_done; \
408 /* Insert implicit MSB of Y. */ \
409 _FP_FRAC_HIGH_##fs (Y) |= _FP_IMPLBIT_SH_##fs; \
411 add1: \
412 /* Shift the mantissa of Y to the right EDIFF steps; \
413 remember to account later for the implicit MSB of X. */ \
414 if (ediff <= _FP_WFRACBITS_##fs) \
415 _FP_FRAC_SRS_##wc (Y, ediff, _FP_WFRACBITS_##fs); \
416 else if (!_FP_FRAC_ZEROP_##wc (Y)) \
417 _FP_FRAC_SET_##wc (Y, _FP_MINFRAC_##wc); \
418 _FP_FRAC_ADD_##wc (R, X, Y); \
420 else if (ediff < 0) \
422 ediff = -ediff; \
423 R##_e = Y##_e; \
424 if (X##_e == 0) \
426 /* X is zero or denormalized. */ \
427 if (_FP_FRAC_ZEROP_##wc (X)) \
429 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
430 _FP_FRAC_COPY_##wc (R, Y); \
431 goto add_done; \
433 else \
435 FP_SET_EXCEPTION (FP_EX_DENORM); \
436 ediff--; \
437 if (ediff == 0) \
439 _FP_FRAC_ADD_##wc (R, Y, X); \
440 goto add3; \
442 if (Y##_e == _FP_EXPMAX_##fs) \
444 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
445 _FP_FRAC_COPY_##wc (R, Y); \
446 goto add_done; \
448 goto add2; \
451 else if (Y##_e == _FP_EXPMAX_##fs) \
453 /* Y is NaN or Inf, X is normal. */ \
454 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
455 _FP_FRAC_COPY_##wc (R, Y); \
456 goto add_done; \
459 /* Insert implicit MSB of X. */ \
460 _FP_FRAC_HIGH_##fs (X) |= _FP_IMPLBIT_SH_##fs; \
462 add2: \
463 /* Shift the mantissa of X to the right EDIFF steps; \
464 remember to account later for the implicit MSB of Y. */ \
465 if (ediff <= _FP_WFRACBITS_##fs) \
466 _FP_FRAC_SRS_##wc (X, ediff, _FP_WFRACBITS_##fs); \
467 else if (!_FP_FRAC_ZEROP_##wc (X)) \
468 _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc); \
469 _FP_FRAC_ADD_##wc (R, Y, X); \
471 else \
473 /* ediff == 0. */ \
474 if (!_FP_EXP_NORMAL (fs, wc, X)) \
476 if (X##_e == 0) \
478 /* X and Y are zero or denormalized. */ \
479 R##_e = 0; \
480 if (_FP_FRAC_ZEROP_##wc (X)) \
482 if (!_FP_FRAC_ZEROP_##wc (Y)) \
483 FP_SET_EXCEPTION (FP_EX_DENORM); \
484 _FP_FRAC_COPY_##wc (R, Y); \
485 goto add_done; \
487 else if (_FP_FRAC_ZEROP_##wc (Y)) \
489 FP_SET_EXCEPTION (FP_EX_DENORM); \
490 _FP_FRAC_COPY_##wc (R, X); \
491 goto add_done; \
493 else \
495 FP_SET_EXCEPTION (FP_EX_DENORM); \
496 _FP_FRAC_ADD_##wc (R, X, Y); \
497 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
499 /* Normalized result. */ \
500 _FP_FRAC_HIGH_##fs (R) \
501 &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
502 R##_e = 1; \
504 goto add_done; \
507 else \
509 /* X and Y are NaN or Inf. */ \
510 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
511 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
512 R##_e = _FP_EXPMAX_##fs; \
513 if (_FP_FRAC_ZEROP_##wc (X)) \
514 _FP_FRAC_COPY_##wc (R, Y); \
515 else if (_FP_FRAC_ZEROP_##wc (Y)) \
516 _FP_FRAC_COPY_##wc (R, X); \
517 else \
518 _FP_CHOOSENAN_SEMIRAW (fs, wc, R, X, Y, OP); \
519 goto add_done; \
522 /* The exponents of X and Y, both normal, are equal. The \
523 implicit MSBs will always add to increase the \
524 exponent. */ \
525 _FP_FRAC_ADD_##wc (R, X, Y); \
526 R##_e = X##_e + 1; \
527 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs); \
528 if (R##_e == _FP_EXPMAX_##fs) \
529 /* Overflow to infinity (depending on rounding mode). */ \
530 _FP_OVERFLOW_SEMIRAW (fs, wc, R); \
531 goto add_done; \
533 add3: \
534 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
536 /* Overflow. */ \
537 _FP_FRAC_HIGH_##fs (R) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
538 R##_e++; \
539 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs); \
540 if (R##_e == _FP_EXPMAX_##fs) \
541 /* Overflow to infinity (depending on rounding mode). */ \
542 _FP_OVERFLOW_SEMIRAW (fs, wc, R); \
544 add_done: ; \
546 else \
548 /* Subtraction. */ \
549 int ediff = X##_e - Y##_e; \
550 if (ediff > 0) \
552 R##_e = X##_e; \
553 R##_s = X##_s; \
554 if (Y##_e == 0) \
556 /* Y is zero or denormalized. */ \
557 if (_FP_FRAC_ZEROP_##wc (Y)) \
559 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
560 _FP_FRAC_COPY_##wc (R, X); \
561 goto sub_done; \
563 else \
565 FP_SET_EXCEPTION (FP_EX_DENORM); \
566 ediff--; \
567 if (ediff == 0) \
569 _FP_FRAC_SUB_##wc (R, X, Y); \
570 goto sub3; \
572 if (X##_e == _FP_EXPMAX_##fs) \
574 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
575 _FP_FRAC_COPY_##wc (R, X); \
576 goto sub_done; \
578 goto sub1; \
581 else if (X##_e == _FP_EXPMAX_##fs) \
583 /* X is NaN or Inf, Y is normal. */ \
584 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
585 _FP_FRAC_COPY_##wc (R, X); \
586 goto sub_done; \
589 /* Insert implicit MSB of Y. */ \
590 _FP_FRAC_HIGH_##fs (Y) |= _FP_IMPLBIT_SH_##fs; \
592 sub1: \
593 /* Shift the mantissa of Y to the right EDIFF steps; \
594 remember to account later for the implicit MSB of X. */ \
595 if (ediff <= _FP_WFRACBITS_##fs) \
596 _FP_FRAC_SRS_##wc (Y, ediff, _FP_WFRACBITS_##fs); \
597 else if (!_FP_FRAC_ZEROP_##wc (Y)) \
598 _FP_FRAC_SET_##wc (Y, _FP_MINFRAC_##wc); \
599 _FP_FRAC_SUB_##wc (R, X, Y); \
601 else if (ediff < 0) \
603 ediff = -ediff; \
604 R##_e = Y##_e; \
605 R##_s = Y##_s; \
606 if (X##_e == 0) \
608 /* X is zero or denormalized. */ \
609 if (_FP_FRAC_ZEROP_##wc (X)) \
611 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
612 _FP_FRAC_COPY_##wc (R, Y); \
613 goto sub_done; \
615 else \
617 FP_SET_EXCEPTION (FP_EX_DENORM); \
618 ediff--; \
619 if (ediff == 0) \
621 _FP_FRAC_SUB_##wc (R, Y, X); \
622 goto sub3; \
624 if (Y##_e == _FP_EXPMAX_##fs) \
626 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
627 _FP_FRAC_COPY_##wc (R, Y); \
628 goto sub_done; \
630 goto sub2; \
633 else if (Y##_e == _FP_EXPMAX_##fs) \
635 /* Y is NaN or Inf, X is normal. */ \
636 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
637 _FP_FRAC_COPY_##wc (R, Y); \
638 goto sub_done; \
641 /* Insert implicit MSB of X. */ \
642 _FP_FRAC_HIGH_##fs (X) |= _FP_IMPLBIT_SH_##fs; \
644 sub2: \
645 /* Shift the mantissa of X to the right EDIFF steps; \
646 remember to account later for the implicit MSB of Y. */ \
647 if (ediff <= _FP_WFRACBITS_##fs) \
648 _FP_FRAC_SRS_##wc (X, ediff, _FP_WFRACBITS_##fs); \
649 else if (!_FP_FRAC_ZEROP_##wc (X)) \
650 _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc); \
651 _FP_FRAC_SUB_##wc (R, Y, X); \
653 else \
655 /* ediff == 0. */ \
656 if (!_FP_EXP_NORMAL (fs, wc, X)) \
658 if (X##_e == 0) \
660 /* X and Y are zero or denormalized. */ \
661 R##_e = 0; \
662 if (_FP_FRAC_ZEROP_##wc (X)) \
664 _FP_FRAC_COPY_##wc (R, Y); \
665 if (_FP_FRAC_ZEROP_##wc (Y)) \
666 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
667 else \
669 FP_SET_EXCEPTION (FP_EX_DENORM); \
670 R##_s = Y##_s; \
672 goto sub_done; \
674 else if (_FP_FRAC_ZEROP_##wc (Y)) \
676 FP_SET_EXCEPTION (FP_EX_DENORM); \
677 _FP_FRAC_COPY_##wc (R, X); \
678 R##_s = X##_s; \
679 goto sub_done; \
681 else \
683 FP_SET_EXCEPTION (FP_EX_DENORM); \
684 _FP_FRAC_SUB_##wc (R, X, Y); \
685 R##_s = X##_s; \
686 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
688 /* |X| < |Y|, negate result. */ \
689 _FP_FRAC_SUB_##wc (R, Y, X); \
690 R##_s = Y##_s; \
692 else if (_FP_FRAC_ZEROP_##wc (R)) \
693 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
694 goto sub_done; \
697 else \
699 /* X and Y are NaN or Inf, of opposite signs. */ \
700 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
701 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
702 R##_e = _FP_EXPMAX_##fs; \
703 if (_FP_FRAC_ZEROP_##wc (X)) \
705 if (_FP_FRAC_ZEROP_##wc (Y)) \
707 /* Inf - Inf. */ \
708 R##_s = _FP_NANSIGN_##fs; \
709 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
710 _FP_FRAC_SLL_##wc (R, _FP_WORKBITS); \
711 FP_SET_EXCEPTION (FP_EX_INVALID); \
713 else \
715 /* Inf - NaN. */ \
716 R##_s = Y##_s; \
717 _FP_FRAC_COPY_##wc (R, Y); \
720 else \
722 if (_FP_FRAC_ZEROP_##wc (Y)) \
724 /* NaN - Inf. */ \
725 R##_s = X##_s; \
726 _FP_FRAC_COPY_##wc (R, X); \
728 else \
730 /* NaN - NaN. */ \
731 _FP_CHOOSENAN_SEMIRAW (fs, wc, R, X, Y, OP); \
734 goto sub_done; \
737 /* The exponents of X and Y, both normal, are equal. The \
738 implicit MSBs cancel. */ \
739 R##_e = X##_e; \
740 _FP_FRAC_SUB_##wc (R, X, Y); \
741 R##_s = X##_s; \
742 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
744 /* |X| < |Y|, negate result. */ \
745 _FP_FRAC_SUB_##wc (R, Y, X); \
746 R##_s = Y##_s; \
748 else if (_FP_FRAC_ZEROP_##wc (R)) \
750 R##_e = 0; \
751 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
752 goto sub_done; \
754 goto norm; \
756 sub3: \
757 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
759 int diff; \
760 /* Carry into most significant bit of larger one of X and Y, \
761 canceling it; renormalize. */ \
762 _FP_FRAC_HIGH_##fs (R) &= _FP_IMPLBIT_SH_##fs - 1; \
763 norm: \
764 _FP_FRAC_CLZ_##wc (diff, R); \
765 diff -= _FP_WFRACXBITS_##fs; \
766 _FP_FRAC_SLL_##wc (R, diff); \
767 if (R##_e <= diff) \
769 /* R is denormalized. */ \
770 diff = diff - R##_e + 1; \
771 _FP_FRAC_SRS_##wc (R, diff, _FP_WFRACBITS_##fs); \
772 R##_e = 0; \
774 else \
776 R##_e -= diff; \
777 _FP_FRAC_HIGH_##fs (R) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
780 sub_done: ; \
783 while (0)
785 #define _FP_ADD(fs, wc, R, X, Y) _FP_ADD_INTERNAL (fs, wc, R, X, Y, '+')
786 #define _FP_SUB(fs, wc, R, X, Y) \
787 do \
789 if (!(Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))) \
790 Y##_s ^= 1; \
791 _FP_ADD_INTERNAL (fs, wc, R, X, Y, '-'); \
793 while (0)
797 * Main negation routine. The input value is raw.
800 #define _FP_NEG(fs, wc, R, X) \
801 do \
803 _FP_FRAC_COPY_##wc (R, X); \
804 R##_e = X##_e; \
805 R##_s = 1 ^ X##_s; \
807 while (0)
811 * Main multiplication routine. The input values should be cooked.
814 #define _FP_MUL(fs, wc, R, X, Y) \
815 do \
817 R##_s = X##_s ^ Y##_s; \
818 R##_e = X##_e + Y##_e + 1; \
819 switch (_FP_CLS_COMBINE (X##_c, Y##_c)) \
821 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL): \
822 R##_c = FP_CLS_NORMAL; \
824 _FP_MUL_MEAT_##fs (R, X, Y); \
826 if (_FP_FRAC_OVERP_##wc (fs, R)) \
827 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs); \
828 else \
829 R##_e--; \
830 break; \
832 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN): \
833 _FP_CHOOSENAN (fs, wc, R, X, Y, '*'); \
834 break; \
836 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL): \
837 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF): \
838 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO): \
839 R##_s = X##_s; \
841 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF): \
842 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL): \
843 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL): \
844 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO): \
845 _FP_FRAC_COPY_##wc (R, X); \
846 R##_c = X##_c; \
847 break; \
849 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN): \
850 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN): \
851 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN): \
852 R##_s = Y##_s; \
854 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF): \
855 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO): \
856 _FP_FRAC_COPY_##wc (R, Y); \
857 R##_c = Y##_c; \
858 break; \
860 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO): \
861 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF): \
862 R##_s = _FP_NANSIGN_##fs; \
863 R##_c = FP_CLS_NAN; \
864 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
865 FP_SET_EXCEPTION (FP_EX_INVALID); \
866 break; \
868 default: \
869 abort (); \
872 while (0)
875 /* Fused multiply-add. The input values should be cooked. */
877 #define _FP_FMA(fs, wc, dwc, R, X, Y, Z) \
878 do \
880 FP_DECL_##fs (T); \
881 T##_s = X##_s ^ Y##_s; \
882 T##_e = X##_e + Y##_e + 1; \
883 switch (_FP_CLS_COMBINE (X##_c, Y##_c)) \
885 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL): \
886 switch (Z##_c) \
888 case FP_CLS_INF: \
889 case FP_CLS_NAN: \
890 R##_s = Z##_s; \
891 _FP_FRAC_COPY_##wc (R, Z); \
892 R##_c = Z##_c; \
893 break; \
895 case FP_CLS_ZERO: \
896 R##_c = FP_CLS_NORMAL; \
897 R##_s = T##_s; \
898 R##_e = T##_e; \
900 _FP_MUL_MEAT_##fs (R, X, Y); \
902 if (_FP_FRAC_OVERP_##wc (fs, R)) \
903 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs); \
904 else \
905 R##_e--; \
906 break; \
908 case FP_CLS_NORMAL:; \
909 _FP_FRAC_DECL_##dwc (TD); \
910 _FP_FRAC_DECL_##dwc (ZD); \
911 _FP_FRAC_DECL_##dwc (RD); \
912 _FP_MUL_MEAT_DW_##fs (TD, X, Y); \
913 R##_e = T##_e; \
914 int tsh = _FP_FRAC_HIGHBIT_DW_##dwc (fs, TD) == 0; \
915 T##_e -= tsh; \
916 int ediff = T##_e - Z##_e; \
917 if (ediff >= 0) \
919 int shift = _FP_WFRACBITS_##fs - tsh - ediff; \
920 if (shift <= -_FP_WFRACBITS_##fs) \
921 _FP_FRAC_SET_##dwc (ZD, _FP_MINFRAC_##dwc); \
922 else \
924 _FP_FRAC_COPY_##dwc##_##wc (ZD, Z); \
925 if (shift < 0) \
926 _FP_FRAC_SRS_##dwc (ZD, -shift, \
927 _FP_WFRACBITS_DW_##fs); \
928 else if (shift > 0) \
929 _FP_FRAC_SLL_##dwc (ZD, shift); \
931 R##_s = T##_s; \
932 if (T##_s == Z##_s) \
933 _FP_FRAC_ADD_##dwc (RD, TD, ZD); \
934 else \
936 _FP_FRAC_SUB_##dwc (RD, TD, ZD); \
937 if (_FP_FRAC_NEGP_##dwc (RD)) \
939 R##_s = Z##_s; \
940 _FP_FRAC_SUB_##dwc (RD, ZD, TD); \
944 else \
946 R##_e = Z##_e; \
947 R##_s = Z##_s; \
948 _FP_FRAC_COPY_##dwc##_##wc (ZD, Z); \
949 _FP_FRAC_SLL_##dwc (ZD, _FP_WFRACBITS_##fs); \
950 int shift = -ediff - tsh; \
951 if (shift >= _FP_WFRACBITS_DW_##fs) \
952 _FP_FRAC_SET_##dwc (TD, _FP_MINFRAC_##dwc); \
953 else if (shift > 0) \
954 _FP_FRAC_SRS_##dwc (TD, shift, \
955 _FP_WFRACBITS_DW_##fs); \
956 if (Z##_s == T##_s) \
957 _FP_FRAC_ADD_##dwc (RD, ZD, TD); \
958 else \
959 _FP_FRAC_SUB_##dwc (RD, ZD, TD); \
961 if (_FP_FRAC_ZEROP_##dwc (RD)) \
963 if (T##_s == Z##_s) \
964 R##_s = Z##_s; \
965 else \
966 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
967 _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc); \
968 R##_c = FP_CLS_ZERO; \
970 else \
972 int rlz; \
973 _FP_FRAC_CLZ_##dwc (rlz, RD); \
974 rlz -= _FP_WFRACXBITS_DW_##fs; \
975 R##_e -= rlz; \
976 int shift = _FP_WFRACBITS_##fs - rlz; \
977 if (shift > 0) \
978 _FP_FRAC_SRS_##dwc (RD, shift, \
979 _FP_WFRACBITS_DW_##fs); \
980 else if (shift < 0) \
981 _FP_FRAC_SLL_##dwc (RD, -shift); \
982 _FP_FRAC_COPY_##wc##_##dwc (R, RD); \
983 R##_c = FP_CLS_NORMAL; \
985 break; \
987 goto done_fma; \
989 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN): \
990 _FP_CHOOSENAN (fs, wc, T, X, Y, '*'); \
991 break; \
993 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL): \
994 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF): \
995 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO): \
996 T##_s = X##_s; \
998 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF): \
999 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL): \
1000 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL): \
1001 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO): \
1002 _FP_FRAC_COPY_##wc (T, X); \
1003 T##_c = X##_c; \
1004 break; \
1006 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN): \
1007 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN): \
1008 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN): \
1009 T##_s = Y##_s; \
1011 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF): \
1012 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO): \
1013 _FP_FRAC_COPY_##wc (T, Y); \
1014 T##_c = Y##_c; \
1015 break; \
1017 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO): \
1018 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF): \
1019 T##_s = _FP_NANSIGN_##fs; \
1020 T##_c = FP_CLS_NAN; \
1021 _FP_FRAC_SET_##wc (T, _FP_NANFRAC_##fs); \
1022 FP_SET_EXCEPTION (FP_EX_INVALID); \
1023 break; \
1025 default: \
1026 abort (); \
1029 /* T = X * Y is zero, infinity or NaN. */ \
1030 switch (_FP_CLS_COMBINE (T##_c, Z##_c)) \
1032 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN): \
1033 _FP_CHOOSENAN (fs, wc, R, T, Z, '+'); \
1034 break; \
1036 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL): \
1037 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF): \
1038 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO): \
1039 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL): \
1040 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO): \
1041 R##_s = T##_s; \
1042 _FP_FRAC_COPY_##wc (R, T); \
1043 R##_c = T##_c; \
1044 break; \
1046 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN): \
1047 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN): \
1048 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL): \
1049 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF): \
1050 R##_s = Z##_s; \
1051 _FP_FRAC_COPY_##wc (R, Z); \
1052 R##_c = Z##_c; \
1053 break; \
1055 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF): \
1056 if (T##_s == Z##_s) \
1058 R##_s = Z##_s; \
1059 _FP_FRAC_COPY_##wc (R, Z); \
1060 R##_c = Z##_c; \
1062 else \
1064 R##_s = _FP_NANSIGN_##fs; \
1065 R##_c = FP_CLS_NAN; \
1066 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
1067 FP_SET_EXCEPTION (FP_EX_INVALID); \
1069 break; \
1071 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO): \
1072 if (T##_s == Z##_s) \
1073 R##_s = Z##_s; \
1074 else \
1075 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
1076 _FP_FRAC_COPY_##wc (R, Z); \
1077 R##_c = Z##_c; \
1078 break; \
1080 default: \
1081 abort (); \
1083 done_fma: ; \
1085 while (0)
1089 * Main division routine. The input values should be cooked.
1092 #define _FP_DIV(fs, wc, R, X, Y) \
1093 do \
1095 R##_s = X##_s ^ Y##_s; \
1096 R##_e = X##_e - Y##_e; \
1097 switch (_FP_CLS_COMBINE (X##_c, Y##_c)) \
1099 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL): \
1100 R##_c = FP_CLS_NORMAL; \
1102 _FP_DIV_MEAT_##fs (R, X, Y); \
1103 break; \
1105 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN): \
1106 _FP_CHOOSENAN (fs, wc, R, X, Y, '/'); \
1107 break; \
1109 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL): \
1110 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF): \
1111 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO): \
1112 R##_s = X##_s; \
1113 _FP_FRAC_COPY_##wc (R, X); \
1114 R##_c = X##_c; \
1115 break; \
1117 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN): \
1118 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN): \
1119 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN): \
1120 R##_s = Y##_s; \
1121 _FP_FRAC_COPY_##wc (R, Y); \
1122 R##_c = Y##_c; \
1123 break; \
1125 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF): \
1126 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF): \
1127 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL): \
1128 R##_c = FP_CLS_ZERO; \
1129 break; \
1131 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO): \
1132 FP_SET_EXCEPTION (FP_EX_DIVZERO); \
1133 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO): \
1134 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL): \
1135 R##_c = FP_CLS_INF; \
1136 break; \
1138 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF): \
1139 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO): \
1140 R##_s = _FP_NANSIGN_##fs; \
1141 R##_c = FP_CLS_NAN; \
1142 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
1143 FP_SET_EXCEPTION (FP_EX_INVALID); \
1144 break; \
1146 default: \
1147 abort (); \
1150 while (0)
1154 * Main differential comparison routine. The inputs should be raw not
1155 * cooked. The return is -1,0,1 for normal values, 2 otherwise.
1158 #define _FP_CMP(fs, wc, ret, X, Y, un) \
1159 do \
1161 /* NANs are unordered */ \
1162 if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
1163 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))) \
1165 ret = un; \
1167 else \
1169 int __is_zero_x; \
1170 int __is_zero_y; \
1172 __is_zero_x = (!X##_e && _FP_FRAC_ZEROP_##wc (X)) ? 1 : 0; \
1173 __is_zero_y = (!Y##_e && _FP_FRAC_ZEROP_##wc (Y)) ? 1 : 0; \
1175 if (__is_zero_x && __is_zero_y) \
1176 ret = 0; \
1177 else if (__is_zero_x) \
1178 ret = Y##_s ? 1 : -1; \
1179 else if (__is_zero_y) \
1180 ret = X##_s ? -1 : 1; \
1181 else if (X##_s != Y##_s) \
1182 ret = X##_s ? -1 : 1; \
1183 else if (X##_e > Y##_e) \
1184 ret = X##_s ? -1 : 1; \
1185 else if (X##_e < Y##_e) \
1186 ret = X##_s ? 1 : -1; \
1187 else if (_FP_FRAC_GT_##wc (X, Y)) \
1188 ret = X##_s ? -1 : 1; \
1189 else if (_FP_FRAC_GT_##wc (Y, X)) \
1190 ret = X##_s ? 1 : -1; \
1191 else \
1192 ret = 0; \
1195 while (0)
1198 /* Simplification for strict equality. */
1200 #define _FP_CMP_EQ(fs, wc, ret, X, Y) \
1201 do \
1203 /* NANs are unordered */ \
1204 if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
1205 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))) \
1207 ret = 1; \
1209 else \
1211 ret = !(X##_e == Y##_e \
1212 && _FP_FRAC_EQ_##wc (X, Y) \
1213 && (X##_s == Y##_s || (!X##_e && _FP_FRAC_ZEROP_##wc (X)))); \
1216 while (0)
1218 /* Version to test unordered. */
1220 #define _FP_CMP_UNORD(fs, wc, ret, X, Y) \
1221 do \
1223 ret = ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
1224 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))); \
1226 while (0)
1229 * Main square root routine. The input value should be cooked.
1232 #define _FP_SQRT(fs, wc, R, X) \
1233 do \
1235 _FP_FRAC_DECL_##wc (T); \
1236 _FP_FRAC_DECL_##wc (S); \
1237 _FP_W_TYPE q; \
1238 switch (X##_c) \
1240 case FP_CLS_NAN: \
1241 _FP_FRAC_COPY_##wc (R, X); \
1242 R##_s = X##_s; \
1243 R##_c = FP_CLS_NAN; \
1244 break; \
1245 case FP_CLS_INF: \
1246 if (X##_s) \
1248 R##_s = _FP_NANSIGN_##fs; \
1249 R##_c = FP_CLS_NAN; /* NAN */ \
1250 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
1251 FP_SET_EXCEPTION (FP_EX_INVALID); \
1253 else \
1255 R##_s = 0; \
1256 R##_c = FP_CLS_INF; /* sqrt(+inf) = +inf */ \
1258 break; \
1259 case FP_CLS_ZERO: \
1260 R##_s = X##_s; \
1261 R##_c = FP_CLS_ZERO; /* sqrt(+-0) = +-0 */ \
1262 break; \
1263 case FP_CLS_NORMAL: \
1264 R##_s = 0; \
1265 if (X##_s) \
1267 R##_c = FP_CLS_NAN; /* NAN */ \
1268 R##_s = _FP_NANSIGN_##fs; \
1269 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
1270 FP_SET_EXCEPTION (FP_EX_INVALID); \
1271 break; \
1273 R##_c = FP_CLS_NORMAL; \
1274 if (X##_e & 1) \
1275 _FP_FRAC_SLL_##wc (X, 1); \
1276 R##_e = X##_e >> 1; \
1277 _FP_FRAC_SET_##wc (S, _FP_ZEROFRAC_##wc); \
1278 _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc); \
1279 q = _FP_OVERFLOW_##fs >> 1; \
1280 _FP_SQRT_MEAT_##wc (R, S, T, X, q); \
1283 while (0)
1286 * Convert from FP to integer. Input is raw.
1289 /* RSIGNED can have following values:
1290 * 0: the number is required to be 0..(2^rsize)-1, if not, NV is set plus
1291 * the result is either 0 or (2^rsize)-1 depending on the sign in such
1292 * case.
1293 * 1: the number is required to be -(2^(rsize-1))..(2^(rsize-1))-1, if not,
1294 * NV is set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
1295 * depending on the sign in such case.
1296 * -1: the number is required to be -(2^(rsize-1))..(2^rsize)-1, if not, NV is
1297 * set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
1298 * depending on the sign in such case.
1300 #define _FP_TO_INT(fs, wc, r, X, rsize, rsigned) \
1301 do \
1303 if (X##_e < _FP_EXPBIAS_##fs) \
1305 r = 0; \
1306 if (X##_e == 0) \
1308 if (!_FP_FRAC_ZEROP_##wc (X)) \
1310 FP_SET_EXCEPTION (FP_EX_INEXACT); \
1311 FP_SET_EXCEPTION (FP_EX_DENORM); \
1314 else \
1315 FP_SET_EXCEPTION (FP_EX_INEXACT); \
1317 else if (X##_e >= _FP_EXPBIAS_##fs + rsize - (rsigned > 0 || X##_s) \
1318 || (!rsigned && X##_s)) \
1320 /* Overflow or converting to the most negative integer. */ \
1321 if (rsigned) \
1323 r = 1; \
1324 r <<= rsize - 1; \
1325 r -= 1 - X##_s; \
1326 } else { \
1327 r = 0; \
1328 if (!X##_s) \
1329 r = ~r; \
1332 if (rsigned && X##_s && X##_e == _FP_EXPBIAS_##fs + rsize - 1) \
1334 /* Possibly converting to most negative integer; check the \
1335 mantissa. */ \
1336 int inexact = 0; \
1337 (void) ((_FP_FRACBITS_##fs > rsize) \
1338 ? ({ \
1339 _FP_FRAC_SRST_##wc (X, inexact, \
1340 _FP_FRACBITS_##fs - rsize, \
1341 _FP_FRACBITS_##fs); \
1342 0; \
1343 }) \
1344 : 0); \
1345 if (!_FP_FRAC_ZEROP_##wc (X)) \
1346 FP_SET_EXCEPTION (FP_EX_INVALID); \
1347 else if (inexact) \
1348 FP_SET_EXCEPTION (FP_EX_INEXACT); \
1350 else \
1351 FP_SET_EXCEPTION (FP_EX_INVALID); \
1353 else \
1355 _FP_FRAC_HIGH_RAW_##fs (X) |= _FP_IMPLBIT_##fs; \
1356 if (X##_e >= _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs - 1) \
1358 _FP_FRAC_ASSEMBLE_##wc (r, X, rsize); \
1359 r <<= X##_e - _FP_EXPBIAS_##fs - _FP_FRACBITS_##fs + 1; \
1361 else \
1363 int inexact; \
1364 _FP_FRAC_SRST_##wc (X, inexact, \
1365 (_FP_FRACBITS_##fs + _FP_EXPBIAS_##fs - 1 \
1366 - X##_e), \
1367 _FP_FRACBITS_##fs); \
1368 if (inexact) \
1369 FP_SET_EXCEPTION (FP_EX_INEXACT); \
1370 _FP_FRAC_ASSEMBLE_##wc (r, X, rsize); \
1372 if (rsigned && X##_s) \
1373 r = -r; \
1376 while (0)
1378 /* Convert integer to fp. Output is raw. RTYPE is unsigned even if
1379 input is signed. */
1380 #define _FP_FROM_INT(fs, wc, X, r, rsize, rtype) \
1381 do \
1383 if (r) \
1385 rtype ur_; \
1387 if ((X##_s = (r < 0))) \
1388 r = -(rtype) r; \
1390 ur_ = (rtype) r; \
1391 (void) ((rsize <= _FP_W_TYPE_SIZE) \
1392 ? ({ \
1393 int lz_; \
1394 __FP_CLZ (lz_, (_FP_W_TYPE) ur_); \
1395 X##_e = _FP_EXPBIAS_##fs + _FP_W_TYPE_SIZE - 1 - lz_; \
1396 }) \
1397 : ((rsize <= 2 * _FP_W_TYPE_SIZE) \
1398 ? ({ \
1399 int lz_; \
1400 __FP_CLZ_2 (lz_, \
1401 (_FP_W_TYPE) (ur_ >> _FP_W_TYPE_SIZE), \
1402 (_FP_W_TYPE) ur_); \
1403 X##_e = (_FP_EXPBIAS_##fs + 2 * _FP_W_TYPE_SIZE - 1 \
1404 - lz_); \
1405 }) \
1406 : (abort (), 0))); \
1408 if (rsize - 1 + _FP_EXPBIAS_##fs >= _FP_EXPMAX_##fs \
1409 && X##_e >= _FP_EXPMAX_##fs) \
1411 /* Exponent too big; overflow to infinity. (May also \
1412 happen after rounding below.) */ \
1413 _FP_OVERFLOW_SEMIRAW (fs, wc, X); \
1414 goto pack_semiraw; \
1417 if (rsize <= _FP_FRACBITS_##fs \
1418 || X##_e < _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs) \
1420 /* Exactly representable; shift left. */ \
1421 _FP_FRAC_DISASSEMBLE_##wc (X, ur_, rsize); \
1422 if (_FP_EXPBIAS_##fs + _FP_FRACBITS_##fs - 1 - X##_e > 0) \
1423 _FP_FRAC_SLL_##wc (X, (_FP_EXPBIAS_##fs \
1424 + _FP_FRACBITS_##fs - 1 - X##_e)); \
1426 else \
1428 /* More bits in integer than in floating type; need to \
1429 round. */ \
1430 if (_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 < X##_e) \
1431 ur_ = ((ur_ >> (X##_e - _FP_EXPBIAS_##fs \
1432 - _FP_WFRACBITS_##fs + 1)) \
1433 | ((ur_ << (rsize - (X##_e - _FP_EXPBIAS_##fs \
1434 - _FP_WFRACBITS_##fs + 1))) \
1435 != 0)); \
1436 _FP_FRAC_DISASSEMBLE_##wc (X, ur_, rsize); \
1437 if ((_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 - X##_e) > 0) \
1438 _FP_FRAC_SLL_##wc (X, (_FP_EXPBIAS_##fs \
1439 + _FP_WFRACBITS_##fs - 1 - X##_e)); \
1440 _FP_FRAC_HIGH_##fs (X) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
1441 pack_semiraw: \
1442 _FP_PACK_SEMIRAW (fs, wc, X); \
1445 else \
1447 X##_s = 0; \
1448 X##_e = 0; \
1449 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
1452 while (0)
1455 /* Extend from a narrower floating-point format to a wider one. Input
1456 and output are raw. */
1457 #define FP_EXTEND(dfs, sfs, dwc, swc, D, S) \
1458 do \
1460 if (_FP_FRACBITS_##dfs < _FP_FRACBITS_##sfs \
1461 || (_FP_EXPMAX_##dfs - _FP_EXPBIAS_##dfs \
1462 < _FP_EXPMAX_##sfs - _FP_EXPBIAS_##sfs) \
1463 || (_FP_EXPBIAS_##dfs < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1 \
1464 && _FP_EXPBIAS_##dfs != _FP_EXPBIAS_##sfs)) \
1465 abort (); \
1466 D##_s = S##_s; \
1467 _FP_FRAC_COPY_##dwc##_##swc (D, S); \
1468 if (_FP_EXP_NORMAL (sfs, swc, S)) \
1470 D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs; \
1471 _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs - _FP_FRACBITS_##sfs)); \
1473 else \
1475 if (S##_e == 0) \
1477 if (_FP_FRAC_ZEROP_##swc (S)) \
1478 D##_e = 0; \
1479 else if (_FP_EXPBIAS_##dfs \
1480 < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1) \
1482 FP_SET_EXCEPTION (FP_EX_DENORM); \
1483 _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs \
1484 - _FP_FRACBITS_##sfs)); \
1485 D##_e = 0; \
1487 else \
1489 int _lz; \
1490 FP_SET_EXCEPTION (FP_EX_DENORM); \
1491 _FP_FRAC_CLZ_##swc (_lz, S); \
1492 _FP_FRAC_SLL_##dwc (D, \
1493 _lz + _FP_FRACBITS_##dfs \
1494 - _FP_FRACTBITS_##sfs); \
1495 D##_e = (_FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs + 1 \
1496 + _FP_FRACXBITS_##sfs - _lz); \
1499 else \
1501 D##_e = _FP_EXPMAX_##dfs; \
1502 if (!_FP_FRAC_ZEROP_##swc (S)) \
1504 if (_FP_FRAC_SNANP (sfs, S)) \
1505 FP_SET_EXCEPTION (FP_EX_INVALID); \
1506 _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs \
1507 - _FP_FRACBITS_##sfs)); \
1508 _FP_SETQNAN (dfs, dwc, D); \
1513 while (0)
1515 /* Truncate from a wider floating-point format to a narrower one.
1516 Input and output are semi-raw. */
1517 #define FP_TRUNC(dfs, sfs, dwc, swc, D, S) \
1518 do \
1520 if (_FP_FRACBITS_##sfs < _FP_FRACBITS_##dfs \
1521 || (_FP_EXPBIAS_##sfs < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1 \
1522 && _FP_EXPBIAS_##sfs != _FP_EXPBIAS_##dfs)) \
1523 abort (); \
1524 D##_s = S##_s; \
1525 if (_FP_EXP_NORMAL (sfs, swc, S)) \
1527 D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs; \
1528 if (D##_e >= _FP_EXPMAX_##dfs) \
1529 _FP_OVERFLOW_SEMIRAW (dfs, dwc, D); \
1530 else \
1532 if (D##_e <= 0) \
1534 if (D##_e < 1 - _FP_FRACBITS_##dfs) \
1536 _FP_FRAC_SET_##swc (S, _FP_ZEROFRAC_##swc); \
1537 _FP_FRAC_LOW_##swc (S) |= 1; \
1539 else \
1541 _FP_FRAC_HIGH_##sfs (S) |= _FP_IMPLBIT_SH_##sfs; \
1542 _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs \
1543 - _FP_WFRACBITS_##dfs \
1544 + 1 - D##_e), \
1545 _FP_WFRACBITS_##sfs); \
1547 D##_e = 0; \
1549 else \
1550 _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs \
1551 - _FP_WFRACBITS_##dfs), \
1552 _FP_WFRACBITS_##sfs); \
1553 _FP_FRAC_COPY_##dwc##_##swc (D, S); \
1556 else \
1558 if (S##_e == 0) \
1560 D##_e = 0; \
1561 if (_FP_FRAC_ZEROP_##swc (S)) \
1562 _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc); \
1563 else \
1565 FP_SET_EXCEPTION (FP_EX_DENORM); \
1566 if (_FP_EXPBIAS_##sfs \
1567 < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1) \
1569 _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs \
1570 - _FP_WFRACBITS_##dfs), \
1571 _FP_WFRACBITS_##sfs); \
1572 _FP_FRAC_COPY_##dwc##_##swc (D, S); \
1574 else \
1576 _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc); \
1577 _FP_FRAC_LOW_##dwc (D) |= 1; \
1581 else \
1583 D##_e = _FP_EXPMAX_##dfs; \
1584 if (_FP_FRAC_ZEROP_##swc (S)) \
1585 _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc); \
1586 else \
1588 _FP_CHECK_SIGNAN_SEMIRAW (sfs, swc, S); \
1589 _FP_FRAC_SRL_##swc (S, (_FP_WFRACBITS_##sfs \
1590 - _FP_WFRACBITS_##dfs)); \
1591 _FP_FRAC_COPY_##dwc##_##swc (D, S); \
1592 /* Semi-raw NaN must have all workbits cleared. */ \
1593 _FP_FRAC_LOW_##dwc (D) \
1594 &= ~(_FP_W_TYPE) ((1 << _FP_WORKBITS) - 1); \
1595 _FP_SETQNAN_SEMIRAW (dfs, dwc, D); \
1600 while (0)
1603 * Helper primitives.
1606 /* Count leading zeros in a word. */
1608 #ifndef __FP_CLZ
1609 /* GCC 3.4 and later provide the builtins for us. */
1610 # define __FP_CLZ(r, x) \
1611 do \
1613 if (sizeof (_FP_W_TYPE) == sizeof (unsigned int)) \
1614 r = __builtin_clz (x); \
1615 else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long)) \
1616 r = __builtin_clzl (x); \
1617 else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long long)) \
1618 r = __builtin_clzll (x); \
1619 else \
1620 abort (); \
1622 while (0)
1623 #endif /* ndef __FP_CLZ */
1625 #define _FP_DIV_HELP_imm(q, r, n, d) \
1626 do \
1628 q = n / d, r = n % d; \
1630 while (0)
1633 /* A restoring bit-by-bit division primitive. */
1635 #define _FP_DIV_MEAT_N_loop(fs, wc, R, X, Y) \
1636 do \
1638 int count = _FP_WFRACBITS_##fs; \
1639 _FP_FRAC_DECL_##wc (u); \
1640 _FP_FRAC_DECL_##wc (v); \
1641 _FP_FRAC_COPY_##wc (u, X); \
1642 _FP_FRAC_COPY_##wc (v, Y); \
1643 _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc); \
1644 /* Normalize U and V. */ \
1645 _FP_FRAC_SLL_##wc (u, _FP_WFRACXBITS_##fs); \
1646 _FP_FRAC_SLL_##wc (v, _FP_WFRACXBITS_##fs); \
1647 /* First round. Since the operands are normalized, either the \
1648 first or second bit will be set in the fraction. Produce a \
1649 normalized result by checking which and adjusting the loop \
1650 count and exponent accordingly. */ \
1651 if (_FP_FRAC_GE_1 (u, v)) \
1653 _FP_FRAC_SUB_##wc (u, u, v); \
1654 _FP_FRAC_LOW_##wc (R) |= 1; \
1655 count--; \
1657 else \
1658 R##_e--; \
1659 /* Subsequent rounds. */ \
1660 do \
1662 int msb = (_FP_WS_TYPE) _FP_FRAC_HIGH_##wc (u) < 0; \
1663 _FP_FRAC_SLL_##wc (u, 1); \
1664 _FP_FRAC_SLL_##wc (R, 1); \
1665 if (msb || _FP_FRAC_GE_1 (u, v)) \
1667 _FP_FRAC_SUB_##wc (u, u, v); \
1668 _FP_FRAC_LOW_##wc (R) |= 1; \
1671 while (--count > 0); \
1672 /* If there's anything left in U, the result is inexact. */ \
1673 _FP_FRAC_LOW_##wc (R) |= !_FP_FRAC_ZEROP_##wc (u); \
1675 while (0)
1677 #define _FP_DIV_MEAT_1_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 1, R, X, Y)
1678 #define _FP_DIV_MEAT_2_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 2, R, X, Y)
1679 #define _FP_DIV_MEAT_4_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 4, R, X, Y)