soft-fp/op-common.h

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