soft-fp: Correct _FP_TO_INT formatting.
[platform/upstream/glibc.git] / soft-fp / op-common.h
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).
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, see
30    <http://www.gnu.org/licenses/>.  */
31
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)
37
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))
47
48 /*
49  * Finish truly unpacking a native fp value by classifying the kind
50  * of fp value and normalizing both the exponent and the fraction.
51  */
52
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;                                                \
64                                                                 \
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 _FP_UNPACK_CANONICAL_shift;            \
72               _FP_FRAC_CLZ_##wc (_FP_UNPACK_CANONICAL_shift,    \
73                                  X);                            \
74               _FP_UNPACK_CANONICAL_shift -= _FP_FRACXBITS_##fs; \
75               _FP_FRAC_SLL_##wc (X, (_FP_UNPACK_CANONICAL_shift \
76                                      + _FP_WORKBITS));          \
77               X##_e -= (_FP_EXPBIAS_##fs - 1                    \
78                         + _FP_UNPACK_CANONICAL_shift);          \
79               X##_c = FP_CLS_NORMAL;                            \
80               FP_SET_EXCEPTION (FP_EX_DENORM);                  \
81             }                                                   \
82           break;                                                \
83                                                                 \
84         case _FP_EXPMAX_##fs:                                   \
85           if (_FP_FRAC_ZEROP_##wc (X))                          \
86             X##_c = FP_CLS_INF;                                 \
87           else                                                  \
88             {                                                   \
89               X##_c = FP_CLS_NAN;                               \
90               /* Check for signaling NaN */                     \
91               if (_FP_FRAC_SNANP (fs, X))                       \
92                 FP_SET_EXCEPTION (FP_EX_INVALID);               \
93             }                                                   \
94           break;                                                \
95         }                                                       \
96     }                                                           \
97   while (0)
98
99 /* Finish unpacking an fp value in semi-raw mode: the mantissa is
100    shifted by _FP_WORKBITS but the implicit MSB is not inserted and
101    other classification is not done.  */
102 #define _FP_UNPACK_SEMIRAW(fs, wc, X)   _FP_FRAC_SLL_##wc (X, _FP_WORKBITS)
103
104 /* A semi-raw value has overflowed to infinity.  Adjust the mantissa
105    and exponent appropriately.  */
106 #define _FP_OVERFLOW_SEMIRAW(fs, wc, X)                 \
107   do                                                    \
108     {                                                   \
109       if (FP_ROUNDMODE == FP_RND_NEAREST                \
110           || (FP_ROUNDMODE == FP_RND_PINF && !X##_s)    \
111           || (FP_ROUNDMODE == FP_RND_MINF && X##_s))    \
112         {                                               \
113           X##_e = _FP_EXPMAX_##fs;                      \
114           _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc);     \
115         }                                               \
116       else                                              \
117         {                                               \
118           X##_e = _FP_EXPMAX_##fs - 1;                  \
119           _FP_FRAC_SET_##wc (X, _FP_MAXFRAC_##wc);      \
120         }                                               \
121       FP_SET_EXCEPTION (FP_EX_INEXACT);                 \
122       FP_SET_EXCEPTION (FP_EX_OVERFLOW);                \
123     }                                                   \
124   while (0)
125
126 /* Check for a semi-raw value being a signaling NaN and raise the
127    invalid exception if so.  */
128 #define _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X)     \
129   do                                            \
130     {                                           \
131       if (X##_e == _FP_EXPMAX_##fs              \
132           && !_FP_FRAC_ZEROP_##wc (X)           \
133           && _FP_FRAC_SNANP_SEMIRAW (fs, X))    \
134         FP_SET_EXCEPTION (FP_EX_INVALID);       \
135     }                                           \
136   while (0)
137
138 /* Choose a NaN result from an operation on two semi-raw NaN
139    values.  */
140 #define _FP_CHOOSENAN_SEMIRAW(fs, wc, R, X, Y, OP)                      \
141   do                                                                    \
142     {                                                                   \
143       /* _FP_CHOOSENAN expects raw values, so shift as required.  */    \
144       _FP_FRAC_SRL_##wc (X, _FP_WORKBITS);                              \
145       _FP_FRAC_SRL_##wc (Y, _FP_WORKBITS);                              \
146       _FP_CHOOSENAN (fs, wc, R, X, Y, OP);                              \
147       _FP_FRAC_SLL_##wc (R, _FP_WORKBITS);                              \
148     }                                                                   \
149   while (0)
150
151 /* Make the fractional part a quiet NaN, preserving the payload
152    if possible, otherwise make it the canonical quiet NaN and set
153    the sign bit accordingly.  */
154 #define _FP_SETQNAN(fs, wc, X)                                  \
155   do                                                            \
156     {                                                           \
157       if (_FP_QNANNEGATEDP)                                     \
158         {                                                       \
159           _FP_FRAC_HIGH_RAW_##fs (X) &= _FP_QNANBIT_##fs - 1;   \
160           if (_FP_FRAC_ZEROP_##wc (X))                          \
161             {                                                   \
162               X##_s = _FP_NANSIGN_##fs;                         \
163               _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs);          \
164             }                                                   \
165         }                                                       \
166       else                                                      \
167         _FP_FRAC_HIGH_RAW_##fs (X) |= _FP_QNANBIT_##fs;         \
168     }                                                           \
169   while (0)
170 #define _FP_SETQNAN_SEMIRAW(fs, wc, X)                          \
171   do                                                            \
172     {                                                           \
173       if (_FP_QNANNEGATEDP)                                     \
174         {                                                       \
175           _FP_FRAC_HIGH_##fs (X) &= _FP_QNANBIT_SH_##fs - 1;    \
176           if (_FP_FRAC_ZEROP_##wc (X))                          \
177             {                                                   \
178               X##_s = _FP_NANSIGN_##fs;                         \
179               _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs);          \
180               _FP_FRAC_SLL_##wc (X, _FP_WORKBITS);              \
181             }                                                   \
182         }                                                       \
183       else                                                      \
184         _FP_FRAC_HIGH_##fs (X) |= _FP_QNANBIT_SH_##fs;          \
185     }                                                           \
186   while (0)
187
188 /* Test whether a biased exponent is normal (not zero or maximum).  */
189 #define _FP_EXP_NORMAL(fs, wc, X)       (((X##_e + 1) & _FP_EXPMAX_##fs) > 1)
190
191 /* Prepare to pack an fp value in semi-raw mode: the mantissa is
192    rounded and shifted right, with the rounding possibly increasing
193    the exponent (including changing a finite value to infinity).  */
194 #define _FP_PACK_SEMIRAW(fs, wc, X)                             \
195   do                                                            \
196     {                                                           \
197       int _FP_PACK_SEMIRAW_is_tiny                              \
198         = X##_e == 0 && !_FP_FRAC_ZEROP_##wc (X);               \
199       if (_FP_TININESS_AFTER_ROUNDING                           \
200           && _FP_PACK_SEMIRAW_is_tiny)                          \
201         {                                                       \
202           FP_DECL_##fs (_FP_PACK_SEMIRAW_T);                    \
203           _FP_FRAC_COPY_##wc (_FP_PACK_SEMIRAW_T, X);           \
204           _FP_PACK_SEMIRAW_T##_s = X##_s;                       \
205           _FP_PACK_SEMIRAW_T##_e = X##_e;                       \
206           _FP_FRAC_SLL_##wc (_FP_PACK_SEMIRAW_T, 1);            \
207           _FP_ROUND (wc, _FP_PACK_SEMIRAW_T);                   \
208           if (_FP_FRAC_OVERP_##wc (fs, _FP_PACK_SEMIRAW_T))     \
209             _FP_PACK_SEMIRAW_is_tiny = 0;                       \
210         }                                                       \
211       _FP_ROUND (wc, X);                                        \
212       if (_FP_PACK_SEMIRAW_is_tiny)                             \
213         {                                                       \
214           if ((FP_CUR_EXCEPTIONS & FP_EX_INEXACT)               \
215               || (FP_TRAPPING_EXCEPTIONS & FP_EX_UNDERFLOW))    \
216             FP_SET_EXCEPTION (FP_EX_UNDERFLOW);                 \
217         }                                                       \
218       if (_FP_FRAC_HIGH_##fs (X)                                \
219           & (_FP_OVERFLOW_##fs >> 1))                           \
220         {                                                       \
221           _FP_FRAC_HIGH_##fs (X) &= ~(_FP_OVERFLOW_##fs >> 1);  \
222           X##_e++;                                              \
223           if (X##_e == _FP_EXPMAX_##fs)                         \
224             _FP_OVERFLOW_SEMIRAW (fs, wc, X);                   \
225         }                                                       \
226       _FP_FRAC_SRL_##wc (X, _FP_WORKBITS);                      \
227       if (X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X)) \
228         {                                                       \
229           if (!_FP_KEEPNANFRACP)                                \
230             {                                                   \
231               _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs);          \
232               X##_s = _FP_NANSIGN_##fs;                         \
233             }                                                   \
234           else                                                  \
235             _FP_SETQNAN (fs, wc, X);                            \
236         }                                                       \
237     }                                                           \
238   while (0)
239
240 /*
241  * Before packing the bits back into the native fp result, take care
242  * of such mundane things as rounding and overflow.  Also, for some
243  * kinds of fp values, the original parts may not have been fully
244  * extracted -- but that is ok, we can regenerate them now.
245  */
246
247 #define _FP_PACK_CANONICAL(fs, wc, X)                                   \
248   do                                                                    \
249     {                                                                   \
250       switch (X##_c)                                                    \
251         {                                                               \
252         case FP_CLS_NORMAL:                                             \
253           X##_e += _FP_EXPBIAS_##fs;                                    \
254           if (X##_e > 0)                                                \
255             {                                                           \
256               _FP_ROUND (wc, X);                                        \
257               if (_FP_FRAC_OVERP_##wc (fs, X))                          \
258                 {                                                       \
259                   _FP_FRAC_CLEAR_OVERP_##wc (fs, X);                    \
260                   X##_e++;                                              \
261                 }                                                       \
262               _FP_FRAC_SRL_##wc (X, _FP_WORKBITS);                      \
263               if (X##_e >= _FP_EXPMAX_##fs)                             \
264                 {                                                       \
265                   /* overflow */                                        \
266                   switch (FP_ROUNDMODE)                                 \
267                     {                                                   \
268                     case FP_RND_NEAREST:                                \
269                       X##_c = FP_CLS_INF;                               \
270                       break;                                            \
271                     case FP_RND_PINF:                                   \
272                       if (!X##_s)                                       \
273                         X##_c = FP_CLS_INF;                             \
274                       break;                                            \
275                     case FP_RND_MINF:                                   \
276                       if (X##_s)                                        \
277                         X##_c = FP_CLS_INF;                             \
278                       break;                                            \
279                     }                                                   \
280                   if (X##_c == FP_CLS_INF)                              \
281                     {                                                   \
282                       /* Overflow to infinity */                        \
283                       X##_e = _FP_EXPMAX_##fs;                          \
284                       _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc);         \
285                     }                                                   \
286                   else                                                  \
287                     {                                                   \
288                       /* Overflow to maximum normal */                  \
289                       X##_e = _FP_EXPMAX_##fs - 1;                      \
290                       _FP_FRAC_SET_##wc (X, _FP_MAXFRAC_##wc);          \
291                     }                                                   \
292                   FP_SET_EXCEPTION (FP_EX_OVERFLOW);                    \
293                   FP_SET_EXCEPTION (FP_EX_INEXACT);                     \
294                 }                                                       \
295             }                                                           \
296           else                                                          \
297             {                                                           \
298               /* we've got a denormalized number */                     \
299               int _FP_PACK_CANONICAL_is_tiny = 1;                       \
300               if (_FP_TININESS_AFTER_ROUNDING && X##_e == 0)            \
301                 {                                                       \
302                   FP_DECL_##fs (_FP_PACK_CANONICAL_T);                  \
303                   _FP_FRAC_COPY_##wc (_FP_PACK_CANONICAL_T, X);         \
304                   _FP_PACK_CANONICAL_T##_s = X##_s;                     \
305                   _FP_PACK_CANONICAL_T##_e = X##_e;                     \
306                   _FP_ROUND (wc, _FP_PACK_CANONICAL_T);                 \
307                   if (_FP_FRAC_OVERP_##wc (fs, _FP_PACK_CANONICAL_T))   \
308                     _FP_PACK_CANONICAL_is_tiny = 0;                     \
309                 }                                                       \
310               X##_e = -X##_e + 1;                                       \
311               if (X##_e <= _FP_WFRACBITS_##fs)                          \
312                 {                                                       \
313                   _FP_FRAC_SRS_##wc (X, X##_e, _FP_WFRACBITS_##fs);     \
314                   _FP_ROUND (wc, X);                                    \
315                   if (_FP_FRAC_HIGH_##fs (X)                            \
316                       & (_FP_OVERFLOW_##fs >> 1))                       \
317                     {                                                   \
318                       X##_e = 1;                                        \
319                       _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc);         \
320                       FP_SET_EXCEPTION (FP_EX_INEXACT);                 \
321                     }                                                   \
322                   else                                                  \
323                     {                                                   \
324                       X##_e = 0;                                        \
325                       _FP_FRAC_SRL_##wc (X, _FP_WORKBITS);              \
326                     }                                                   \
327                   if (_FP_PACK_CANONICAL_is_tiny                        \
328                       && ((FP_CUR_EXCEPTIONS & FP_EX_INEXACT)           \
329                           || (FP_TRAPPING_EXCEPTIONS                    \
330                               & FP_EX_UNDERFLOW)))                      \
331                     FP_SET_EXCEPTION (FP_EX_UNDERFLOW);                 \
332                 }                                                       \
333               else                                                      \
334                 {                                                       \
335                   /* underflow to zero */                               \
336                   X##_e = 0;                                            \
337                   if (!_FP_FRAC_ZEROP_##wc (X))                         \
338                     {                                                   \
339                       _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc);          \
340                       _FP_ROUND (wc, X);                                \
341                       _FP_FRAC_LOW_##wc (X) >>= (_FP_WORKBITS);         \
342                     }                                                   \
343                   FP_SET_EXCEPTION (FP_EX_UNDERFLOW);                   \
344                 }                                                       \
345             }                                                           \
346           break;                                                        \
347                                                                         \
348         case FP_CLS_ZERO:                                               \
349           X##_e = 0;                                                    \
350           _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc);                     \
351           break;                                                        \
352                                                                         \
353         case FP_CLS_INF:                                                \
354           X##_e = _FP_EXPMAX_##fs;                                      \
355           _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc);                     \
356           break;                                                        \
357                                                                         \
358         case FP_CLS_NAN:                                                \
359           X##_e = _FP_EXPMAX_##fs;                                      \
360           if (!_FP_KEEPNANFRACP)                                        \
361             {                                                           \
362               _FP_FRAC_SET_##wc (X, _FP_NANFRAC_##fs);                  \
363               X##_s = _FP_NANSIGN_##fs;                                 \
364             }                                                           \
365           else                                                          \
366             _FP_SETQNAN (fs, wc, X);                                    \
367           break;                                                        \
368         }                                                               \
369     }                                                                   \
370   while (0)
371
372 /* This one accepts raw argument and not cooked,  returns
373  * 1 if X is a signaling NaN.
374  */
375 #define _FP_ISSIGNAN(fs, wc, X)                 \
376   ({                                            \
377     int _FP_ISSIGNAN_ret = 0;                   \
378     if (X##_e == _FP_EXPMAX_##fs)               \
379       {                                         \
380         if (!_FP_FRAC_ZEROP_##wc (X)            \
381             && _FP_FRAC_SNANP (fs, X))          \
382           _FP_ISSIGNAN_ret = 1;                 \
383       }                                         \
384     _FP_ISSIGNAN_ret;                           \
385   })
386
387
388
389
390
391 /* Addition on semi-raw values.  */
392 #define _FP_ADD_INTERNAL(fs, wc, R, X, Y, OP)                           \
393   do                                                                    \
394     {                                                                   \
395       if (X##_s == Y##_s)                                               \
396         {                                                               \
397           /* Addition.  */                                              \
398           R##_s = X##_s;                                                \
399           int _FP_ADD_INTERNAL_ediff = X##_e - Y##_e;                   \
400           if (_FP_ADD_INTERNAL_ediff > 0)                               \
401             {                                                           \
402               R##_e = X##_e;                                            \
403               if (Y##_e == 0)                                           \
404                 {                                                       \
405                   /* Y is zero or denormalized.  */                     \
406                   if (_FP_FRAC_ZEROP_##wc (Y))                          \
407                     {                                                   \
408                       _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X);             \
409                       _FP_FRAC_COPY_##wc (R, X);                        \
410                       goto add_done;                                    \
411                     }                                                   \
412                   else                                                  \
413                     {                                                   \
414                       FP_SET_EXCEPTION (FP_EX_DENORM);                  \
415                       _FP_ADD_INTERNAL_ediff--;                         \
416                       if (_FP_ADD_INTERNAL_ediff == 0)                  \
417                         {                                               \
418                           _FP_FRAC_ADD_##wc (R, X, Y);                  \
419                           goto add3;                                    \
420                         }                                               \
421                       if (X##_e == _FP_EXPMAX_##fs)                     \
422                         {                                               \
423                           _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X);         \
424                           _FP_FRAC_COPY_##wc (R, X);                    \
425                           goto add_done;                                \
426                         }                                               \
427                       goto add1;                                        \
428                     }                                                   \
429                 }                                                       \
430               else if (X##_e == _FP_EXPMAX_##fs)                        \
431                 {                                                       \
432                   /* X is NaN or Inf, Y is normal.  */                  \
433                   _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X);                 \
434                   _FP_FRAC_COPY_##wc (R, X);                            \
435                   goto add_done;                                        \
436                 }                                                       \
437                                                                         \
438               /* Insert implicit MSB of Y.  */                          \
439               _FP_FRAC_HIGH_##fs (Y) |= _FP_IMPLBIT_SH_##fs;            \
440                                                                         \
441             add1:                                                       \
442               /* Shift the mantissa of Y to the right                   \
443                  _FP_ADD_INTERNAL_EDIFF steps; remember to account      \
444                  later for the implicit MSB of X.  */                   \
445               if (_FP_ADD_INTERNAL_ediff <= _FP_WFRACBITS_##fs)         \
446                 _FP_FRAC_SRS_##wc (Y, _FP_ADD_INTERNAL_ediff,           \
447                                    _FP_WFRACBITS_##fs);                 \
448               else if (!_FP_FRAC_ZEROP_##wc (Y))                        \
449                 _FP_FRAC_SET_##wc (Y, _FP_MINFRAC_##wc);                \
450               _FP_FRAC_ADD_##wc (R, X, Y);                              \
451             }                                                           \
452           else if (_FP_ADD_INTERNAL_ediff < 0)                          \
453             {                                                           \
454               _FP_ADD_INTERNAL_ediff = -_FP_ADD_INTERNAL_ediff;         \
455               R##_e = Y##_e;                                            \
456               if (X##_e == 0)                                           \
457                 {                                                       \
458                   /* X is zero or denormalized.  */                     \
459                   if (_FP_FRAC_ZEROP_##wc (X))                          \
460                     {                                                   \
461                       _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y);             \
462                       _FP_FRAC_COPY_##wc (R, Y);                        \
463                       goto add_done;                                    \
464                     }                                                   \
465                   else                                                  \
466                     {                                                   \
467                       FP_SET_EXCEPTION (FP_EX_DENORM);                  \
468                       _FP_ADD_INTERNAL_ediff--;                         \
469                       if (_FP_ADD_INTERNAL_ediff == 0)                  \
470                         {                                               \
471                           _FP_FRAC_ADD_##wc (R, Y, X);                  \
472                           goto add3;                                    \
473                         }                                               \
474                       if (Y##_e == _FP_EXPMAX_##fs)                     \
475                         {                                               \
476                           _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y);         \
477                           _FP_FRAC_COPY_##wc (R, Y);                    \
478                           goto add_done;                                \
479                         }                                               \
480                       goto add2;                                        \
481                     }                                                   \
482                 }                                                       \
483               else if (Y##_e == _FP_EXPMAX_##fs)                        \
484                 {                                                       \
485                   /* Y is NaN or Inf, X is normal.  */                  \
486                   _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y);                 \
487                   _FP_FRAC_COPY_##wc (R, Y);                            \
488                   goto add_done;                                        \
489                 }                                                       \
490                                                                         \
491               /* Insert implicit MSB of X.  */                          \
492               _FP_FRAC_HIGH_##fs (X) |= _FP_IMPLBIT_SH_##fs;            \
493                                                                         \
494             add2:                                                       \
495               /* Shift the mantissa of X to the right                   \
496                  _FP_ADD_INTERNAL_EDIFF steps; remember to account      \
497                  later for the implicit MSB of Y.  */                   \
498               if (_FP_ADD_INTERNAL_ediff <= _FP_WFRACBITS_##fs)         \
499                 _FP_FRAC_SRS_##wc (X, _FP_ADD_INTERNAL_ediff,           \
500                                    _FP_WFRACBITS_##fs);                 \
501               else if (!_FP_FRAC_ZEROP_##wc (X))                        \
502                 _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc);                \
503               _FP_FRAC_ADD_##wc (R, Y, X);                              \
504             }                                                           \
505           else                                                          \
506             {                                                           \
507               /* _FP_ADD_INTERNAL_ediff == 0.  */                       \
508               if (!_FP_EXP_NORMAL (fs, wc, X))                          \
509                 {                                                       \
510                   if (X##_e == 0)                                       \
511                     {                                                   \
512                       /* X and Y are zero or denormalized.  */          \
513                       R##_e = 0;                                        \
514                       if (_FP_FRAC_ZEROP_##wc (X))                      \
515                         {                                               \
516                           if (!_FP_FRAC_ZEROP_##wc (Y))                 \
517                             FP_SET_EXCEPTION (FP_EX_DENORM);            \
518                           _FP_FRAC_COPY_##wc (R, Y);                    \
519                           goto add_done;                                \
520                         }                                               \
521                       else if (_FP_FRAC_ZEROP_##wc (Y))                 \
522                         {                                               \
523                           FP_SET_EXCEPTION (FP_EX_DENORM);              \
524                           _FP_FRAC_COPY_##wc (R, X);                    \
525                           goto add_done;                                \
526                         }                                               \
527                       else                                              \
528                         {                                               \
529                           FP_SET_EXCEPTION (FP_EX_DENORM);              \
530                           _FP_FRAC_ADD_##wc (R, X, Y);                  \
531                           if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
532                             {                                           \
533                               /* Normalized result.  */                 \
534                               _FP_FRAC_HIGH_##fs (R)                    \
535                                 &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs;   \
536                               R##_e = 1;                                \
537                             }                                           \
538                           goto add_done;                                \
539                         }                                               \
540                     }                                                   \
541                   else                                                  \
542                     {                                                   \
543                       /* X and Y are NaN or Inf.  */                    \
544                       _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X);             \
545                       _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y);             \
546                       R##_e = _FP_EXPMAX_##fs;                          \
547                       if (_FP_FRAC_ZEROP_##wc (X))                      \
548                         _FP_FRAC_COPY_##wc (R, Y);                      \
549                       else if (_FP_FRAC_ZEROP_##wc (Y))                 \
550                         _FP_FRAC_COPY_##wc (R, X);                      \
551                       else                                              \
552                         _FP_CHOOSENAN_SEMIRAW (fs, wc, R, X, Y, OP);    \
553                       goto add_done;                                    \
554                     }                                                   \
555                 }                                                       \
556               /* The exponents of X and Y, both normal, are equal.  The \
557                  implicit MSBs will always add to increase the          \
558                  exponent.  */                                          \
559               _FP_FRAC_ADD_##wc (R, X, Y);                              \
560               R##_e = X##_e + 1;                                        \
561               _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs);             \
562               if (R##_e == _FP_EXPMAX_##fs)                             \
563                 /* Overflow to infinity (depending on rounding mode).  */ \
564                 _FP_OVERFLOW_SEMIRAW (fs, wc, R);                       \
565               goto add_done;                                            \
566             }                                                           \
567         add3:                                                           \
568           if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs)             \
569             {                                                           \
570               /* Overflow.  */                                          \
571               _FP_FRAC_HIGH_##fs (R) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
572               R##_e++;                                                  \
573               _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs);             \
574               if (R##_e == _FP_EXPMAX_##fs)                             \
575                 /* Overflow to infinity (depending on rounding mode).  */ \
576                 _FP_OVERFLOW_SEMIRAW (fs, wc, R);                       \
577             }                                                           \
578         add_done: ;                                                     \
579         }                                                               \
580       else                                                              \
581         {                                                               \
582           /* Subtraction.  */                                           \
583           int _FP_ADD_INTERNAL_ediff = X##_e - Y##_e;                   \
584           if (_FP_ADD_INTERNAL_ediff > 0)                               \
585             {                                                           \
586               R##_e = X##_e;                                            \
587               R##_s = X##_s;                                            \
588               if (Y##_e == 0)                                           \
589                 {                                                       \
590                   /* Y is zero or denormalized.  */                     \
591                   if (_FP_FRAC_ZEROP_##wc (Y))                          \
592                     {                                                   \
593                       _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X);             \
594                       _FP_FRAC_COPY_##wc (R, X);                        \
595                       goto sub_done;                                    \
596                     }                                                   \
597                   else                                                  \
598                     {                                                   \
599                       FP_SET_EXCEPTION (FP_EX_DENORM);                  \
600                       _FP_ADD_INTERNAL_ediff--;                         \
601                       if (_FP_ADD_INTERNAL_ediff == 0)                  \
602                         {                                               \
603                           _FP_FRAC_SUB_##wc (R, X, Y);                  \
604                           goto sub3;                                    \
605                         }                                               \
606                       if (X##_e == _FP_EXPMAX_##fs)                     \
607                         {                                               \
608                           _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X);         \
609                           _FP_FRAC_COPY_##wc (R, X);                    \
610                           goto sub_done;                                \
611                         }                                               \
612                       goto sub1;                                        \
613                     }                                                   \
614                 }                                                       \
615               else if (X##_e == _FP_EXPMAX_##fs)                        \
616                 {                                                       \
617                   /* X is NaN or Inf, Y is normal.  */                  \
618                   _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X);                 \
619                   _FP_FRAC_COPY_##wc (R, X);                            \
620                   goto sub_done;                                        \
621                 }                                                       \
622                                                                         \
623               /* Insert implicit MSB of Y.  */                          \
624               _FP_FRAC_HIGH_##fs (Y) |= _FP_IMPLBIT_SH_##fs;            \
625                                                                         \
626             sub1:                                                       \
627               /* Shift the mantissa of Y to the right                   \
628                  _FP_ADD_INTERNAL_EDIFF steps; remember to account      \
629                  later for the implicit MSB of X.  */                   \
630               if (_FP_ADD_INTERNAL_ediff <= _FP_WFRACBITS_##fs)         \
631                 _FP_FRAC_SRS_##wc (Y, _FP_ADD_INTERNAL_ediff,           \
632                                    _FP_WFRACBITS_##fs);                 \
633               else if (!_FP_FRAC_ZEROP_##wc (Y))                        \
634                 _FP_FRAC_SET_##wc (Y, _FP_MINFRAC_##wc);                \
635               _FP_FRAC_SUB_##wc (R, X, Y);                              \
636             }                                                           \
637           else if (_FP_ADD_INTERNAL_ediff < 0)                          \
638             {                                                           \
639               _FP_ADD_INTERNAL_ediff = -_FP_ADD_INTERNAL_ediff;         \
640               R##_e = Y##_e;                                            \
641               R##_s = Y##_s;                                            \
642               if (X##_e == 0)                                           \
643                 {                                                       \
644                   /* X is zero or denormalized.  */                     \
645                   if (_FP_FRAC_ZEROP_##wc (X))                          \
646                     {                                                   \
647                       _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y);             \
648                       _FP_FRAC_COPY_##wc (R, Y);                        \
649                       goto sub_done;                                    \
650                     }                                                   \
651                   else                                                  \
652                     {                                                   \
653                       FP_SET_EXCEPTION (FP_EX_DENORM);                  \
654                       _FP_ADD_INTERNAL_ediff--;                         \
655                       if (_FP_ADD_INTERNAL_ediff == 0)                  \
656                         {                                               \
657                           _FP_FRAC_SUB_##wc (R, Y, X);                  \
658                           goto sub3;                                    \
659                         }                                               \
660                       if (Y##_e == _FP_EXPMAX_##fs)                     \
661                         {                                               \
662                           _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y);         \
663                           _FP_FRAC_COPY_##wc (R, Y);                    \
664                           goto sub_done;                                \
665                         }                                               \
666                       goto sub2;                                        \
667                     }                                                   \
668                 }                                                       \
669               else if (Y##_e == _FP_EXPMAX_##fs)                        \
670                 {                                                       \
671                   /* Y is NaN or Inf, X is normal.  */                  \
672                   _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y);                 \
673                   _FP_FRAC_COPY_##wc (R, Y);                            \
674                   goto sub_done;                                        \
675                 }                                                       \
676                                                                         \
677               /* Insert implicit MSB of X.  */                          \
678               _FP_FRAC_HIGH_##fs (X) |= _FP_IMPLBIT_SH_##fs;            \
679                                                                         \
680             sub2:                                                       \
681               /* Shift the mantissa of X to the right                   \
682                  _FP_ADD_INTERNAL_EDIFF steps; remember to account      \
683                  later for the implicit MSB of Y.  */                   \
684               if (_FP_ADD_INTERNAL_ediff <= _FP_WFRACBITS_##fs)         \
685                 _FP_FRAC_SRS_##wc (X, _FP_ADD_INTERNAL_ediff,           \
686                                    _FP_WFRACBITS_##fs);                 \
687               else if (!_FP_FRAC_ZEROP_##wc (X))                        \
688                 _FP_FRAC_SET_##wc (X, _FP_MINFRAC_##wc);                \
689               _FP_FRAC_SUB_##wc (R, Y, X);                              \
690             }                                                           \
691           else                                                          \
692             {                                                           \
693               /* ediff == 0.  */                                        \
694               if (!_FP_EXP_NORMAL (fs, wc, X))                          \
695                 {                                                       \
696                   if (X##_e == 0)                                       \
697                     {                                                   \
698                       /* X and Y are zero or denormalized.  */          \
699                       R##_e = 0;                                        \
700                       if (_FP_FRAC_ZEROP_##wc (X))                      \
701                         {                                               \
702                           _FP_FRAC_COPY_##wc (R, Y);                    \
703                           if (_FP_FRAC_ZEROP_##wc (Y))                  \
704                             R##_s = (FP_ROUNDMODE == FP_RND_MINF);      \
705                           else                                          \
706                             {                                           \
707                               FP_SET_EXCEPTION (FP_EX_DENORM);          \
708                               R##_s = Y##_s;                            \
709                             }                                           \
710                           goto sub_done;                                \
711                         }                                               \
712                       else if (_FP_FRAC_ZEROP_##wc (Y))                 \
713                         {                                               \
714                           FP_SET_EXCEPTION (FP_EX_DENORM);              \
715                           _FP_FRAC_COPY_##wc (R, X);                    \
716                           R##_s = X##_s;                                \
717                           goto sub_done;                                \
718                         }                                               \
719                       else                                              \
720                         {                                               \
721                           FP_SET_EXCEPTION (FP_EX_DENORM);              \
722                           _FP_FRAC_SUB_##wc (R, X, Y);                  \
723                           R##_s = X##_s;                                \
724                           if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs) \
725                             {                                           \
726                               /* |X| < |Y|, negate result.  */          \
727                               _FP_FRAC_SUB_##wc (R, Y, X);              \
728                               R##_s = Y##_s;                            \
729                             }                                           \
730                           else if (_FP_FRAC_ZEROP_##wc (R))             \
731                             R##_s = (FP_ROUNDMODE == FP_RND_MINF);      \
732                           goto sub_done;                                \
733                         }                                               \
734                     }                                                   \
735                   else                                                  \
736                     {                                                   \
737                       /* X and Y are NaN or Inf, of opposite signs.  */ \
738                       _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, X);             \
739                       _FP_CHECK_SIGNAN_SEMIRAW (fs, wc, Y);             \
740                       R##_e = _FP_EXPMAX_##fs;                          \
741                       if (_FP_FRAC_ZEROP_##wc (X))                      \
742                         {                                               \
743                           if (_FP_FRAC_ZEROP_##wc (Y))                  \
744                             {                                           \
745                               /* Inf - Inf.  */                         \
746                               R##_s = _FP_NANSIGN_##fs;                 \
747                               _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs);  \
748                               _FP_FRAC_SLL_##wc (R, _FP_WORKBITS);      \
749                               FP_SET_EXCEPTION (FP_EX_INVALID);         \
750                             }                                           \
751                           else                                          \
752                             {                                           \
753                               /* Inf - NaN.  */                         \
754                               R##_s = Y##_s;                            \
755                               _FP_FRAC_COPY_##wc (R, Y);                \
756                             }                                           \
757                         }                                               \
758                       else                                              \
759                         {                                               \
760                           if (_FP_FRAC_ZEROP_##wc (Y))                  \
761                             {                                           \
762                               /* NaN - Inf.  */                         \
763                               R##_s = X##_s;                            \
764                               _FP_FRAC_COPY_##wc (R, X);                \
765                             }                                           \
766                           else                                          \
767                             {                                           \
768                               /* NaN - NaN.  */                         \
769                               _FP_CHOOSENAN_SEMIRAW (fs, wc, R, X, Y, OP); \
770                             }                                           \
771                         }                                               \
772                       goto sub_done;                                    \
773                     }                                                   \
774                 }                                                       \
775               /* The exponents of X and Y, both normal, are equal.  The \
776                  implicit MSBs cancel.  */                              \
777               R##_e = X##_e;                                            \
778               _FP_FRAC_SUB_##wc (R, X, Y);                              \
779               R##_s = X##_s;                                            \
780               if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs)         \
781                 {                                                       \
782                   /* |X| < |Y|, negate result.  */                      \
783                   _FP_FRAC_SUB_##wc (R, Y, X);                          \
784                   R##_s = Y##_s;                                        \
785                 }                                                       \
786               else if (_FP_FRAC_ZEROP_##wc (R))                         \
787                 {                                                       \
788                   R##_e = 0;                                            \
789                   R##_s = (FP_ROUNDMODE == FP_RND_MINF);                \
790                   goto sub_done;                                        \
791                 }                                                       \
792               goto norm;                                                \
793             }                                                           \
794         sub3:                                                           \
795           if (_FP_FRAC_HIGH_##fs (R) & _FP_IMPLBIT_SH_##fs)             \
796             {                                                           \
797               int _FP_ADD_INTERNAL_diff;                                \
798               /* Carry into most significant bit of larger one of X and Y, \
799                  canceling it; renormalize.  */                         \
800               _FP_FRAC_HIGH_##fs (R) &= _FP_IMPLBIT_SH_##fs - 1;        \
801             norm:                                                       \
802               _FP_FRAC_CLZ_##wc (_FP_ADD_INTERNAL_diff, R);             \
803               _FP_ADD_INTERNAL_diff -= _FP_WFRACXBITS_##fs;             \
804               _FP_FRAC_SLL_##wc (R, _FP_ADD_INTERNAL_diff);             \
805               if (R##_e <= _FP_ADD_INTERNAL_diff)                       \
806                 {                                                       \
807                   /* R is denormalized.  */                             \
808                   _FP_ADD_INTERNAL_diff                                 \
809                     = _FP_ADD_INTERNAL_diff - R##_e + 1;                \
810                   _FP_FRAC_SRS_##wc (R, _FP_ADD_INTERNAL_diff,          \
811                                      _FP_WFRACBITS_##fs);               \
812                   R##_e = 0;                                            \
813                 }                                                       \
814               else                                                      \
815                 {                                                       \
816                   R##_e -= _FP_ADD_INTERNAL_diff;                       \
817                   _FP_FRAC_HIGH_##fs (R) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
818                 }                                                       \
819             }                                                           \
820         sub_done: ;                                                     \
821         }                                                               \
822     }                                                                   \
823   while (0)
824
825 #define _FP_ADD(fs, wc, R, X, Y) _FP_ADD_INTERNAL (fs, wc, R, X, Y, '+')
826 #define _FP_SUB(fs, wc, R, X, Y)                                        \
827   do                                                                    \
828     {                                                                   \
829       if (!(Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y)))      \
830         Y##_s ^= 1;                                                     \
831       _FP_ADD_INTERNAL (fs, wc, R, X, Y, '-');                          \
832     }                                                                   \
833   while (0)
834
835
836 /*
837  * Main negation routine.  The input value is raw.
838  */
839
840 #define _FP_NEG(fs, wc, R, X)                   \
841   do                                            \
842     {                                           \
843       _FP_FRAC_COPY_##wc (R, X);                \
844       R##_e = X##_e;                            \
845       R##_s = 1 ^ X##_s;                        \
846     }                                           \
847   while (0)
848
849
850 /*
851  * Main multiplication routine.  The input values should be cooked.
852  */
853
854 #define _FP_MUL(fs, wc, R, X, Y)                                \
855   do                                                            \
856     {                                                           \
857       R##_s = X##_s ^ Y##_s;                                    \
858       R##_e = X##_e + Y##_e + 1;                                \
859       switch (_FP_CLS_COMBINE (X##_c, Y##_c))                   \
860         {                                                       \
861         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL):    \
862           R##_c = FP_CLS_NORMAL;                                \
863                                                                 \
864           _FP_MUL_MEAT_##fs (R, X, Y);                          \
865                                                                 \
866           if (_FP_FRAC_OVERP_##wc (fs, R))                      \
867             _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs);       \
868           else                                                  \
869             R##_e--;                                            \
870           break;                                                \
871                                                                 \
872         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN):          \
873           _FP_CHOOSENAN (fs, wc, R, X, Y, '*');                 \
874           break;                                                \
875                                                                 \
876         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL):       \
877         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF):          \
878         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO):         \
879           R##_s = X##_s;                                        \
880                                                                 \
881         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF):          \
882         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL):       \
883         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL):      \
884         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO):        \
885           _FP_FRAC_COPY_##wc (R, X);                            \
886           R##_c = X##_c;                                        \
887           break;                                                \
888                                                                 \
889         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN):       \
890         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN):          \
891         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN):         \
892           R##_s = Y##_s;                                        \
893                                                                 \
894         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF):       \
895         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO):      \
896           _FP_FRAC_COPY_##wc (R, Y);                            \
897           R##_c = Y##_c;                                        \
898           break;                                                \
899                                                                 \
900         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO):         \
901         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF):         \
902           R##_s = _FP_NANSIGN_##fs;                             \
903           R##_c = FP_CLS_NAN;                                   \
904           _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs);              \
905           FP_SET_EXCEPTION (FP_EX_INVALID);                     \
906           break;                                                \
907                                                                 \
908         default:                                                \
909           abort ();                                             \
910         }                                                       \
911     }                                                           \
912   while (0)
913
914
915 /* Fused multiply-add.  The input values should be cooked.  */
916
917 #define _FP_FMA(fs, wc, dwc, R, X, Y, Z)                                \
918   do                                                                    \
919     {                                                                   \
920       FP_DECL_##fs (_FP_FMA_T);                                         \
921       _FP_FMA_T##_s = X##_s ^ Y##_s;                                    \
922       _FP_FMA_T##_e = X##_e + Y##_e + 1;                                \
923       switch (_FP_CLS_COMBINE (X##_c, Y##_c))                           \
924         {                                                               \
925         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL):            \
926           switch (Z##_c)                                                \
927             {                                                           \
928             case FP_CLS_INF:                                            \
929             case FP_CLS_NAN:                                            \
930               R##_s = Z##_s;                                            \
931               _FP_FRAC_COPY_##wc (R, Z);                                \
932               R##_c = Z##_c;                                            \
933               break;                                                    \
934                                                                         \
935             case FP_CLS_ZERO:                                           \
936               R##_c = FP_CLS_NORMAL;                                    \
937               R##_s = _FP_FMA_T##_s;                                    \
938               R##_e = _FP_FMA_T##_e;                                    \
939                                                                         \
940               _FP_MUL_MEAT_##fs (R, X, Y);                              \
941                                                                         \
942               if (_FP_FRAC_OVERP_##wc (fs, R))                          \
943                 _FP_FRAC_SRS_##wc (R, 1, _FP_WFRACBITS_##fs);           \
944               else                                                      \
945                 R##_e--;                                                \
946               break;                                                    \
947                                                                         \
948             case FP_CLS_NORMAL:;                                        \
949               _FP_FRAC_DECL_##dwc (_FP_FMA_TD);                         \
950               _FP_FRAC_DECL_##dwc (_FP_FMA_ZD);                         \
951               _FP_FRAC_DECL_##dwc (_FP_FMA_RD);                         \
952               _FP_MUL_MEAT_DW_##fs (_FP_FMA_TD, X, Y);                  \
953               R##_e = _FP_FMA_T##_e;                                    \
954               int _FP_FMA_tsh                                           \
955                 = _FP_FRAC_HIGHBIT_DW_##dwc (fs, _FP_FMA_TD) == 0;      \
956               _FP_FMA_T##_e -= _FP_FMA_tsh;                             \
957               int _FP_FMA_ediff = _FP_FMA_T##_e - Z##_e;                \
958               if (_FP_FMA_ediff >= 0)                                   \
959                 {                                                       \
960                   int _FP_FMA_shift                                     \
961                     = _FP_WFRACBITS_##fs - _FP_FMA_tsh - _FP_FMA_ediff; \
962                   if (_FP_FMA_shift <= -_FP_WFRACBITS_##fs)             \
963                     _FP_FRAC_SET_##dwc (_FP_FMA_ZD, _FP_MINFRAC_##dwc); \
964                   else                                                  \
965                     {                                                   \
966                       _FP_FRAC_COPY_##dwc##_##wc (_FP_FMA_ZD, Z);       \
967                       if (_FP_FMA_shift < 0)                            \
968                         _FP_FRAC_SRS_##dwc (_FP_FMA_ZD, -_FP_FMA_shift, \
969                                             _FP_WFRACBITS_DW_##fs);     \
970                       else if (_FP_FMA_shift > 0)                       \
971                         _FP_FRAC_SLL_##dwc (_FP_FMA_ZD, _FP_FMA_shift); \
972                     }                                                   \
973                   R##_s = _FP_FMA_T##_s;                                \
974                   if (_FP_FMA_T##_s == Z##_s)                           \
975                     _FP_FRAC_ADD_##dwc (_FP_FMA_RD, _FP_FMA_TD,         \
976                                         _FP_FMA_ZD);                    \
977                   else                                                  \
978                     {                                                   \
979                       _FP_FRAC_SUB_##dwc (_FP_FMA_RD, _FP_FMA_TD,       \
980                                           _FP_FMA_ZD);                  \
981                       if (_FP_FRAC_NEGP_##dwc (_FP_FMA_RD))             \
982                         {                                               \
983                           R##_s = Z##_s;                                \
984                           _FP_FRAC_SUB_##dwc (_FP_FMA_RD, _FP_FMA_ZD,   \
985                                               _FP_FMA_TD);              \
986                         }                                               \
987                     }                                                   \
988                 }                                                       \
989               else                                                      \
990                 {                                                       \
991                   R##_e = Z##_e;                                        \
992                   R##_s = Z##_s;                                        \
993                   _FP_FRAC_COPY_##dwc##_##wc (_FP_FMA_ZD, Z);           \
994                   _FP_FRAC_SLL_##dwc (_FP_FMA_ZD, _FP_WFRACBITS_##fs);  \
995                   int _FP_FMA_shift = -_FP_FMA_ediff - _FP_FMA_tsh;     \
996                   if (_FP_FMA_shift >= _FP_WFRACBITS_DW_##fs)           \
997                     _FP_FRAC_SET_##dwc (_FP_FMA_TD, _FP_MINFRAC_##dwc); \
998                   else if (_FP_FMA_shift > 0)                           \
999                     _FP_FRAC_SRS_##dwc (_FP_FMA_TD, _FP_FMA_shift,      \
1000                                         _FP_WFRACBITS_DW_##fs);         \
1001                   if (Z##_s == _FP_FMA_T##_s)                           \
1002                     _FP_FRAC_ADD_##dwc (_FP_FMA_RD, _FP_FMA_ZD,         \
1003                                         _FP_FMA_TD);                    \
1004                   else                                                  \
1005                     _FP_FRAC_SUB_##dwc (_FP_FMA_RD, _FP_FMA_ZD,         \
1006                                         _FP_FMA_TD);                    \
1007                 }                                                       \
1008               if (_FP_FRAC_ZEROP_##dwc (_FP_FMA_RD))                    \
1009                 {                                                       \
1010                   if (_FP_FMA_T##_s == Z##_s)                           \
1011                     R##_s = Z##_s;                                      \
1012                   else                                                  \
1013                     R##_s = (FP_ROUNDMODE == FP_RND_MINF);              \
1014                   _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc);             \
1015                   R##_c = FP_CLS_ZERO;                                  \
1016                 }                                                       \
1017               else                                                      \
1018                 {                                                       \
1019                   int _FP_FMA_rlz;                                      \
1020                   _FP_FRAC_CLZ_##dwc (_FP_FMA_rlz, _FP_FMA_RD);         \
1021                   _FP_FMA_rlz -= _FP_WFRACXBITS_DW_##fs;                \
1022                   R##_e -= _FP_FMA_rlz;                                 \
1023                   int _FP_FMA_shift = _FP_WFRACBITS_##fs - _FP_FMA_rlz; \
1024                   if (_FP_FMA_shift > 0)                                \
1025                     _FP_FRAC_SRS_##dwc (_FP_FMA_RD, _FP_FMA_shift,      \
1026                                         _FP_WFRACBITS_DW_##fs);         \
1027                   else if (_FP_FMA_shift < 0)                           \
1028                     _FP_FRAC_SLL_##dwc (_FP_FMA_RD, -_FP_FMA_shift);    \
1029                   _FP_FRAC_COPY_##wc##_##dwc (R, _FP_FMA_RD);           \
1030                   R##_c = FP_CLS_NORMAL;                                \
1031                 }                                                       \
1032               break;                                                    \
1033             }                                                           \
1034           goto done_fma;                                                \
1035                                                                         \
1036         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN):                  \
1037           _FP_CHOOSENAN (fs, wc, _FP_FMA_T, X, Y, '*');                 \
1038           break;                                                        \
1039                                                                         \
1040         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL):               \
1041         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF):                  \
1042         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO):                 \
1043           _FP_FMA_T##_s = X##_s;                                        \
1044                                                                         \
1045         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF):                  \
1046         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL):               \
1047         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL):              \
1048         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO):                \
1049           _FP_FRAC_COPY_##wc (_FP_FMA_T, X);                            \
1050           _FP_FMA_T##_c = X##_c;                                        \
1051           break;                                                        \
1052                                                                         \
1053         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN):               \
1054         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN):                  \
1055         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN):                 \
1056           _FP_FMA_T##_s = Y##_s;                                        \
1057                                                                         \
1058         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF):               \
1059         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO):              \
1060           _FP_FRAC_COPY_##wc (_FP_FMA_T, Y);                            \
1061           _FP_FMA_T##_c = Y##_c;                                        \
1062           break;                                                        \
1063                                                                         \
1064         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO):                 \
1065         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF):                 \
1066           _FP_FMA_T##_s = _FP_NANSIGN_##fs;                             \
1067           _FP_FMA_T##_c = FP_CLS_NAN;                                   \
1068           _FP_FRAC_SET_##wc (_FP_FMA_T, _FP_NANFRAC_##fs);              \
1069           FP_SET_EXCEPTION (FP_EX_INVALID);                             \
1070           break;                                                        \
1071                                                                         \
1072         default:                                                        \
1073           abort ();                                                     \
1074         }                                                               \
1075                                                                         \
1076       /* T = X * Y is zero, infinity or NaN.  */                        \
1077       switch (_FP_CLS_COMBINE (_FP_FMA_T##_c, Z##_c))                   \
1078         {                                                               \
1079         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN):                  \
1080           _FP_CHOOSENAN (fs, wc, R, _FP_FMA_T, Z, '+');                 \
1081           break;                                                        \
1082                                                                         \
1083         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL):               \
1084         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF):                  \
1085         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO):                 \
1086         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL):               \
1087         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO):                 \
1088           R##_s = _FP_FMA_T##_s;                                        \
1089           _FP_FRAC_COPY_##wc (R, _FP_FMA_T);                            \
1090           R##_c = _FP_FMA_T##_c;                                        \
1091           break;                                                        \
1092                                                                         \
1093         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN):                  \
1094         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN):                 \
1095         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL):              \
1096         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF):                 \
1097           R##_s = Z##_s;                                                \
1098           _FP_FRAC_COPY_##wc (R, Z);                                    \
1099           R##_c = Z##_c;                                                \
1100           break;                                                        \
1101                                                                         \
1102         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF):                  \
1103           if (_FP_FMA_T##_s == Z##_s)                                   \
1104             {                                                           \
1105               R##_s = Z##_s;                                            \
1106               _FP_FRAC_COPY_##wc (R, Z);                                \
1107               R##_c = Z##_c;                                            \
1108             }                                                           \
1109           else                                                          \
1110             {                                                           \
1111               R##_s = _FP_NANSIGN_##fs;                                 \
1112               R##_c = FP_CLS_NAN;                                       \
1113               _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs);                  \
1114               FP_SET_EXCEPTION (FP_EX_INVALID);                         \
1115             }                                                           \
1116           break;                                                        \
1117                                                                         \
1118         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO):                \
1119           if (_FP_FMA_T##_s == Z##_s)                                   \
1120             R##_s = Z##_s;                                              \
1121           else                                                          \
1122             R##_s = (FP_ROUNDMODE == FP_RND_MINF);                      \
1123           _FP_FRAC_COPY_##wc (R, Z);                                    \
1124           R##_c = Z##_c;                                                \
1125           break;                                                        \
1126                                                                         \
1127         default:                                                        \
1128           abort ();                                                     \
1129         }                                                               \
1130     done_fma: ;                                                         \
1131     }                                                                   \
1132   while (0)
1133
1134
1135 /*
1136  * Main division routine.  The input values should be cooked.
1137  */
1138
1139 #define _FP_DIV(fs, wc, R, X, Y)                                \
1140   do                                                            \
1141     {                                                           \
1142       R##_s = X##_s ^ Y##_s;                                    \
1143       R##_e = X##_e - Y##_e;                                    \
1144       switch (_FP_CLS_COMBINE (X##_c, Y##_c))                   \
1145         {                                                       \
1146         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NORMAL):    \
1147           R##_c = FP_CLS_NORMAL;                                \
1148                                                                 \
1149           _FP_DIV_MEAT_##fs (R, X, Y);                          \
1150           break;                                                \
1151                                                                 \
1152         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NAN):          \
1153           _FP_CHOOSENAN (fs, wc, R, X, Y, '/');                 \
1154           break;                                                \
1155                                                                 \
1156         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_NORMAL):       \
1157         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_INF):          \
1158         case _FP_CLS_COMBINE (FP_CLS_NAN, FP_CLS_ZERO):         \
1159           R##_s = X##_s;                                        \
1160           _FP_FRAC_COPY_##wc (R, X);                            \
1161           R##_c = X##_c;                                        \
1162           break;                                                \
1163                                                                 \
1164         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_NAN):       \
1165         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NAN):          \
1166         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NAN):         \
1167           R##_s = Y##_s;                                        \
1168           _FP_FRAC_COPY_##wc (R, Y);                            \
1169           R##_c = Y##_c;                                        \
1170           break;                                                \
1171                                                                 \
1172         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_INF):       \
1173         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_INF):         \
1174         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_NORMAL):      \
1175           R##_c = FP_CLS_ZERO;                                  \
1176           break;                                                \
1177                                                                 \
1178         case _FP_CLS_COMBINE (FP_CLS_NORMAL, FP_CLS_ZERO):      \
1179           FP_SET_EXCEPTION (FP_EX_DIVZERO);                     \
1180         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_ZERO):         \
1181         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_NORMAL):       \
1182           R##_c = FP_CLS_INF;                                   \
1183           break;                                                \
1184                                                                 \
1185         case _FP_CLS_COMBINE (FP_CLS_INF, FP_CLS_INF):          \
1186         case _FP_CLS_COMBINE (FP_CLS_ZERO, FP_CLS_ZERO):        \
1187           R##_s = _FP_NANSIGN_##fs;                             \
1188           R##_c = FP_CLS_NAN;                                   \
1189           _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs);              \
1190           FP_SET_EXCEPTION (FP_EX_INVALID);                     \
1191           break;                                                \
1192                                                                 \
1193         default:                                                \
1194           abort ();                                             \
1195         }                                                       \
1196     }                                                           \
1197   while (0)
1198
1199
1200 /*
1201  * Main differential comparison routine.  The inputs should be raw not
1202  * cooked.  The return is -1,0,1 for normal values, 2 otherwise.
1203  */
1204
1205 #define _FP_CMP(fs, wc, ret, X, Y, un)                                  \
1206   do                                                                    \
1207     {                                                                   \
1208       /* NANs are unordered */                                          \
1209       if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X))        \
1210           || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y)))    \
1211         {                                                               \
1212           ret = un;                                                     \
1213         }                                                               \
1214       else                                                              \
1215         {                                                               \
1216           int _FP_CMP_is_zero_x;                                        \
1217           int _FP_CMP_is_zero_y;                                        \
1218                                                                         \
1219           _FP_CMP_is_zero_x                                             \
1220             = (!X##_e && _FP_FRAC_ZEROP_##wc (X)) ? 1 : 0;              \
1221           _FP_CMP_is_zero_y                                             \
1222             = (!Y##_e && _FP_FRAC_ZEROP_##wc (Y)) ? 1 : 0;              \
1223                                                                         \
1224           if (_FP_CMP_is_zero_x && _FP_CMP_is_zero_y)                   \
1225             ret = 0;                                                    \
1226           else if (_FP_CMP_is_zero_x)                                   \
1227             ret = Y##_s ? 1 : -1;                                       \
1228           else if (_FP_CMP_is_zero_y)                                   \
1229             ret = X##_s ? -1 : 1;                                       \
1230           else if (X##_s != Y##_s)                                      \
1231             ret = X##_s ? -1 : 1;                                       \
1232           else if (X##_e > Y##_e)                                       \
1233             ret = X##_s ? -1 : 1;                                       \
1234           else if (X##_e < Y##_e)                                       \
1235             ret = X##_s ? 1 : -1;                                       \
1236           else if (_FP_FRAC_GT_##wc (X, Y))                             \
1237             ret = X##_s ? -1 : 1;                                       \
1238           else if (_FP_FRAC_GT_##wc (Y, X))                             \
1239             ret = X##_s ? 1 : -1;                                       \
1240           else                                                          \
1241             ret = 0;                                                    \
1242         }                                                               \
1243     }                                                                   \
1244   while (0)
1245
1246
1247 /* Simplification for strict equality.  */
1248
1249 #define _FP_CMP_EQ(fs, wc, ret, X, Y)                                   \
1250   do                                                                    \
1251     {                                                                   \
1252       /* NANs are unordered */                                          \
1253       if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X))        \
1254           || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y)))    \
1255         {                                                               \
1256           ret = 1;                                                      \
1257         }                                                               \
1258       else                                                              \
1259         {                                                               \
1260           ret = !(X##_e == Y##_e                                        \
1261                   && _FP_FRAC_EQ_##wc (X, Y)                            \
1262                   && (X##_s == Y##_s || (!X##_e && _FP_FRAC_ZEROP_##wc (X)))); \
1263         }                                                               \
1264     }                                                                   \
1265   while (0)
1266
1267 /* Version to test unordered.  */
1268
1269 #define _FP_CMP_UNORD(fs, wc, ret, X, Y)                                \
1270   do                                                                    \
1271     {                                                                   \
1272       ret = ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (X))     \
1273              || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc (Y))); \
1274     }                                                                   \
1275   while (0)
1276
1277 /*
1278  * Main square root routine.  The input value should be cooked.
1279  */
1280
1281 #define _FP_SQRT(fs, wc, R, X)                                  \
1282   do                                                            \
1283     {                                                           \
1284       _FP_FRAC_DECL_##wc (_FP_SQRT_T);                          \
1285       _FP_FRAC_DECL_##wc (_FP_SQRT_S);                          \
1286       _FP_W_TYPE _FP_SQRT_q;                                    \
1287       switch (X##_c)                                            \
1288         {                                                       \
1289         case FP_CLS_NAN:                                        \
1290           _FP_FRAC_COPY_##wc (R, X);                            \
1291           R##_s = X##_s;                                        \
1292           R##_c = FP_CLS_NAN;                                   \
1293           break;                                                \
1294         case FP_CLS_INF:                                        \
1295           if (X##_s)                                            \
1296             {                                                   \
1297               R##_s = _FP_NANSIGN_##fs;                         \
1298               R##_c = FP_CLS_NAN; /* NAN */                     \
1299               _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs);          \
1300               FP_SET_EXCEPTION (FP_EX_INVALID);                 \
1301             }                                                   \
1302           else                                                  \
1303             {                                                   \
1304               R##_s = 0;                                        \
1305               R##_c = FP_CLS_INF; /* sqrt(+inf) = +inf */       \
1306             }                                                   \
1307           break;                                                \
1308         case FP_CLS_ZERO:                                       \
1309           R##_s = X##_s;                                        \
1310           R##_c = FP_CLS_ZERO; /* sqrt(+-0) = +-0 */            \
1311           break;                                                \
1312         case FP_CLS_NORMAL:                                     \
1313           R##_s = 0;                                            \
1314           if (X##_s)                                            \
1315             {                                                   \
1316               R##_c = FP_CLS_NAN; /* NAN */                     \
1317               R##_s = _FP_NANSIGN_##fs;                         \
1318               _FP_FRAC_SET_##wc (R, _FP_NANFRAC_##fs);          \
1319               FP_SET_EXCEPTION (FP_EX_INVALID);                 \
1320               break;                                            \
1321             }                                                   \
1322           R##_c = FP_CLS_NORMAL;                                \
1323           if (X##_e & 1)                                        \
1324             _FP_FRAC_SLL_##wc (X, 1);                           \
1325           R##_e = X##_e >> 1;                                   \
1326           _FP_FRAC_SET_##wc (_FP_SQRT_S, _FP_ZEROFRAC_##wc);    \
1327           _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc);             \
1328           _FP_SQRT_q = _FP_OVERFLOW_##fs >> 1;                  \
1329           _FP_SQRT_MEAT_##wc (R, _FP_SQRT_S, _FP_SQRT_T, X,     \
1330                               _FP_SQRT_q);                      \
1331         }                                                       \
1332     }                                                           \
1333   while (0)
1334
1335 /*
1336  * Convert from FP to integer.  Input is raw.
1337  */
1338
1339 /* RSIGNED can have following values:
1340  * 0:  the number is required to be 0..(2^rsize)-1, if not, NV is set plus
1341  *     the result is either 0 or (2^rsize)-1 depending on the sign in such
1342  *     case.
1343  * 1:  the number is required to be -(2^(rsize-1))..(2^(rsize-1))-1, if not,
1344  *     NV is set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
1345  *     depending on the sign in such case.
1346  * -1: the number is required to be -(2^(rsize-1))..(2^rsize)-1, if not, NV is
1347  *     set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
1348  *     depending on the sign in such case.
1349  */
1350 #define _FP_TO_INT(fs, wc, r, X, rsize, rsigned)                        \
1351   do                                                                    \
1352     {                                                                   \
1353       if (X##_e < _FP_EXPBIAS_##fs)                                     \
1354         {                                                               \
1355           r = 0;                                                        \
1356           if (X##_e == 0)                                               \
1357             {                                                           \
1358               if (!_FP_FRAC_ZEROP_##wc (X))                             \
1359                 {                                                       \
1360                   FP_SET_EXCEPTION (FP_EX_INEXACT);                     \
1361                   FP_SET_EXCEPTION (FP_EX_DENORM);                      \
1362                 }                                                       \
1363             }                                                           \
1364           else                                                          \
1365             FP_SET_EXCEPTION (FP_EX_INEXACT);                           \
1366         }                                                               \
1367       else if (X##_e >= _FP_EXPBIAS_##fs + rsize - (rsigned > 0 || X##_s) \
1368                || (!rsigned && X##_s))                                  \
1369         {                                                               \
1370           /* Overflow or converting to the most negative integer.  */   \
1371           if (rsigned)                                                  \
1372             {                                                           \
1373               r = 1;                                                    \
1374               r <<= rsize - 1;                                          \
1375               r -= 1 - X##_s;                                           \
1376             }                                                           \
1377           else                                                          \
1378             {                                                           \
1379               r = 0;                                                    \
1380               if (!X##_s)                                               \
1381                 r = ~r;                                                 \
1382             }                                                           \
1383                                                                         \
1384           if (rsigned && X##_s && X##_e == _FP_EXPBIAS_##fs + rsize - 1) \
1385             {                                                           \
1386               /* Possibly converting to most negative integer; check the \
1387                  mantissa.  */                                          \
1388               int _FP_TO_INT_inexact = 0;                               \
1389               (void) ((_FP_FRACBITS_##fs > rsize)                       \
1390                       ? ({                                              \
1391                           _FP_FRAC_SRST_##wc (X, _FP_TO_INT_inexact,    \
1392                                               _FP_FRACBITS_##fs - rsize, \
1393                                               _FP_FRACBITS_##fs);       \
1394                           0;                                            \
1395                         })                                              \
1396                       : 0);                                             \
1397               if (!_FP_FRAC_ZEROP_##wc (X))                             \
1398                 FP_SET_EXCEPTION (FP_EX_INVALID);                       \
1399               else if (_FP_TO_INT_inexact)                              \
1400                 FP_SET_EXCEPTION (FP_EX_INEXACT);                       \
1401             }                                                           \
1402           else                                                          \
1403             FP_SET_EXCEPTION (FP_EX_INVALID);                           \
1404         }                                                               \
1405       else                                                              \
1406         {                                                               \
1407           _FP_FRAC_HIGH_RAW_##fs (X) |= _FP_IMPLBIT_##fs;               \
1408           if (X##_e >= _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs - 1)        \
1409             {                                                           \
1410               _FP_FRAC_ASSEMBLE_##wc (r, X, rsize);                     \
1411               r <<= X##_e - _FP_EXPBIAS_##fs - _FP_FRACBITS_##fs + 1;   \
1412             }                                                           \
1413           else                                                          \
1414             {                                                           \
1415               int _FP_TO_INT_inexact;                                   \
1416               _FP_FRAC_SRST_##wc (X, _FP_TO_INT_inexact,                \
1417                                   (_FP_FRACBITS_##fs + _FP_EXPBIAS_##fs - 1 \
1418                                    - X##_e),                            \
1419                                   _FP_FRACBITS_##fs);                   \
1420               if (_FP_TO_INT_inexact)                                   \
1421                 FP_SET_EXCEPTION (FP_EX_INEXACT);                       \
1422               _FP_FRAC_ASSEMBLE_##wc (r, X, rsize);                     \
1423             }                                                           \
1424           if (rsigned && X##_s)                                         \
1425             r = -r;                                                     \
1426         }                                                               \
1427     }                                                                   \
1428   while (0)
1429
1430 /* Convert integer to fp.  Output is raw.  RTYPE is unsigned even if
1431    input is signed.  */
1432 #define _FP_FROM_INT(fs, wc, X, r, rsize, rtype)                        \
1433   do                                                                    \
1434     {                                                                   \
1435       if (r)                                                            \
1436         {                                                               \
1437           rtype _FP_FROM_INT_ur;                                        \
1438                                                                         \
1439           if ((X##_s = (r < 0)))                                        \
1440             r = -(rtype) r;                                             \
1441                                                                         \
1442           _FP_FROM_INT_ur = (rtype) r;                                  \
1443           (void) ((rsize <= _FP_W_TYPE_SIZE)                            \
1444                   ? ({                                                  \
1445                       int _FP_FROM_INT_lz;                              \
1446                       __FP_CLZ (_FP_FROM_INT_lz,                        \
1447                                 (_FP_W_TYPE) _FP_FROM_INT_ur);          \
1448                       X##_e = (_FP_EXPBIAS_##fs + _FP_W_TYPE_SIZE - 1   \
1449                                - _FP_FROM_INT_lz);                      \
1450                     })                                                  \
1451                   : ((rsize <= 2 * _FP_W_TYPE_SIZE)                     \
1452                      ? ({                                               \
1453                          int _FP_FROM_INT_lz;                           \
1454                          __FP_CLZ_2 (_FP_FROM_INT_lz,                   \
1455                                      (_FP_W_TYPE) (_FP_FROM_INT_ur      \
1456                                                    >> _FP_W_TYPE_SIZE), \
1457                                      (_FP_W_TYPE) _FP_FROM_INT_ur);     \
1458                          X##_e = (_FP_EXPBIAS_##fs + 2 * _FP_W_TYPE_SIZE - 1 \
1459                                   - _FP_FROM_INT_lz);                   \
1460                        })                                               \
1461                      : (abort (), 0)));                                 \
1462                                                                         \
1463           if (rsize - 1 + _FP_EXPBIAS_##fs >= _FP_EXPMAX_##fs           \
1464               && X##_e >= _FP_EXPMAX_##fs)                              \
1465             {                                                           \
1466               /* Exponent too big; overflow to infinity.  (May also     \
1467                  happen after rounding below.)  */                      \
1468               _FP_OVERFLOW_SEMIRAW (fs, wc, X);                         \
1469               goto pack_semiraw;                                        \
1470             }                                                           \
1471                                                                         \
1472           if (rsize <= _FP_FRACBITS_##fs                                \
1473               || X##_e < _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs)          \
1474             {                                                           \
1475               /* Exactly representable; shift left.  */                 \
1476               _FP_FRAC_DISASSEMBLE_##wc (X, _FP_FROM_INT_ur, rsize);    \
1477               if (_FP_EXPBIAS_##fs + _FP_FRACBITS_##fs - 1 - X##_e > 0) \
1478                 _FP_FRAC_SLL_##wc (X, (_FP_EXPBIAS_##fs                 \
1479                                        + _FP_FRACBITS_##fs - 1 - X##_e)); \
1480             }                                                           \
1481           else                                                          \
1482             {                                                           \
1483               /* More bits in integer than in floating type; need to    \
1484                  round.  */                                             \
1485               if (_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 < X##_e)    \
1486                 _FP_FROM_INT_ur                                         \
1487                   = ((_FP_FROM_INT_ur >> (X##_e - _FP_EXPBIAS_##fs      \
1488                                           - _FP_WFRACBITS_##fs + 1))    \
1489                      | ((_FP_FROM_INT_ur                                \
1490                          << (rsize - (X##_e - _FP_EXPBIAS_##fs          \
1491                                       - _FP_WFRACBITS_##fs + 1)))       \
1492                         != 0));                                         \
1493               _FP_FRAC_DISASSEMBLE_##wc (X, _FP_FROM_INT_ur, rsize);    \
1494               if ((_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 - X##_e) > 0) \
1495                 _FP_FRAC_SLL_##wc (X, (_FP_EXPBIAS_##fs                 \
1496                                        + _FP_WFRACBITS_##fs - 1 - X##_e)); \
1497               _FP_FRAC_HIGH_##fs (X) &= ~(_FP_W_TYPE) _FP_IMPLBIT_SH_##fs; \
1498             pack_semiraw:                                               \
1499               _FP_PACK_SEMIRAW (fs, wc, X);                             \
1500             }                                                           \
1501         }                                                               \
1502       else                                                              \
1503         {                                                               \
1504           X##_s = 0;                                                    \
1505           X##_e = 0;                                                    \
1506           _FP_FRAC_SET_##wc (X, _FP_ZEROFRAC_##wc);                     \
1507         }                                                               \
1508     }                                                                   \
1509   while (0)
1510
1511
1512 /* Extend from a narrower floating-point format to a wider one.  Input
1513    and output are raw.  */
1514 #define FP_EXTEND(dfs, sfs, dwc, swc, D, S)                             \
1515   do                                                                    \
1516     {                                                                   \
1517       if (_FP_FRACBITS_##dfs < _FP_FRACBITS_##sfs                       \
1518           || (_FP_EXPMAX_##dfs - _FP_EXPBIAS_##dfs                      \
1519               < _FP_EXPMAX_##sfs - _FP_EXPBIAS_##sfs)                   \
1520           || (_FP_EXPBIAS_##dfs < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1 \
1521               && _FP_EXPBIAS_##dfs != _FP_EXPBIAS_##sfs))               \
1522         abort ();                                                       \
1523       D##_s = S##_s;                                                    \
1524       _FP_FRAC_COPY_##dwc##_##swc (D, S);                               \
1525       if (_FP_EXP_NORMAL (sfs, swc, S))                                 \
1526         {                                                               \
1527           D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs;        \
1528           _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs - _FP_FRACBITS_##sfs)); \
1529         }                                                               \
1530       else                                                              \
1531         {                                                               \
1532           if (S##_e == 0)                                               \
1533             {                                                           \
1534               if (_FP_FRAC_ZEROP_##swc (S))                             \
1535                 D##_e = 0;                                              \
1536               else if (_FP_EXPBIAS_##dfs                                \
1537                        < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1)    \
1538                 {                                                       \
1539                   FP_SET_EXCEPTION (FP_EX_DENORM);                      \
1540                   _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs            \
1541                                           - _FP_FRACBITS_##sfs));       \
1542                   D##_e = 0;                                            \
1543                 }                                                       \
1544               else                                                      \
1545                 {                                                       \
1546                   int FP_EXTEND_lz;                                     \
1547                   FP_SET_EXCEPTION (FP_EX_DENORM);                      \
1548                   _FP_FRAC_CLZ_##swc (FP_EXTEND_lz, S);                 \
1549                   _FP_FRAC_SLL_##dwc (D,                                \
1550                                       FP_EXTEND_lz + _FP_FRACBITS_##dfs \
1551                                       - _FP_FRACTBITS_##sfs);           \
1552                   D##_e = (_FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs + 1    \
1553                            + _FP_FRACXBITS_##sfs - FP_EXTEND_lz);       \
1554                 }                                                       \
1555             }                                                           \
1556           else                                                          \
1557             {                                                           \
1558               D##_e = _FP_EXPMAX_##dfs;                                 \
1559               if (!_FP_FRAC_ZEROP_##swc (S))                            \
1560                 {                                                       \
1561                   if (_FP_FRAC_SNANP (sfs, S))                          \
1562                     FP_SET_EXCEPTION (FP_EX_INVALID);                   \
1563                   _FP_FRAC_SLL_##dwc (D, (_FP_FRACBITS_##dfs            \
1564                                           - _FP_FRACBITS_##sfs));       \
1565                   _FP_SETQNAN (dfs, dwc, D);                            \
1566                 }                                                       \
1567             }                                                           \
1568         }                                                               \
1569     }                                                                   \
1570   while (0)
1571
1572 /* Truncate from a wider floating-point format to a narrower one.
1573    Input and output are semi-raw.  */
1574 #define FP_TRUNC(dfs, sfs, dwc, swc, D, S)                              \
1575   do                                                                    \
1576     {                                                                   \
1577       if (_FP_FRACBITS_##sfs < _FP_FRACBITS_##dfs                       \
1578           || (_FP_EXPBIAS_##sfs < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1 \
1579               && _FP_EXPBIAS_##sfs != _FP_EXPBIAS_##dfs))               \
1580         abort ();                                                       \
1581       D##_s = S##_s;                                                    \
1582       if (_FP_EXP_NORMAL (sfs, swc, S))                                 \
1583         {                                                               \
1584           D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs;        \
1585           if (D##_e >= _FP_EXPMAX_##dfs)                                \
1586             _FP_OVERFLOW_SEMIRAW (dfs, dwc, D);                         \
1587           else                                                          \
1588             {                                                           \
1589               if (D##_e <= 0)                                           \
1590                 {                                                       \
1591                   if (D##_e < 1 - _FP_FRACBITS_##dfs)                   \
1592                     {                                                   \
1593                       _FP_FRAC_SET_##swc (S, _FP_ZEROFRAC_##swc);       \
1594                       _FP_FRAC_LOW_##swc (S) |= 1;                      \
1595                     }                                                   \
1596                   else                                                  \
1597                     {                                                   \
1598                       _FP_FRAC_HIGH_##sfs (S) |= _FP_IMPLBIT_SH_##sfs;  \
1599                       _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs       \
1600                                               - _FP_WFRACBITS_##dfs     \
1601                                               + 1 - D##_e),             \
1602                                           _FP_WFRACBITS_##sfs);         \
1603                     }                                                   \
1604                   D##_e = 0;                                            \
1605                 }                                                       \
1606               else                                                      \
1607                 _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs             \
1608                                         - _FP_WFRACBITS_##dfs),         \
1609                                     _FP_WFRACBITS_##sfs);               \
1610               _FP_FRAC_COPY_##dwc##_##swc (D, S);                       \
1611             }                                                           \
1612         }                                                               \
1613       else                                                              \
1614         {                                                               \
1615           if (S##_e == 0)                                               \
1616             {                                                           \
1617               D##_e = 0;                                                \
1618               if (_FP_FRAC_ZEROP_##swc (S))                             \
1619                 _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc);             \
1620               else                                                      \
1621                 {                                                       \
1622                   FP_SET_EXCEPTION (FP_EX_DENORM);                      \
1623                   if (_FP_EXPBIAS_##sfs                                 \
1624                       < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1)     \
1625                     {                                                   \
1626                       _FP_FRAC_SRS_##swc (S, (_FP_WFRACBITS_##sfs       \
1627                                               - _FP_WFRACBITS_##dfs),   \
1628                                           _FP_WFRACBITS_##sfs);         \
1629                       _FP_FRAC_COPY_##dwc##_##swc (D, S);               \
1630                     }                                                   \
1631                   else                                                  \
1632                     {                                                   \
1633                       _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc);       \
1634                       _FP_FRAC_LOW_##dwc (D) |= 1;                      \
1635                     }                                                   \
1636                 }                                                       \
1637             }                                                           \
1638           else                                                          \
1639             {                                                           \
1640               D##_e = _FP_EXPMAX_##dfs;                                 \
1641               if (_FP_FRAC_ZEROP_##swc (S))                             \
1642                 _FP_FRAC_SET_##dwc (D, _FP_ZEROFRAC_##dwc);             \
1643               else                                                      \
1644                 {                                                       \
1645                   _FP_CHECK_SIGNAN_SEMIRAW (sfs, swc, S);               \
1646                   _FP_FRAC_SRL_##swc (S, (_FP_WFRACBITS_##sfs           \
1647                                           - _FP_WFRACBITS_##dfs));      \
1648                   _FP_FRAC_COPY_##dwc##_##swc (D, S);                   \
1649                   /* Semi-raw NaN must have all workbits cleared.  */   \
1650                   _FP_FRAC_LOW_##dwc (D)                                \
1651                     &= ~(_FP_W_TYPE) ((1 << _FP_WORKBITS) - 1);         \
1652                   _FP_SETQNAN_SEMIRAW (dfs, dwc, D);                    \
1653                 }                                                       \
1654             }                                                           \
1655         }                                                               \
1656     }                                                                   \
1657   while (0)
1658
1659 /*
1660  * Helper primitives.
1661  */
1662
1663 /* Count leading zeros in a word.  */
1664
1665 #ifndef __FP_CLZ
1666 /* GCC 3.4 and later provide the builtins for us.  */
1667 # define __FP_CLZ(r, x)                                                 \
1668   do                                                                    \
1669     {                                                                   \
1670       if (sizeof (_FP_W_TYPE) == sizeof (unsigned int))                 \
1671         r = __builtin_clz (x);                                          \
1672       else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long))           \
1673         r = __builtin_clzl (x);                                         \
1674       else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long long))      \
1675         r = __builtin_clzll (x);                                        \
1676       else                                                              \
1677         abort ();                                                       \
1678     }                                                                   \
1679   while (0)
1680 #endif /* ndef __FP_CLZ */
1681
1682 #define _FP_DIV_HELP_imm(q, r, n, d)            \
1683   do                                            \
1684     {                                           \
1685       q = n / d, r = n % d;                     \
1686     }                                           \
1687   while (0)
1688
1689
1690 /* A restoring bit-by-bit division primitive.  */
1691
1692 #define _FP_DIV_MEAT_N_loop(fs, wc, R, X, Y)                            \
1693   do                                                                    \
1694     {                                                                   \
1695       int _FP_DIV_MEAT_N_loop_count = _FP_WFRACBITS_##fs;               \
1696       _FP_FRAC_DECL_##wc (_FP_DIV_MEAT_N_loop_u);                       \
1697       _FP_FRAC_DECL_##wc (_FP_DIV_MEAT_N_loop_v);                       \
1698       _FP_FRAC_COPY_##wc (_FP_DIV_MEAT_N_loop_u, X);                    \
1699       _FP_FRAC_COPY_##wc (_FP_DIV_MEAT_N_loop_v, Y);                    \
1700       _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc);                         \
1701       /* Normalize _FP_DIV_MEAT_N_LOOP_U and _FP_DIV_MEAT_N_LOOP_V.  */ \
1702       _FP_FRAC_SLL_##wc (_FP_DIV_MEAT_N_loop_u, _FP_WFRACXBITS_##fs);   \
1703       _FP_FRAC_SLL_##wc (_FP_DIV_MEAT_N_loop_v, _FP_WFRACXBITS_##fs);   \
1704       /* First round.  Since the operands are normalized, either the    \
1705          first or second bit will be set in the fraction.  Produce a    \
1706          normalized result by checking which and adjusting the loop     \
1707          count and exponent accordingly.  */                            \
1708       if (_FP_FRAC_GE_1 (_FP_DIV_MEAT_N_loop_u, _FP_DIV_MEAT_N_loop_v)) \
1709         {                                                               \
1710           _FP_FRAC_SUB_##wc (_FP_DIV_MEAT_N_loop_u,                     \
1711                              _FP_DIV_MEAT_N_loop_u,                     \
1712                              _FP_DIV_MEAT_N_loop_v);                    \
1713           _FP_FRAC_LOW_##wc (R) |= 1;                                   \
1714           _FP_DIV_MEAT_N_loop_count--;                                  \
1715         }                                                               \
1716       else                                                              \
1717         R##_e--;                                                        \
1718       /* Subsequent rounds.  */                                         \
1719       do                                                                \
1720         {                                                               \
1721           int _FP_DIV_MEAT_N_loop_msb                                   \
1722             = (_FP_WS_TYPE) _FP_FRAC_HIGH_##wc (_FP_DIV_MEAT_N_loop_u) < 0; \
1723           _FP_FRAC_SLL_##wc (_FP_DIV_MEAT_N_loop_u, 1);                 \
1724           _FP_FRAC_SLL_##wc (R, 1);                                     \
1725           if (_FP_DIV_MEAT_N_loop_msb                                   \
1726               || _FP_FRAC_GE_1 (_FP_DIV_MEAT_N_loop_u,                  \
1727                                 _FP_DIV_MEAT_N_loop_v))                 \
1728             {                                                           \
1729               _FP_FRAC_SUB_##wc (_FP_DIV_MEAT_N_loop_u,                 \
1730                                  _FP_DIV_MEAT_N_loop_u,                 \
1731                                  _FP_DIV_MEAT_N_loop_v);                \
1732               _FP_FRAC_LOW_##wc (R) |= 1;                               \
1733             }                                                           \
1734         }                                                               \
1735       while (--_FP_DIV_MEAT_N_loop_count > 0);                          \
1736       /* If there's anything left in _FP_DIV_MEAT_N_LOOP_U, the result  \
1737          is inexact.  */                                                \
1738       _FP_FRAC_LOW_##wc (R)                                             \
1739         |= !_FP_FRAC_ZEROP_##wc (_FP_DIV_MEAT_N_loop_u);                \
1740     }                                                                   \
1741   while (0)
1742
1743 #define _FP_DIV_MEAT_1_loop(fs, R, X, Y)  _FP_DIV_MEAT_N_loop (fs, 1, R, X, Y)
1744 #define _FP_DIV_MEAT_2_loop(fs, R, X, Y)  _FP_DIV_MEAT_N_loop (fs, 2, R, X, Y)
1745 #define _FP_DIV_MEAT_4_loop(fs, R, X, Y)  _FP_DIV_MEAT_N_loop (fs, 4, R, X, Y)