Replace __int128 with __int128_t in bits/link.h
[glibc.git] / soft-fp / op-common.h
blob6f8c1f3cb04c1300a37be433b528c44d5d464103
1 /* Software floating-point emulation. Common operations.
2 Copyright (C) 1997-2014 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 __attribute__ ((unused)); \
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 int _FP_PACK_SEMIRAW_is_tiny \
195 = X##_e == 0 && !_FP_FRAC_ZEROP_##wc (X); \
196 if (_FP_TININESS_AFTER_ROUNDING \
197 && _FP_PACK_SEMIRAW_is_tiny) \
199 FP_DECL_##fs (_FP_PACK_SEMIRAW_T); \
200 _FP_FRAC_COPY_##wc (_FP_PACK_SEMIRAW_T, X); \
201 _FP_PACK_SEMIRAW_T##_s = X##_s; \
202 _FP_PACK_SEMIRAW_T##_e = X##_e; \
203 _FP_FRAC_SLL_##wc (_FP_PACK_SEMIRAW_T, 1); \
204 _FP_ROUND (wc, _FP_PACK_SEMIRAW_T); \
205 if (_FP_FRAC_OVERP_##wc (fs, _FP_PACK_SEMIRAW_T)) \
206 _FP_PACK_SEMIRAW_is_tiny = 0; \
208 _FP_ROUND (wc, X); \
209 if (_FP_PACK_SEMIRAW_is_tiny) \
211 if ((FP_CUR_EXCEPTIONS & FP_EX_INEXACT) \
212 || (FP_TRAPPING_EXCEPTIONS & FP_EX_UNDERFLOW)) \
213 FP_SET_EXCEPTION (FP_EX_UNDERFLOW); \
215 if (_FP_FRAC_HIGH_##fs (X) \
216 & (_FP_OVERFLOW_##fs >> 1)) \
218 _FP_FRAC_HIGH_##fs (X) &= ~(_FP_OVERFLOW_##fs >> 1); \
219 X##_e++; \
220 if (X##_e == _FP_EXPMAX_##fs) \
221 _FP_OVERFLOW_SEMIRAW (fs, wc, X); \
223 _FP_FRAC_SRL_##wc (X, _FP_WORKBITS); \
224 if (X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
226 if (!_FP_KEEPNANFRACP) \
228 _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs); \
229 X##_s = _FP_NANSIGN_##fs; \
231 else \
232 _FP_SETQNAN (fs, wc, X); \
235 while (0)
238 * Before packing the bits back into the native fp result, take care
239 * of such mundane things as rounding and overflow. Also, for some
240 * kinds of fp values, the original parts may not have been fully
241 * extracted -- but that is ok, we can regenerate them now.
244 #define _FP_PACK_CANONICAL(fs, wc, X) \
245 do \
247 switch (X##_c) \
249 case FP_CLS_NORMAL: \
250 X##_e += _FP_EXPBIAS_##fs; \
251 if (X##_e > 0) \
253 _FP_ROUND (wc, X); \
254 if (_FP_FRAC_OVERP_##wc (fs, X)) \
256 _FP_FRAC_CLEAR_OVERP_##wc (fs, X); \
257 X##_e++; \
259 _FP_FRAC_SRL_##wc (X, _FP_WORKBITS); \
260 if (X##_e >= _FP_EXPMAX_##fs) \
262 /* overflow */ \
263 switch (FP_ROUNDMODE) \
265 case FP_RND_NEAREST: \
266 X##_c = FP_CLS_INF; \
267 break; \
268 case FP_RND_PINF: \
269 if (!X##_s) \
270 X##_c = FP_CLS_INF; \
271 break; \
272 case FP_RND_MINF: \
273 if (X##_s) \
274 X##_c = FP_CLS_INF; \
275 break; \
277 if (X##_c == FP_CLS_INF) \
279 /* Overflow to infinity */ \
280 X##_e = _FP_EXPMAX_##fs; \
281 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
283 else \
285 /* Overflow to maximum normal */ \
286 X##_e = _FP_EXPMAX_##fs - 1; \
287 _FP_FRAC_SET_##wc (X, _FP_MAXFRAC_##wc); \
289 FP_SET_EXCEPTION (FP_EX_OVERFLOW); \
290 FP_SET_EXCEPTION (FP_EX_INEXACT); \
293 else \
295 /* we've got a denormalized number */ \
296 int _FP_PACK_CANONICAL_is_tiny = 1; \
297 if (_FP_TININESS_AFTER_ROUNDING && X##_e == 0) \
299 FP_DECL_##fs (_FP_PACK_CANONICAL_T); \
300 _FP_FRAC_COPY_##wc (_FP_PACK_CANONICAL_T, X); \
301 _FP_PACK_CANONICAL_T##_s = X##_s; \
302 _FP_PACK_CANONICAL_T##_e = X##_e; \
303 _FP_ROUND (wc, _FP_PACK_CANONICAL_T); \
304 if (_FP_FRAC_OVERP_##wc (fs, _FP_PACK_CANONICAL_T)) \
305 _FP_PACK_CANONICAL_is_tiny = 0; \
307 X##_e = -X##_e + 1; \
308 if (X##_e <= _FP_WFRACBITS_##fs) \
310 _FP_FRAC_SRS_##wc (X, X##_e, _FP_WFRACBITS_##fs); \
311 _FP_ROUND (wc, X); \
312 if (_FP_FRAC_HIGH_##fs (X) \
313 & (_FP_OVERFLOW_##fs >> 1)) \
315 X##_e = 1; \
316 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
317 FP_SET_EXCEPTION (FP_EX_INEXACT); \
319 else \
321 X##_e = 0; \
322 _FP_FRAC_SRL_##wc (X, _FP_WORKBITS); \
324 if (_FP_PACK_CANONICAL_is_tiny \
325 && ((FP_CUR_EXCEPTIONS & FP_EX_INEXACT) \
326 || (FP_TRAPPING_EXCEPTIONS \
327 & FP_EX_UNDERFLOW))) \
328 FP_SET_EXCEPTION (FP_EX_UNDERFLOW); \
330 else \
332 /* underflow to zero */ \
333 X##_e = 0; \
334 if (!_FP_FRAC_ZEROP_##wc (X)) \
336 _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc); \
337 _FP_ROUND (wc, X); \
338 _FP_FRAC_LOW_##wc (X) >>= (_FP_WORKBITS); \
340 FP_SET_EXCEPTION (FP_EX_UNDERFLOW); \
343 break; \
345 case FP_CLS_ZERO: \
346 X##_e = 0; \
347 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
348 break; \
350 case FP_CLS_INF: \
351 X##_e = _FP_EXPMAX_##fs; \
352 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
353 break; \
355 case FP_CLS_NAN: \
356 X##_e = _FP_EXPMAX_##fs; \
357 if (!_FP_KEEPNANFRACP) \
359 _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs); \
360 X##_s = _FP_NANSIGN_##fs; \
362 else \
363 _FP_SETQNAN (fs, wc, X); \
364 break; \
367 while (0)
369 /* This one accepts raw argument and not cooked, returns
370 * 1 if X is a signaling NaN.
372 #define _FP_ISSIGNAN(fs, wc, X) \
373 ({ \
374 int __ret = 0; \
375 if (X##_e == _FP_EXPMAX_##fs) \
377 if (!_FP_FRAC_ZEROP_##wc (X) \
378 && _FP_FRAC_SNANP (fs, X)) \
379 __ret = 1; \
381 __ret; \
388 /* Addition on semi-raw values. */
389 #define _FP_ADD_INTERNAL(fs, wc, R, X, Y, OP) \
390 do \
392 if (X##_s == Y##_s) \
394 /* Addition. */ \
395 R##_s = X##_s; \
396 int ediff = X##_e - Y##_e; \
397 if (ediff > 0) \
399 R##_e = X##_e; \
400 if (Y##_e == 0) \
402 /* Y is zero or denormalized. */ \
403 if (_FP_FRAC_ZEROP_##wc (Y)) \
405 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
406 _FP_FRAC_COPY_##wc (R, X); \
407 goto add_done; \
409 else \
411 FP_SET_EXCEPTION (FP_EX_DENORM); \
412 ediff--; \
413 if (ediff == 0) \
415 _FP_FRAC_ADD_##wc (R, X, Y); \
416 goto add3; \
418 if (X##_e == _FP_EXPMAX_##fs) \
420 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
421 _FP_FRAC_COPY_##wc (R, X); \
422 goto add_done; \
424 goto add1; \
427 else if (X##_e == _FP_EXPMAX_##fs) \
429 /* X is NaN or Inf, Y is normal. */ \
430 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
431 _FP_FRAC_COPY_##wc (R, X); \
432 goto add_done; \
435 /* Insert implicit MSB of Y. */ \
436 _FP_FRAC_HIGH_##fs (Y) |= _FP_IMPLBIT_SH_##fs; \
438 add1: \
439 /* Shift the mantissa of Y to the right EDIFF steps; \
440 remember to account later for the implicit MSB of X. */ \
441 if (ediff <= _FP_WFRACBITS_##fs) \
442 _FP_FRAC_SRS_##wc (Y, ediff, _FP_WFRACBITS_##fs); \
443 else if (!_FP_FRAC_ZEROP_##wc (Y)) \
444 _FP_FRAC_SET_##wc (Y, _FP_MINFRAC_##wc); \
445 _FP_FRAC_ADD_##wc (R, X, Y); \
447 else if (ediff < 0) \
449 ediff = -ediff; \
450 R##_e = Y##_e; \
451 if (X##_e == 0) \
453 /* X is zero or denormalized. */ \
454 if (_FP_FRAC_ZEROP_##wc (X)) \
456 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
457 _FP_FRAC_COPY_##wc (R, Y); \
458 goto add_done; \
460 else \
462 FP_SET_EXCEPTION (FP_EX_DENORM); \
463 ediff--; \
464 if (ediff == 0) \
466 _FP_FRAC_ADD_##wc (R, Y, X); \
467 goto add3; \
469 if (Y##_e == _FP_EXPMAX_##fs) \
471 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
472 _FP_FRAC_COPY_##wc (R, Y); \
473 goto add_done; \
475 goto add2; \
478 else if (Y##_e == _FP_EXPMAX_##fs) \
480 /* Y is NaN or Inf, X is normal. */ \
481 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
482 _FP_FRAC_COPY_##wc (R, Y); \
483 goto add_done; \
486 /* Insert implicit MSB of X. */ \
487 _FP_FRAC_HIGH_##fs (X) |= _FP_IMPLBIT_SH_##fs; \
489 add2: \
490 /* Shift the mantissa of X to the right EDIFF steps; \
491 remember to account later for the implicit MSB of Y. */ \
492 if (ediff <= _FP_WFRACBITS_##fs) \
493 _FP_FRAC_SRS_##wc (X, ediff, _FP_WFRACBITS_##fs); \
494 else if (!_FP_FRAC_ZEROP_##wc (X)) \
495 _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc); \
496 _FP_FRAC_ADD_##wc (R, Y, X); \
498 else \
500 /* ediff == 0. */ \
501 if (!_FP_EXP_NORMAL (fs, wc, X)) \
503 if (X##_e == 0) \
505 /* X and Y are zero or denormalized. */ \
506 R##_e = 0; \
507 if (_FP_FRAC_ZEROP_##wc (X)) \
509 if (!_FP_FRAC_ZEROP_##wc (Y)) \
510 FP_SET_EXCEPTION (FP_EX_DENORM); \
511 _FP_FRAC_COPY_##wc (R, Y); \
512 goto add_done; \
514 else if (_FP_FRAC_ZEROP_##wc (Y)) \
516 FP_SET_EXCEPTION (FP_EX_DENORM); \
517 _FP_FRAC_COPY_##wc (R, X); \
518 goto add_done; \
520 else \
522 FP_SET_EXCEPTION (FP_EX_DENORM); \
523 _FP_FRAC_ADD_##wc (R, X, Y); \
524 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
526 /* Normalized result. */ \
527 _FP_FRAC_HIGH_##fs (R) \
528 &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
529 R##_e = 1; \
531 goto add_done; \
534 else \
536 /* X and Y are NaN or Inf. */ \
537 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
538 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
539 R##_e = _FP_EXPMAX_##fs; \
540 if (_FP_FRAC_ZEROP_##wc (X)) \
541 _FP_FRAC_COPY_##wc (R, Y); \
542 else if (_FP_FRAC_ZEROP_##wc (Y)) \
543 _FP_FRAC_COPY_##wc (R, X); \
544 else \
545 _FP_CHOOSENAN_SEMIRAW (fs, wc, R, X, Y, OP); \
546 goto add_done; \
549 /* The exponents of X and Y, both normal, are equal. The \
550 implicit MSBs will always add to increase the \
551 exponent. */ \
552 _FP_FRAC_ADD_##wc (R, X, Y); \
553 R##_e = X##_e + 1; \
554 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs); \
555 if (R##_e == _FP_EXPMAX_##fs) \
556 /* Overflow to infinity (depending on rounding mode). */ \
557 _FP_OVERFLOW_SEMIRAW (fs, wc, R); \
558 goto add_done; \
560 add3: \
561 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
563 /* Overflow. */ \
564 _FP_FRAC_HIGH_##fs (R) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
565 R##_e++; \
566 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs); \
567 if (R##_e == _FP_EXPMAX_##fs) \
568 /* Overflow to infinity (depending on rounding mode). */ \
569 _FP_OVERFLOW_SEMIRAW (fs, wc, R); \
571 add_done: ; \
573 else \
575 /* Subtraction. */ \
576 int ediff = X##_e - Y##_e; \
577 if (ediff > 0) \
579 R##_e = X##_e; \
580 R##_s = X##_s; \
581 if (Y##_e == 0) \
583 /* Y is zero or denormalized. */ \
584 if (_FP_FRAC_ZEROP_##wc (Y)) \
586 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
587 _FP_FRAC_COPY_##wc (R, X); \
588 goto sub_done; \
590 else \
592 FP_SET_EXCEPTION (FP_EX_DENORM); \
593 ediff--; \
594 if (ediff == 0) \
596 _FP_FRAC_SUB_##wc (R, X, Y); \
597 goto sub3; \
599 if (X##_e == _FP_EXPMAX_##fs) \
601 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
602 _FP_FRAC_COPY_##wc (R, X); \
603 goto sub_done; \
605 goto sub1; \
608 else if (X##_e == _FP_EXPMAX_##fs) \
610 /* X is NaN or Inf, Y is normal. */ \
611 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
612 _FP_FRAC_COPY_##wc (R, X); \
613 goto sub_done; \
616 /* Insert implicit MSB of Y. */ \
617 _FP_FRAC_HIGH_##fs (Y) |= _FP_IMPLBIT_SH_##fs; \
619 sub1: \
620 /* Shift the mantissa of Y to the right EDIFF steps; \
621 remember to account later for the implicit MSB of X. */ \
622 if (ediff <= _FP_WFRACBITS_##fs) \
623 _FP_FRAC_SRS_##wc (Y, ediff, _FP_WFRACBITS_##fs); \
624 else if (!_FP_FRAC_ZEROP_##wc (Y)) \
625 _FP_FRAC_SET_##wc (Y, _FP_MINFRAC_##wc); \
626 _FP_FRAC_SUB_##wc (R, X, Y); \
628 else if (ediff < 0) \
630 ediff = -ediff; \
631 R##_e = Y##_e; \
632 R##_s = Y##_s; \
633 if (X##_e == 0) \
635 /* X is zero or denormalized. */ \
636 if (_FP_FRAC_ZEROP_##wc (X)) \
638 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
639 _FP_FRAC_COPY_##wc (R, Y); \
640 goto sub_done; \
642 else \
644 FP_SET_EXCEPTION (FP_EX_DENORM); \
645 ediff--; \
646 if (ediff == 0) \
648 _FP_FRAC_SUB_##wc (R, Y, X); \
649 goto sub3; \
651 if (Y##_e == _FP_EXPMAX_##fs) \
653 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
654 _FP_FRAC_COPY_##wc (R, Y); \
655 goto sub_done; \
657 goto sub2; \
660 else if (Y##_e == _FP_EXPMAX_##fs) \
662 /* Y is NaN or Inf, X is normal. */ \
663 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
664 _FP_FRAC_COPY_##wc (R, Y); \
665 goto sub_done; \
668 /* Insert implicit MSB of X. */ \
669 _FP_FRAC_HIGH_##fs (X) |= _FP_IMPLBIT_SH_##fs; \
671 sub2: \
672 /* Shift the mantissa of X to the right EDIFF steps; \
673 remember to account later for the implicit MSB of Y. */ \
674 if (ediff <= _FP_WFRACBITS_##fs) \
675 _FP_FRAC_SRS_##wc (X, ediff, _FP_WFRACBITS_##fs); \
676 else if (!_FP_FRAC_ZEROP_##wc (X)) \
677 _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc); \
678 _FP_FRAC_SUB_##wc (R, Y, X); \
680 else \
682 /* ediff == 0. */ \
683 if (!_FP_EXP_NORMAL (fs, wc, X)) \
685 if (X##_e == 0) \
687 /* X and Y are zero or denormalized. */ \
688 R##_e = 0; \
689 if (_FP_FRAC_ZEROP_##wc (X)) \
691 _FP_FRAC_COPY_##wc (R, Y); \
692 if (_FP_FRAC_ZEROP_##wc (Y)) \
693 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
694 else \
696 FP_SET_EXCEPTION (FP_EX_DENORM); \
697 R##_s = Y##_s; \
699 goto sub_done; \
701 else if (_FP_FRAC_ZEROP_##wc (Y)) \
703 FP_SET_EXCEPTION (FP_EX_DENORM); \
704 _FP_FRAC_COPY_##wc (R, X); \
705 R##_s = X##_s; \
706 goto sub_done; \
708 else \
710 FP_SET_EXCEPTION (FP_EX_DENORM); \
711 _FP_FRAC_SUB_##wc (R, X, Y); \
712 R##_s = X##_s; \
713 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
715 /* |X| < |Y|, negate result. */ \
716 _FP_FRAC_SUB_##wc (R, Y, X); \
717 R##_s = Y##_s; \
719 else if (_FP_FRAC_ZEROP_##wc (R)) \
720 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
721 goto sub_done; \
724 else \
726 /* X and Y are NaN or Inf, of opposite signs. */ \
727 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X); \
728 _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y); \
729 R##_e = _FP_EXPMAX_##fs; \
730 if (_FP_FRAC_ZEROP_##wc (X)) \
732 if (_FP_FRAC_ZEROP_##wc (Y)) \
734 /* Inf - Inf. */ \
735 R##_s = _FP_NANSIGN_##fs; \
736 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
737 _FP_FRAC_SLL_##wc (R, _FP_WORKBITS); \
738 FP_SET_EXCEPTION (FP_EX_INVALID); \
740 else \
742 /* Inf - NaN. */ \
743 R##_s = Y##_s; \
744 _FP_FRAC_COPY_##wc (R, Y); \
747 else \
749 if (_FP_FRAC_ZEROP_##wc (Y)) \
751 /* NaN - Inf. */ \
752 R##_s = X##_s; \
753 _FP_FRAC_COPY_##wc (R, X); \
755 else \
757 /* NaN - NaN. */ \
758 _FP_CHOOSENAN_SEMIRAW (fs, wc, R, X, Y, OP); \
761 goto sub_done; \
764 /* The exponents of X and Y, both normal, are equal. The \
765 implicit MSBs cancel. */ \
766 R##_e = X##_e; \
767 _FP_FRAC_SUB_##wc (R, X, Y); \
768 R##_s = X##_s; \
769 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
771 /* |X| < |Y|, negate result. */ \
772 _FP_FRAC_SUB_##wc (R, Y, X); \
773 R##_s = Y##_s; \
775 else if (_FP_FRAC_ZEROP_##wc (R)) \
777 R##_e = 0; \
778 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
779 goto sub_done; \
781 goto norm; \
783 sub3: \
784 if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
786 int diff; \
787 /* Carry into most significant bit of larger one of X and Y, \
788 canceling it; renormalize. */ \
789 _FP_FRAC_HIGH_##fs (R) &= _FP_IMPLBIT_SH_##fs - 1; \
790 norm: \
791 _FP_FRAC_CLZ_##wc (diff, R); \
792 diff -= _FP_WFRACXBITS_##fs; \
793 _FP_FRAC_SLL_##wc (R, diff); \
794 if (R##_e <= diff) \
796 /* R is denormalized. */ \
797 diff = diff - R##_e + 1; \
798 _FP_FRAC_SRS_##wc (R, diff, _FP_WFRACBITS_##fs); \
799 R##_e = 0; \
801 else \
803 R##_e -= diff; \
804 _FP_FRAC_HIGH_##fs (R) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
807 sub_done: ; \
810 while (0)
812 #define _FP_ADD(fs, wc, R, X, Y) _FP_ADD_INTERNAL (fs, wc, R, X, Y, '+')
813 #define _FP_SUB(fs, wc, R, X, Y) \
814 do \
816 if (!(Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))) \
817 Y##_s ^= 1; \
818 _FP_ADD_INTERNAL (fs, wc, R, X, Y, '-'); \
820 while (0)
824 * Main negation routine. The input value is raw.
827 #define _FP_NEG(fs, wc, R, X) \
828 do \
830 _FP_FRAC_COPY_##wc (R, X); \
831 R##_e = X##_e; \
832 R##_s = 1 ^ X##_s; \
834 while (0)
838 * Main multiplication routine. The input values should be cooked.
841 #define _FP_MUL(fs, wc, R, X, Y) \
842 do \
844 R##_s = X##_s ^ Y##_s; \
845 R##_e = X##_e + Y##_e + 1; \
846 switch (_FP_CLS_COMBINE (X##_c, Y##_c)) \
848 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL): \
849 R##_c = FP_CLS_NORMAL; \
851 _FP_MUL_MEAT_##fs (R, X, Y); \
853 if (_FP_FRAC_OVERP_##wc (fs, R)) \
854 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs); \
855 else \
856 R##_e--; \
857 break; \
859 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN): \
860 _FP_CHOOSENAN (fs, wc, R, X, Y, '*'); \
861 break; \
863 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL): \
864 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF): \
865 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO): \
866 R##_s = X##_s; \
868 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF): \
869 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL): \
870 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL): \
871 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO): \
872 _FP_FRAC_COPY_##wc (R, X); \
873 R##_c = X##_c; \
874 break; \
876 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN): \
877 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN): \
878 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN): \
879 R##_s = Y##_s; \
881 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF): \
882 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO): \
883 _FP_FRAC_COPY_##wc (R, Y); \
884 R##_c = Y##_c; \
885 break; \
887 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO): \
888 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF): \
889 R##_s = _FP_NANSIGN_##fs; \
890 R##_c = FP_CLS_NAN; \
891 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
892 FP_SET_EXCEPTION (FP_EX_INVALID); \
893 break; \
895 default: \
896 abort (); \
899 while (0)
902 /* Fused multiply-add. The input values should be cooked. */
904 #define _FP_FMA(fs, wc, dwc, R, X, Y, Z) \
905 do \
907 FP_DECL_##fs (T); \
908 T##_s = X##_s ^ Y##_s; \
909 T##_e = X##_e + Y##_e + 1; \
910 switch (_FP_CLS_COMBINE (X##_c, Y##_c)) \
912 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL): \
913 switch (Z##_c) \
915 case FP_CLS_INF: \
916 case FP_CLS_NAN: \
917 R##_s = Z##_s; \
918 _FP_FRAC_COPY_##wc (R, Z); \
919 R##_c = Z##_c; \
920 break; \
922 case FP_CLS_ZERO: \
923 R##_c = FP_CLS_NORMAL; \
924 R##_s = T##_s; \
925 R##_e = T##_e; \
927 _FP_MUL_MEAT_##fs (R, X, Y); \
929 if (_FP_FRAC_OVERP_##wc (fs, R)) \
930 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs); \
931 else \
932 R##_e--; \
933 break; \
935 case FP_CLS_NORMAL:; \
936 _FP_FRAC_DECL_##dwc (TD); \
937 _FP_FRAC_DECL_##dwc (ZD); \
938 _FP_FRAC_DECL_##dwc (RD); \
939 _FP_MUL_MEAT_DW_##fs (TD, X, Y); \
940 R##_e = T##_e; \
941 int tsh = _FP_FRAC_HIGHBIT_DW_##dwc (fs, TD) == 0; \
942 T##_e -= tsh; \
943 int ediff = T##_e - Z##_e; \
944 if (ediff >= 0) \
946 int shift = _FP_WFRACBITS_##fs - tsh - ediff; \
947 if (shift <= -_FP_WFRACBITS_##fs) \
948 _FP_FRAC_SET_##dwc (ZD, _FP_MINFRAC_##dwc); \
949 else \
951 _FP_FRAC_COPY_##dwc##_##wc (ZD, Z); \
952 if (shift < 0) \
953 _FP_FRAC_SRS_##dwc (ZD, -shift, \
954 _FP_WFRACBITS_DW_##fs); \
955 else if (shift > 0) \
956 _FP_FRAC_SLL_##dwc (ZD, shift); \
958 R##_s = T##_s; \
959 if (T##_s == Z##_s) \
960 _FP_FRAC_ADD_##dwc (RD, TD, ZD); \
961 else \
963 _FP_FRAC_SUB_##dwc (RD, TD, ZD); \
964 if (_FP_FRAC_NEGP_##dwc (RD)) \
966 R##_s = Z##_s; \
967 _FP_FRAC_SUB_##dwc (RD, ZD, TD); \
971 else \
973 R##_e = Z##_e; \
974 R##_s = Z##_s; \
975 _FP_FRAC_COPY_##dwc##_##wc (ZD, Z); \
976 _FP_FRAC_SLL_##dwc (ZD, _FP_WFRACBITS_##fs); \
977 int shift = -ediff - tsh; \
978 if (shift >= _FP_WFRACBITS_DW_##fs) \
979 _FP_FRAC_SET_##dwc (TD, _FP_MINFRAC_##dwc); \
980 else if (shift > 0) \
981 _FP_FRAC_SRS_##dwc (TD, shift, \
982 _FP_WFRACBITS_DW_##fs); \
983 if (Z##_s == T##_s) \
984 _FP_FRAC_ADD_##dwc (RD, ZD, TD); \
985 else \
986 _FP_FRAC_SUB_##dwc (RD, ZD, TD); \
988 if (_FP_FRAC_ZEROP_##dwc (RD)) \
990 if (T##_s == Z##_s) \
991 R##_s = Z##_s; \
992 else \
993 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
994 _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc); \
995 R##_c = FP_CLS_ZERO; \
997 else \
999 int rlz; \
1000 _FP_FRAC_CLZ_##dwc (rlz, RD); \
1001 rlz -= _FP_WFRACXBITS_DW_##fs; \
1002 R##_e -= rlz; \
1003 int shift = _FP_WFRACBITS_##fs - rlz; \
1004 if (shift > 0) \
1005 _FP_FRAC_SRS_##dwc (RD, shift, \
1006 _FP_WFRACBITS_DW_##fs); \
1007 else if (shift < 0) \
1008 _FP_FRAC_SLL_##dwc (RD, -shift); \
1009 _FP_FRAC_COPY_##wc##_##dwc (R, RD); \
1010 R##_c = FP_CLS_NORMAL; \
1012 break; \
1014 goto done_fma; \
1016 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN): \
1017 _FP_CHOOSENAN (fs, wc, T, X, Y, '*'); \
1018 break; \
1020 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL): \
1021 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF): \
1022 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO): \
1023 T##_s = X##_s; \
1025 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF): \
1026 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL): \
1027 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL): \
1028 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO): \
1029 _FP_FRAC_COPY_##wc (T, X); \
1030 T##_c = X##_c; \
1031 break; \
1033 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN): \
1034 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN): \
1035 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN): \
1036 T##_s = Y##_s; \
1038 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF): \
1039 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO): \
1040 _FP_FRAC_COPY_##wc (T, Y); \
1041 T##_c = Y##_c; \
1042 break; \
1044 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO): \
1045 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF): \
1046 T##_s = _FP_NANSIGN_##fs; \
1047 T##_c = FP_CLS_NAN; \
1048 _FP_FRAC_SET_##wc (T, _FP_NANFRAC_##fs); \
1049 FP_SET_EXCEPTION (FP_EX_INVALID); \
1050 break; \
1052 default: \
1053 abort (); \
1056 /* T = X * Y is zero, infinity or NaN. */ \
1057 switch (_FP_CLS_COMBINE (T##_c, Z##_c)) \
1059 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN): \
1060 _FP_CHOOSENAN (fs, wc, R, T, Z, '+'); \
1061 break; \
1063 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL): \
1064 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF): \
1065 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO): \
1066 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL): \
1067 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO): \
1068 R##_s = T##_s; \
1069 _FP_FRAC_COPY_##wc (R, T); \
1070 R##_c = T##_c; \
1071 break; \
1073 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN): \
1074 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN): \
1075 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL): \
1076 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF): \
1077 R##_s = Z##_s; \
1078 _FP_FRAC_COPY_##wc (R, Z); \
1079 R##_c = Z##_c; \
1080 break; \
1082 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF): \
1083 if (T##_s == Z##_s) \
1085 R##_s = Z##_s; \
1086 _FP_FRAC_COPY_##wc (R, Z); \
1087 R##_c = Z##_c; \
1089 else \
1091 R##_s = _FP_NANSIGN_##fs; \
1092 R##_c = FP_CLS_NAN; \
1093 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
1094 FP_SET_EXCEPTION (FP_EX_INVALID); \
1096 break; \
1098 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO): \
1099 if (T##_s == Z##_s) \
1100 R##_s = Z##_s; \
1101 else \
1102 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
1103 _FP_FRAC_COPY_##wc (R, Z); \
1104 R##_c = Z##_c; \
1105 break; \
1107 default: \
1108 abort (); \
1110 done_fma: ; \
1112 while (0)
1116 * Main division routine. The input values should be cooked.
1119 #define _FP_DIV(fs, wc, R, X, Y) \
1120 do \
1122 R##_s = X##_s ^ Y##_s; \
1123 R##_e = X##_e - Y##_e; \
1124 switch (_FP_CLS_COMBINE (X##_c, Y##_c)) \
1126 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL): \
1127 R##_c = FP_CLS_NORMAL; \
1129 _FP_DIV_MEAT_##fs (R, X, Y); \
1130 break; \
1132 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN): \
1133 _FP_CHOOSENAN (fs, wc, R, X, Y, '/'); \
1134 break; \
1136 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL): \
1137 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF): \
1138 case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO): \
1139 R##_s = X##_s; \
1140 _FP_FRAC_COPY_##wc (R, X); \
1141 R##_c = X##_c; \
1142 break; \
1144 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN): \
1145 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN): \
1146 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN): \
1147 R##_s = Y##_s; \
1148 _FP_FRAC_COPY_##wc (R, Y); \
1149 R##_c = Y##_c; \
1150 break; \
1152 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF): \
1153 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF): \
1154 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL): \
1155 R##_c = FP_CLS_ZERO; \
1156 break; \
1158 case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO): \
1159 FP_SET_EXCEPTION (FP_EX_DIVZERO); \
1160 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO): \
1161 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL): \
1162 R##_c = FP_CLS_INF; \
1163 break; \
1165 case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF): \
1166 case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO): \
1167 R##_s = _FP_NANSIGN_##fs; \
1168 R##_c = FP_CLS_NAN; \
1169 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
1170 FP_SET_EXCEPTION (FP_EX_INVALID); \
1171 break; \
1173 default: \
1174 abort (); \
1177 while (0)
1181 * Main differential comparison routine. The inputs should be raw not
1182 * cooked. The return is -1,0,1 for normal values, 2 otherwise.
1185 #define _FP_CMP(fs, wc, ret, X, Y, un) \
1186 do \
1188 /* NANs are unordered */ \
1189 if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
1190 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))) \
1192 ret = un; \
1194 else \
1196 int __is_zero_x; \
1197 int __is_zero_y; \
1199 __is_zero_x = (!X##_e && _FP_FRAC_ZEROP_##wc (X)) ? 1 : 0; \
1200 __is_zero_y = (!Y##_e && _FP_FRAC_ZEROP_##wc (Y)) ? 1 : 0; \
1202 if (__is_zero_x && __is_zero_y) \
1203 ret = 0; \
1204 else if (__is_zero_x) \
1205 ret = Y##_s ? 1 : -1; \
1206 else if (__is_zero_y) \
1207 ret = X##_s ? -1 : 1; \
1208 else if (X##_s != Y##_s) \
1209 ret = X##_s ? -1 : 1; \
1210 else if (X##_e > Y##_e) \
1211 ret = X##_s ? -1 : 1; \
1212 else if (X##_e < Y##_e) \
1213 ret = X##_s ? 1 : -1; \
1214 else if (_FP_FRAC_GT_##wc (X, Y)) \
1215 ret = X##_s ? -1 : 1; \
1216 else if (_FP_FRAC_GT_##wc (Y, X)) \
1217 ret = X##_s ? 1 : -1; \
1218 else \
1219 ret = 0; \
1222 while (0)
1225 /* Simplification for strict equality. */
1227 #define _FP_CMP_EQ(fs, wc, ret, X, Y) \
1228 do \
1230 /* NANs are unordered */ \
1231 if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
1232 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))) \
1234 ret = 1; \
1236 else \
1238 ret = !(X##_e == Y##_e \
1239 && _FP_FRAC_EQ_##wc (X, Y) \
1240 && (X##_s == Y##_s || (!X##_e && _FP_FRAC_ZEROP_##wc (X)))); \
1243 while (0)
1245 /* Version to test unordered. */
1247 #define _FP_CMP_UNORD(fs, wc, ret, X, Y) \
1248 do \
1250 ret = ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
1251 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))); \
1253 while (0)
1256 * Main square root routine. The input value should be cooked.
1259 #define _FP_SQRT(fs, wc, R, X) \
1260 do \
1262 _FP_FRAC_DECL_##wc (T); \
1263 _FP_FRAC_DECL_##wc (S); \
1264 _FP_W_TYPE q; \
1265 switch (X##_c) \
1267 case FP_CLS_NAN: \
1268 _FP_FRAC_COPY_##wc (R, X); \
1269 R##_s = X##_s; \
1270 R##_c = FP_CLS_NAN; \
1271 break; \
1272 case FP_CLS_INF: \
1273 if (X##_s) \
1275 R##_s = _FP_NANSIGN_##fs; \
1276 R##_c = FP_CLS_NAN; /* NAN */ \
1277 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
1278 FP_SET_EXCEPTION (FP_EX_INVALID); \
1280 else \
1282 R##_s = 0; \
1283 R##_c = FP_CLS_INF; /* sqrt(+inf) = +inf */ \
1285 break; \
1286 case FP_CLS_ZERO: \
1287 R##_s = X##_s; \
1288 R##_c = FP_CLS_ZERO; /* sqrt(+-0) = +-0 */ \
1289 break; \
1290 case FP_CLS_NORMAL: \
1291 R##_s = 0; \
1292 if (X##_s) \
1294 R##_c = FP_CLS_NAN; /* NAN */ \
1295 R##_s = _FP_NANSIGN_##fs; \
1296 _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs); \
1297 FP_SET_EXCEPTION (FP_EX_INVALID); \
1298 break; \
1300 R##_c = FP_CLS_NORMAL; \
1301 if (X##_e & 1) \
1302 _FP_FRAC_SLL_##wc (X, 1); \
1303 R##_e = X##_e >> 1; \
1304 _FP_FRAC_SET_##wc (S, _FP_ZEROFRAC_##wc); \
1305 _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc); \
1306 q = _FP_OVERFLOW_##fs >> 1; \
1307 _FP_SQRT_MEAT_##wc (R, S, T, X, q); \
1310 while (0)
1313 * Convert from FP to integer. Input is raw.
1316 /* RSIGNED can have following values:
1317 * 0: the number is required to be 0..(2^rsize)-1, if not, NV is set plus
1318 * the result is either 0 or (2^rsize)-1 depending on the sign in such
1319 * case.
1320 * 1: the number is required to be -(2^(rsize-1))..(2^(rsize-1))-1, if not,
1321 * NV is set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
1322 * depending on the sign in such case.
1323 * -1: the number is required to be -(2^(rsize-1))..(2^rsize)-1, if not, NV is
1324 * set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
1325 * depending on the sign in such case.
1327 #define _FP_TO_INT(fs, wc, r, X, rsize, rsigned) \
1328 do \
1330 if (X##_e < _FP_EXPBIAS_##fs) \
1332 r = 0; \
1333 if (X##_e == 0) \
1335 if (!_FP_FRAC_ZEROP_##wc (X)) \
1337 FP_SET_EXCEPTION (FP_EX_INEXACT); \
1338 FP_SET_EXCEPTION (FP_EX_DENORM); \
1341 else \
1342 FP_SET_EXCEPTION (FP_EX_INEXACT); \
1344 else if (X##_e >= _FP_EXPBIAS_##fs + rsize - (rsigned > 0 || X##_s) \
1345 || (!rsigned && X##_s)) \
1347 /* Overflow or converting to the most negative integer. */ \
1348 if (rsigned) \
1350 r = 1; \
1351 r <<= rsize - 1; \
1352 r -= 1 - X##_s; \
1353 } else { \
1354 r = 0; \
1355 if (!X##_s) \
1356 r = ~r; \
1359 if (rsigned && X##_s && X##_e == _FP_EXPBIAS_##fs + rsize - 1) \
1361 /* Possibly converting to most negative integer; check the \
1362 mantissa. */ \
1363 int inexact = 0; \
1364 (void) ((_FP_FRACBITS_##fs > rsize) \
1365 ? ({ \
1366 _FP_FRAC_SRST_##wc (X, inexact, \
1367 _FP_FRACBITS_##fs - rsize, \
1368 _FP_FRACBITS_##fs); \
1369 0; \
1370 }) \
1371 : 0); \
1372 if (!_FP_FRAC_ZEROP_##wc (X)) \
1373 FP_SET_EXCEPTION (FP_EX_INVALID); \
1374 else if (inexact) \
1375 FP_SET_EXCEPTION (FP_EX_INEXACT); \
1377 else \
1378 FP_SET_EXCEPTION (FP_EX_INVALID); \
1380 else \
1382 _FP_FRAC_HIGH_RAW_##fs (X) |= _FP_IMPLBIT_##fs; \
1383 if (X##_e >= _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs - 1) \
1385 _FP_FRAC_ASSEMBLE_##wc (r, X, rsize); \
1386 r <<= X##_e - _FP_EXPBIAS_##fs - _FP_FRACBITS_##fs + 1; \
1388 else \
1390 int inexact; \
1391 _FP_FRAC_SRST_##wc (X, inexact, \
1392 (_FP_FRACBITS_##fs + _FP_EXPBIAS_##fs - 1 \
1393 - X##_e), \
1394 _FP_FRACBITS_##fs); \
1395 if (inexact) \
1396 FP_SET_EXCEPTION (FP_EX_INEXACT); \
1397 _FP_FRAC_ASSEMBLE_##wc (r, X, rsize); \
1399 if (rsigned && X##_s) \
1400 r = -r; \
1403 while (0)
1405 /* Convert integer to fp. Output is raw. RTYPE is unsigned even if
1406 input is signed. */
1407 #define _FP_FROM_INT(fs, wc, X, r, rsize, rtype) \
1408 do \
1410 if (r) \
1412 rtype ur_; \
1414 if ((X##_s = (r < 0))) \
1415 r = -(rtype) r; \
1417 ur_ = (rtype) r; \
1418 (void) ((rsize <= _FP_W_TYPE_SIZE) \
1419 ? ({ \
1420 int lz_; \
1421 __FP_CLZ (lz_, (_FP_W_TYPE) ur_); \
1422 X##_e = _FP_EXPBIAS_##fs + _FP_W_TYPE_SIZE - 1 - lz_; \
1423 }) \
1424 : ((rsize <= 2 * _FP_W_TYPE_SIZE) \
1425 ? ({ \
1426 int lz_; \
1427 __FP_CLZ_2 (lz_, \
1428 (_FP_W_TYPE) (ur_ >> _FP_W_TYPE_SIZE), \
1429 (_FP_W_TYPE) ur_); \
1430 X##_e = (_FP_EXPBIAS_##fs + 2 * _FP_W_TYPE_SIZE - 1 \
1431 - lz_); \
1432 }) \
1433 : (abort (), 0))); \
1435 if (rsize - 1 + _FP_EXPBIAS_##fs >= _FP_EXPMAX_##fs \
1436 && X##_e >= _FP_EXPMAX_##fs) \
1438 /* Exponent too big; overflow to infinity. (May also \
1439 happen after rounding below.) */ \
1440 _FP_OVERFLOW_SEMIRAW (fs, wc, X); \
1441 goto pack_semiraw; \
1444 if (rsize <= _FP_FRACBITS_##fs \
1445 || X##_e < _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs) \
1447 /* Exactly representable; shift left. */ \
1448 _FP_FRAC_DISASSEMBLE_##wc (X, ur_, rsize); \
1449 if (_FP_EXPBIAS_##fs + _FP_FRACBITS_##fs - 1 - X##_e > 0) \
1450 _FP_FRAC_SLL_##wc (X, (_FP_EXPBIAS_##fs \
1451 + _FP_FRACBITS_##fs - 1 - X##_e)); \
1453 else \
1455 /* More bits in integer than in floating type; need to \
1456 round. */ \
1457 if (_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 < X##_e) \
1458 ur_ = ((ur_ >> (X##_e - _FP_EXPBIAS_##fs \
1459 - _FP_WFRACBITS_##fs + 1)) \
1460 | ((ur_ << (rsize - (X##_e - _FP_EXPBIAS_##fs \
1461 - _FP_WFRACBITS_##fs + 1))) \
1462 != 0)); \
1463 _FP_FRAC_DISASSEMBLE_##wc (X, ur_, rsize); \
1464 if ((_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 - X##_e) > 0) \
1465 _FP_FRAC_SLL_##wc (X, (_FP_EXPBIAS_##fs \
1466 + _FP_WFRACBITS_##fs - 1 - X##_e)); \
1467 _FP_FRAC_HIGH_##fs (X) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
1468 pack_semiraw: \
1469 _FP_PACK_SEMIRAW (fs, wc, X); \
1472 else \
1474 X##_s = 0; \
1475 X##_e = 0; \
1476 _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc); \
1479 while (0)
1482 /* Extend from a narrower floating-point format to a wider one. Input
1483 and output are raw. */
1484 #define FP_EXTEND(dfs, sfs, dwc, swc, D, S) \
1485 do \
1487 if (_FP_FRACBITS_##dfs < _FP_FRACBITS_##sfs \
1488 || (_FP_EXPMAX_##dfs - _FP_EXPBIAS_##dfs \
1489 < _FP_EXPMAX_##sfs - _FP_EXPBIAS_##sfs) \
1490 || (_FP_EXPBIAS_##dfs < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1 \
1491 && _FP_EXPBIAS_##dfs != _FP_EXPBIAS_##sfs)) \
1492 abort (); \
1493 D##_s = S##_s; \
1494 _FP_FRAC_COPY_##dwc##_##swc (D, S); \
1495 if (_FP_EXP_NORMAL (sfs, swc, S)) \
1497 D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs; \
1498 _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs - _FP_FRACBITS_##sfs)); \
1500 else \
1502 if (S##_e == 0) \
1504 if (_FP_FRAC_ZEROP_##swc (S)) \
1505 D##_e = 0; \
1506 else if (_FP_EXPBIAS_##dfs \
1507 < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1) \
1509 FP_SET_EXCEPTION (FP_EX_DENORM); \
1510 _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs \
1511 - _FP_FRACBITS_##sfs)); \
1512 D##_e = 0; \
1514 else \
1516 int _lz; \
1517 FP_SET_EXCEPTION (FP_EX_DENORM); \
1518 _FP_FRAC_CLZ_##swc (_lz, S); \
1519 _FP_FRAC_SLL_##dwc (D, \
1520 _lz + _FP_FRACBITS_##dfs \
1521 - _FP_FRACTBITS_##sfs); \
1522 D##_e = (_FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs + 1 \
1523 + _FP_FRACXBITS_##sfs - _lz); \
1526 else \
1528 D##_e = _FP_EXPMAX_##dfs; \
1529 if (!_FP_FRAC_ZEROP_##swc (S)) \
1531 if (_FP_FRAC_SNANP (sfs, S)) \
1532 FP_SET_EXCEPTION (FP_EX_INVALID); \
1533 _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs \
1534 - _FP_FRACBITS_##sfs)); \
1535 _FP_SETQNAN (dfs, dwc, D); \
1540 while (0)
1542 /* Truncate from a wider floating-point format to a narrower one.
1543 Input and output are semi-raw. */
1544 #define FP_TRUNC(dfs, sfs, dwc, swc, D, S) \
1545 do \
1547 if (_FP_FRACBITS_##sfs < _FP_FRACBITS_##dfs \
1548 || (_FP_EXPBIAS_##sfs < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1 \
1549 && _FP_EXPBIAS_##sfs != _FP_EXPBIAS_##dfs)) \
1550 abort (); \
1551 D##_s = S##_s; \
1552 if (_FP_EXP_NORMAL (sfs, swc, S)) \
1554 D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs; \
1555 if (D##_e >= _FP_EXPMAX_##dfs) \
1556 _FP_OVERFLOW_SEMIRAW (dfs, dwc, D); \
1557 else \
1559 if (D##_e <= 0) \
1561 if (D##_e < 1 - _FP_FRACBITS_##dfs) \
1563 _FP_FRAC_SET_##swc (S, _FP_ZEROFRAC_##swc); \
1564 _FP_FRAC_LOW_##swc (S) |= 1; \
1566 else \
1568 _FP_FRAC_HIGH_##sfs (S) |= _FP_IMPLBIT_SH_##sfs; \
1569 _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs \
1570 - _FP_WFRACBITS_##dfs \
1571 + 1 - D##_e), \
1572 _FP_WFRACBITS_##sfs); \
1574 D##_e = 0; \
1576 else \
1577 _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs \
1578 - _FP_WFRACBITS_##dfs), \
1579 _FP_WFRACBITS_##sfs); \
1580 _FP_FRAC_COPY_##dwc##_##swc (D, S); \
1583 else \
1585 if (S##_e == 0) \
1587 D##_e = 0; \
1588 if (_FP_FRAC_ZEROP_##swc (S)) \
1589 _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc); \
1590 else \
1592 FP_SET_EXCEPTION (FP_EX_DENORM); \
1593 if (_FP_EXPBIAS_##sfs \
1594 < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1) \
1596 _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs \
1597 - _FP_WFRACBITS_##dfs), \
1598 _FP_WFRACBITS_##sfs); \
1599 _FP_FRAC_COPY_##dwc##_##swc (D, S); \
1601 else \
1603 _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc); \
1604 _FP_FRAC_LOW_##dwc (D) |= 1; \
1608 else \
1610 D##_e = _FP_EXPMAX_##dfs; \
1611 if (_FP_FRAC_ZEROP_##swc (S)) \
1612 _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc); \
1613 else \
1615 _FP_CHECK_SIGNAN_SEMIRAW (sfs, swc, S); \
1616 _FP_FRAC_SRL_##swc (S, (_FP_WFRACBITS_##sfs \
1617 - _FP_WFRACBITS_##dfs)); \
1618 _FP_FRAC_COPY_##dwc##_##swc (D, S); \
1619 /* Semi-raw NaN must have all workbits cleared. */ \
1620 _FP_FRAC_LOW_##dwc (D) \
1621 &= ~(_FP_W_TYPE) ((1 << _FP_WORKBITS) - 1); \
1622 _FP_SETQNAN_SEMIRAW (dfs, dwc, D); \
1627 while (0)
1630 * Helper primitives.
1633 /* Count leading zeros in a word. */
1635 #ifndef __FP_CLZ
1636 /* GCC 3.4 and later provide the builtins for us. */
1637 # define __FP_CLZ(r, x) \
1638 do \
1640 if (sizeof (_FP_W_TYPE) == sizeof (unsigned int)) \
1641 r = __builtin_clz (x); \
1642 else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long)) \
1643 r = __builtin_clzl (x); \
1644 else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long long)) \
1645 r = __builtin_clzll (x); \
1646 else \
1647 abort (); \
1649 while (0)
1650 #endif /* ndef __FP_CLZ */
1652 #define _FP_DIV_HELP_imm(q, r, n, d) \
1653 do \
1655 q = n / d, r = n % d; \
1657 while (0)
1660 /* A restoring bit-by-bit division primitive. */
1662 #define _FP_DIV_MEAT_N_loop(fs, wc, R, X, Y) \
1663 do \
1665 int count = _FP_WFRACBITS_##fs; \
1666 _FP_FRAC_DECL_##wc (u); \
1667 _FP_FRAC_DECL_##wc (v); \
1668 _FP_FRAC_COPY_##wc (u, X); \
1669 _FP_FRAC_COPY_##wc (v, Y); \
1670 _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc); \
1671 /* Normalize U and V. */ \
1672 _FP_FRAC_SLL_##wc (u, _FP_WFRACXBITS_##fs); \
1673 _FP_FRAC_SLL_##wc (v, _FP_WFRACXBITS_##fs); \
1674 /* First round. Since the operands are normalized, either the \
1675 first or second bit will be set in the fraction. Produce a \
1676 normalized result by checking which and adjusting the loop \
1677 count and exponent accordingly. */ \
1678 if (_FP_FRAC_GE_1 (u, v)) \
1680 _FP_FRAC_SUB_##wc (u, u, v); \
1681 _FP_FRAC_LOW_##wc (R) |= 1; \
1682 count--; \
1684 else \
1685 R##_e--; \
1686 /* Subsequent rounds. */ \
1687 do \
1689 int msb = (_FP_WS_TYPE) _FP_FRAC_HIGH_##wc (u) < 0; \
1690 _FP_FRAC_SLL_##wc (u, 1); \
1691 _FP_FRAC_SLL_##wc (R, 1); \
1692 if (msb || _FP_FRAC_GE_1 (u, v)) \
1694 _FP_FRAC_SUB_##wc (u, u, v); \
1695 _FP_FRAC_LOW_##wc (R) |= 1; \
1698 while (--count > 0); \
1699 /* If there's anything left in U, the result is inexact. */ \
1700 _FP_FRAC_LOW_##wc (R) |= !_FP_FRAC_ZEROP_##wc (u); \
1702 while (0)
1704 #define _FP_DIV_MEAT_1_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 1, R, X, Y)
1705 #define _FP_DIV_MEAT_2_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 2, R, X, Y)
1706 #define _FP_DIV_MEAT_4_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 4, R, X, Y)