soft-fp/op-common.h

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