soft-fp/op-common.h

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