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