libgcc/config/libbid/bid128_to_int32.c

   1 /* Copyright (C) 2007  Free Software Foundation, Inc.
   2
   3 This file is part of GCC.
   4
   5 GCC is free software; you can redistribute it and/or modify it under
   6 the terms of the GNU General Public License as published by the Free
   7 Software Foundation; either version 2, or (at your option) any later
   8 version.
   9
  10 In addition to the permissions in the GNU General Public License, the
  11 Free Software Foundation gives you unlimited permission to link the
  12 compiled version of this file into combinations with other programs,
  13 and to distribute those combinations without any restriction coming
  14 from the use of this file.  (The General Public License restrictions
  15 do apply in other respects; for example, they cover modification of
  16 the file, and distribution when not linked into a combine
  17 executable.)
  18
  19 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
  20 WARRANTY; without even the implied warranty of MERCHANTABILITY or
  21 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  22 for more details.
  23
  24 You should have received a copy of the GNU General Public License
  25 along with GCC; see the file COPYING.  If not, write to the Free
  26 Software Foundation, 51 Franklin Street, Fifth Floor, Boston, MA
  27 02110-1301, USA.  */
  28
  29 #include "bid_internal.h"
  30
  31 /*****************************************************************************
  32  *  BID128_to_int32_rnint
  33  ****************************************************************************/
  34
  35 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_rnint, x)
  36
  37   int res;
  38   UINT64 x_sign;
  39   UINT64 x_exp;
  40   int exp;      // unbiased exponent
  41   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
  42   UINT64 tmp64;
  43   BID_UI64DOUBLE tmp1;
  44   unsigned int x_nr_bits;
  45   int q, ind, shift;
  46   UINT128 C1, C;
  47   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
  48   UINT256 fstar;
  49   UINT256 P256;
  50
  51   // unpack x
  52   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
  53   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
  54   C1.w[1] = x.w[1] & MASK_COEFF;
  55   C1.w[0] = x.w[0];
  56
  57   // check for NaN or Infinity
  58   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
  59     // x is special
  60     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
  61       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
  62         // set invalid flag
  63         *pfpsf |= INVALID_EXCEPTION;
  64         // return Integer Indefinite
  65         res = 0x80000000;
  66       } else { // x is QNaN
  67         // set invalid flag
  68         *pfpsf |= INVALID_EXCEPTION;
  69         // return Integer Indefinite
  70         res = 0x80000000;
  71       }
  72       BID_RETURN (res);
  73     } else { // x is not a NaN, so it must be infinity
  74       if (!x_sign) { // x is +inf
  75         // set invalid flag
  76         *pfpsf |= INVALID_EXCEPTION;
  77         // return Integer Indefinite
  78         res = 0x80000000;
  79       } else { // x is -inf
  80         // set invalid flag
  81         *pfpsf |= INVALID_EXCEPTION;
  82         // return Integer Indefinite
  83         res = 0x80000000;
  84       }
  85       BID_RETURN (res);
  86     }
  87   }
  88   // check for non-canonical values (after the check for special values)
  89   if ((C1.w[1] > 0x0001ed09bead87c0ull)
  90       || (C1.w[1] == 0x0001ed09bead87c0ull
  91           && (C1.w[0] > 0x378d8e63ffffffffull))
  92       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
  93     res = 0x00000000;
  94     BID_RETURN (res);
  95   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
  96     // x is 0
  97     res = 0x00000000;
  98     BID_RETURN (res);
  99   } else { // x is not special and is not zero
 100
 101     // q = nr. of decimal digits in x
 102     //  determine first the nr. of bits in x
 103     if (C1.w[1] == 0) {
 104       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
 105         // split the 64-bit value in two 32-bit halves to avoid rounding errors
 106         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
 107           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
 108           x_nr_bits =
 109             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 110         } else { // x < 2^32
 111           tmp1.d = (double) (C1.w[0]); // exact conversion
 112           x_nr_bits =
 113             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 114         }
 115       } else { // if x < 2^53
 116         tmp1.d = (double) C1.w[0]; // exact conversion
 117         x_nr_bits =
 118           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 119       }
 120     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
 121       tmp1.d = (double) C1.w[1]; // exact conversion
 122       x_nr_bits =
 123         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 124     }
 125     q = __bid_nr_digits[x_nr_bits - 1].digits;
 126     if (q == 0) {
 127       q = __bid_nr_digits[x_nr_bits - 1].digits1;
 128       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
 129           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
 130           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
 131         q++;
 132     }
 133     exp = (x_exp >> 49) - 6176;
 134     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
 135       // set invalid flag
 136       *pfpsf |= INVALID_EXCEPTION;
 137       // return Integer Indefinite
 138       res = 0x80000000;
 139       BID_RETURN (res);
 140     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
 141       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
 142       // so x rounded to an integer may or may not fit in a signed 32-bit int
 143       // the cases that do not fit are identified here; the ones that fit
 144       // fall through and will be handled with other cases further,
 145       // under '1 <= q + exp <= 10'
 146       if (x_sign) { // if n < 0 and q + exp = 10
 147         // if n < -2^31 - 1/2 then n is too large
 148         // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31+1/2
 149         // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000005, 1<=q<=34
 150         if (q <= 11) {
 151           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
 152           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
 153           if (tmp64 > 0x500000005ull) {
 154             // set invalid flag
 155             *pfpsf |= INVALID_EXCEPTION;
 156             // return Integer Indefinite
 157             res = 0x80000000;
 158             BID_RETURN (res);
 159           }
 160           // else cases that can be rounded to a 32-bit int fall through
 161           // to '1 <= q + exp <= 10'
 162         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
 163           // 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000005 <=>
 164           // C > 0x500000005 * 10^(q-11) where 1 <= q - 11 <= 23
 165           // (scale 2^31+1/2 up)
 166           tmp64 = 0x500000005ull;
 167           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
 168             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
 169           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
 170             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
 171           }
 172           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] > C.w[0])) {
 173             // set invalid flag
 174             *pfpsf |= INVALID_EXCEPTION;
 175             // return Integer Indefinite
 176             res = 0x80000000;
 177             BID_RETURN (res);
 178           }
 179           // else cases that can be rounded to a 32-bit int fall through
 180           // to '1 <= q + exp <= 10'
 181         }
 182       } else { // if n > 0 and q + exp = 10
 183         // if n >= 2^31 - 1/2 then n is too large
 184         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31-1/2
 185         // too large if 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb, 1<=q<=34
 186         if (q <= 11) {
 187           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
 188           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
 189           if (tmp64 >= 0x4fffffffbull) {
 190             // set invalid flag
 191             *pfpsf |= INVALID_EXCEPTION;
 192             // return Integer Indefinite
 193             res = 0x80000000;
 194             BID_RETURN (res);
 195           }
 196           // else cases that can be rounded to a 32-bit int fall through
 197           // to '1 <= q + exp <= 10'
 198         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
 199           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb <=>
 200           // C >= 0x4fffffffb * 10^(q-11) where 1 <= q - 11 <= 23
 201           // (scale 2^31-1/2 up)
 202           tmp64 = 0x4fffffffbull;
 203           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
 204             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
 205           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
 206             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
 207           }
 208           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1]
 209               && C1.w[0] >= C.w[0])) {
 210             // set invalid flag
 211             *pfpsf |= INVALID_EXCEPTION;
 212             // return Integer Indefinite
 213             res = 0x80000000;
 214             BID_RETURN (res);
 215           }
 216           // else cases that can be rounded to a 32-bit int fall through
 217           // to '1 <= q + exp <= 10'
 218         }
 219       }
 220     }
 221     // n is not too large to be converted to int32: -2^31 - 1/2 < n < 2^31 - 1/2
 222     // Note: some of the cases tested for above fall through to this point
 223     if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1)
 224       // return 0
 225       res = 0x00000000;
 226       BID_RETURN (res);
 227     } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1)
 228       // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1)
 229       //   res = 0
 230       // else
 231       //   res = +/-1
 232       ind = q - 1;
 233       if (ind <= 18) { // 0 <= ind <= 18
 234         if ((C1.w[1] == 0) && (C1.w[0] <= __bid_midpoint64[ind])) {
 235           res = 0x00000000; // return 0
 236         } else if (x_sign) { // n < 0
 237           res = 0xffffffff; // return -1
 238         } else { // n > 0
 239           res = 0x00000001; // return +1
 240         }
 241       } else { // 19 <= ind <= 33
 242         if ((C1.w[1] < __bid_midpoint128[ind - 19].w[1])
 243             || ((C1.w[1] == __bid_midpoint128[ind - 19].w[1])
 244             && (C1.w[0] <= __bid_midpoint128[ind - 19].w[0]))) {
 245           res = 0x00000000; // return 0
 246         } else if (x_sign) { // n < 0
 247           res = 0xffffffff; // return -1
 248         } else { // n > 0
 249           res = 0x00000001; // return +1
 250         }
 251       }
 252     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
 253       // -2^31-1/2 <= x <= -1 or 1 <= x < 2^31-1/2 so x can be rounded
 254       // to nearest to a 32-bit signed integer
 255       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
 256         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
 257         // chop off ind digits from the lower part of C1
 258         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
 259         tmp64 = C1.w[0];
 260         if (ind <= 19) {
 261           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
 262         } else {
 263           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
 264           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
 265         }
 266         if (C1.w[0] < tmp64)
 267           C1.w[1]++;
 268         // calculate C* and f*
 269         // C* is actually floor(C*) in this case
 270         // C* and f* need shifting and masking, as shown by
 271         // __bid_shiftright128[] and __bid_maskhigh128[]
 272         // 1 <= x <= 33
 273         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
 274         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
 275         // the approximation of 10^(-x) was rounded up to 118 bits
 276         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
 277         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
 278           Cstar.w[1] = P256.w[3];
 279           Cstar.w[0] = P256.w[2];
 280           fstar.w[3] = 0;
 281           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
 282           fstar.w[1] = P256.w[1];
 283           fstar.w[0] = P256.w[0];
 284         } else { // 22 <= ind - 1 <= 33
 285           Cstar.w[1] = 0;
 286           Cstar.w[0] = P256.w[3];
 287           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
 288           fstar.w[2] = P256.w[2];
 289           fstar.w[1] = P256.w[1];
 290           fstar.w[0] = P256.w[0];
 291         }
 292         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
 293         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
 294         // if (0 < f* < 10^(-x)) then the result is a midpoint
 295         //   if floor(C*) is even then C* = floor(C*) - logical right
 296         //       shift; C* has p decimal digits, correct by Prop. 1)
 297         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
 298         //       shift; C* has p decimal digits, correct by Pr. 1)
 299         // else
 300         //   C* = floor(C*) (logical right shift; C has p decimal digits,
 301         //       correct by Property 1)
 302         // n = C* * 10^(e+x)
 303
 304         // shift right C* by Ex-128 = __bid_shiftright128[ind]
 305         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
 306         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
 307           Cstar.w[0] =
 308             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
 309           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
 310         } else { // 22 <= ind - 1 <= 33
 311           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
 312         }
 313         // if the result was a midpoint it was rounded away from zero, so
 314         // it will need a correction
 315         // check for midpoints
 316         if ((fstar.w[3] == 0) && (fstar.w[2] == 0)
 317             && (fstar.w[1] || fstar.w[0])
 318             && (fstar.w[1] < __bid_ten2mk128trunc[ind - 1].w[1]
 319                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
 320                 && fstar.w[0] <= __bid_ten2mk128trunc[ind - 1].w[0]))) {
 321           // the result is a midpoint; round to nearest
 322           if (Cstar.w[0] & 0x01) { // Cstar.w[0] is odd; MP in [EVEN, ODD]
 323             // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
 324             Cstar.w[0]--; // Cstar.w[0] is now even
 325           } // else MP in [ODD, EVEN]
 326         }
 327         if (x_sign)
 328           res = -Cstar.w[0];
 329         else
 330           res = Cstar.w[0];
 331       } else if (exp == 0) {
 332         // 1 <= q <= 10
 333         // res = +/-C (exact)
 334         if (x_sign)
 335           res = -C1.w[0];
 336         else
 337           res = C1.w[0];
 338       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
 339         // res = +/-C * 10^exp (exact)
 340         if (x_sign)
 341           res = -C1.w[0] * __bid_ten2k64[exp];
 342         else
 343           res = C1.w[0] * __bid_ten2k64[exp];
 344       }
 345     }
 346   }
 347   BID_RETURN (res);
 348 }
 349
 350 /*****************************************************************************
 351  *  BID128_to_int32_xrnint
 352  ****************************************************************************/
 353
 354 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_xrnint, x)
 355
 356   int res;
 357   UINT64 x_sign;
 358   UINT64 x_exp;
 359   int exp;      // unbiased exponent
 360   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
 361   UINT64 tmp64, tmp64A;
 362   BID_UI64DOUBLE tmp1;
 363   unsigned int x_nr_bits;
 364   int q, ind, shift;
 365   UINT128 C1, C;
 366   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
 367   UINT256 fstar;
 368   UINT256 P256;
 369
 370   // unpack x
 371   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
 372   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
 373   C1.w[1] = x.w[1] & MASK_COEFF;
 374   C1.w[0] = x.w[0];
 375
 376   // check for NaN or Infinity
 377   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
 378     // x is special
 379     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
 380       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
 381         // set invalid flag
 382         *pfpsf |= INVALID_EXCEPTION;
 383         // return Integer Indefinite
 384         res = 0x80000000;
 385       } else { // x is QNaN
 386         // set invalid flag
 387         *pfpsf |= INVALID_EXCEPTION;
 388         // return Integer Indefinite
 389         res = 0x80000000;
 390       }
 391       BID_RETURN (res);
 392     } else { // x is not a NaN, so it must be infinity
 393       if (!x_sign) { // x is +inf
 394         // set invalid flag
 395         *pfpsf |= INVALID_EXCEPTION;
 396         // return Integer Indefinite
 397         res = 0x80000000;
 398       } else { // x is -inf
 399         // set invalid flag
 400         *pfpsf |= INVALID_EXCEPTION;
 401         // return Integer Indefinite
 402         res = 0x80000000;
 403       }
 404       BID_RETURN (res);
 405     }
 406   }
 407   // check for non-canonical values (after the check for special values)
 408   if ((C1.w[1] > 0x0001ed09bead87c0ull)
 409       || (C1.w[1] == 0x0001ed09bead87c0ull
 410           && (C1.w[0] > 0x378d8e63ffffffffull))
 411       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
 412     res = 0x00000000;
 413     BID_RETURN (res);
 414   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
 415     // x is 0
 416     res = 0x00000000;
 417     BID_RETURN (res);
 418   } else { // x is not special and is not zero
 419
 420     // q = nr. of decimal digits in x
 421     //  determine first the nr. of bits in x
 422     if (C1.w[1] == 0) {
 423       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
 424         // split the 64-bit value in two 32-bit halves to avoid rounding errors
 425         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
 426           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
 427           x_nr_bits =
 428             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 429         } else { // x < 2^32
 430           tmp1.d = (double) (C1.w[0]); // exact conversion
 431           x_nr_bits =
 432             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 433         }
 434       } else { // if x < 2^53
 435         tmp1.d = (double) C1.w[0]; // exact conversion
 436         x_nr_bits =
 437           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 438       }
 439     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
 440       tmp1.d = (double) C1.w[1]; // exact conversion
 441       x_nr_bits =
 442         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 443     }
 444     q = __bid_nr_digits[x_nr_bits - 1].digits;
 445     if (q == 0) {
 446       q = __bid_nr_digits[x_nr_bits - 1].digits1;
 447       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
 448           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
 449           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
 450         q++;
 451     }
 452     exp = (x_exp >> 49) - 6176;
 453     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
 454       // set invalid flag
 455       *pfpsf |= INVALID_EXCEPTION;
 456       // return Integer Indefinite
 457       res = 0x80000000;
 458       BID_RETURN (res);
 459     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
 460       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
 461       // so x rounded to an integer may or may not fit in a signed 32-bit int
 462       // the cases that do not fit are identified here; the ones that fit
 463       // fall through and will be handled with other cases further,
 464       // under '1 <= q + exp <= 10'
 465       if (x_sign) { // if n < 0 and q + exp = 10
 466         // if n < -2^31 - 1/2 then n is too large
 467         // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31+1/2
 468         // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000005, 1<=q<=34
 469         if (q <= 11) {
 470           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
 471           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
 472           if (tmp64 > 0x500000005ull) {
 473             // set invalid flag
 474             *pfpsf |= INVALID_EXCEPTION;
 475             // return Integer Indefinite
 476             res = 0x80000000;
 477             BID_RETURN (res);
 478           }
 479           // else cases that can be rounded to a 32-bit int fall through
 480           // to '1 <= q + exp <= 10'
 481         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
 482           // 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000005 <=>
 483           // C > 0x500000005 * 10^(q-11) where 1 <= q - 11 <= 23
 484           // (scale 2^31+1/2 up)
 485           tmp64 = 0x500000005ull;
 486           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
 487             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
 488           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
 489             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
 490           }
 491           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] > C.w[0])) {
 492             // set invalid flag
 493             *pfpsf |= INVALID_EXCEPTION;
 494             // return Integer Indefinite
 495             res = 0x80000000;
 496             BID_RETURN (res);
 497           }
 498           // else cases that can be rounded to a 32-bit int fall through
 499           // to '1 <= q + exp <= 10'
 500         }
 501       } else { // if n > 0 and q + exp = 10
 502         // if n >= 2^31 - 1/2 then n is too large
 503         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31-1/2
 504         // too large if 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb, 1<=q<=34
 505         if (q <= 11) {
 506           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
 507           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
 508           if (tmp64 >= 0x4fffffffbull) {
 509             // set invalid flag
 510             *pfpsf |= INVALID_EXCEPTION;
 511             // return Integer Indefinite
 512             res = 0x80000000;
 513             BID_RETURN (res);
 514           }
 515           // else cases that can be rounded to a 32-bit int fall through
 516           // to '1 <= q + exp <= 10'
 517         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
 518           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb <=>
 519           // C >= 0x4fffffffb * 10^(q-11) where 1 <= q - 11 <= 23
 520           // (scale 2^31-1/2 up)
 521           tmp64 = 0x4fffffffbull;
 522           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
 523             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
 524           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
 525             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
 526           }
 527           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
 528             // set invalid flag
 529             *pfpsf |= INVALID_EXCEPTION;
 530             // return Integer Indefinite
 531             res = 0x80000000;
 532             BID_RETURN (res);
 533           }
 534           // else cases that can be rounded to a 32-bit int fall through
 535           // to '1 <= q + exp <= 10'
 536         }
 537       }
 538     }
 539     // n is not too large to be converted to int32: -2^31 - 1/2 < n < 2^31 - 1/2
 540     // Note: some of the cases tested for above fall through to this point
 541     if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1)
 542       // set inexact flag
 543       *pfpsf |= INEXACT_EXCEPTION;
 544       // return 0
 545       res = 0x00000000;
 546       BID_RETURN (res);
 547     } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1)
 548       // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1)
 549       //   res = 0
 550       // else
 551       //   res = +/-1
 552       ind = q - 1;
 553       if (ind <= 18) { // 0 <= ind <= 18
 554         if ((C1.w[1] == 0) && (C1.w[0] <= __bid_midpoint64[ind])) {
 555           res = 0x00000000; // return 0
 556         } else if (x_sign) { // n < 0
 557           res = 0xffffffff; // return -1
 558         } else { // n > 0
 559           res = 0x00000001; // return +1
 560         }
 561       } else { // 19 <= ind <= 33
 562         if ((C1.w[1] < __bid_midpoint128[ind - 19].w[1])
 563             || ((C1.w[1] == __bid_midpoint128[ind - 19].w[1])
 564             && (C1.w[0] <= __bid_midpoint128[ind - 19].w[0]))) {
 565           res = 0x00000000; // return 0
 566         } else if (x_sign) { // n < 0
 567           res = 0xffffffff; // return -1
 568         } else { // n > 0
 569           res = 0x00000001; // return +1
 570         }
 571       }
 572       // set inexact flag
 573       *pfpsf |= INEXACT_EXCEPTION;
 574     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
 575       // -2^31-1/2 <= x <= -1 or 1 <= x < 2^31-1/2 so x can be rounded
 576       // to nearest to a 32-bit signed integer
 577       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
 578         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
 579         // chop off ind digits from the lower part of C1
 580         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
 581         tmp64 = C1.w[0];
 582         if (ind <= 19) {
 583           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
 584         } else {
 585           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
 586           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
 587         }
 588         if (C1.w[0] < tmp64)
 589           C1.w[1]++;
 590         // calculate C* and f*
 591         // C* is actually floor(C*) in this case
 592         // C* and f* need shifting and masking, as shown by
 593         // __bid_shiftright128[] and __bid_maskhigh128[]
 594         // 1 <= x <= 33
 595         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
 596         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
 597         // the approximation of 10^(-x) was rounded up to 118 bits
 598         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
 599         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
 600           Cstar.w[1] = P256.w[3];
 601           Cstar.w[0] = P256.w[2];
 602           fstar.w[3] = 0;
 603           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
 604           fstar.w[1] = P256.w[1];
 605           fstar.w[0] = P256.w[0];
 606         } else { // 22 <= ind - 1 <= 33
 607           Cstar.w[1] = 0;
 608           Cstar.w[0] = P256.w[3];
 609           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
 610           fstar.w[2] = P256.w[2];
 611           fstar.w[1] = P256.w[1];
 612           fstar.w[0] = P256.w[0];
 613         }
 614         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
 615         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
 616         // if (0 < f* < 10^(-x)) then the result is a midpoint
 617         //   if floor(C*) is even then C* = floor(C*) - logical right
 618         //       shift; C* has p decimal digits, correct by Prop. 1)
 619         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
 620         //       shift; C* has p decimal digits, correct by Pr. 1)
 621         // else
 622         //   C* = floor(C*) (logical right shift; C has p decimal digits,
 623         //       correct by Property 1)
 624         // n = C* * 10^(e+x)
 625
 626         // shift right C* by Ex-128 = __bid_shiftright128[ind]
 627         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
 628         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
 629           Cstar.w[0] =
 630             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
 631           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
 632         } else { // 22 <= ind - 1 <= 33
 633           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
 634         }
 635         // determine inexactness of the rounding of C*
 636         // if (0 < f* - 1/2 < 10^(-x)) then
 637         //   the result is exact
 638         // else // if (f* - 1/2 > T*) then
 639         //   the result is inexact
 640         if (ind - 1 <= 2) {
 641           if (fstar.w[1] > 0x8000000000000000ull ||
 642               (fstar.w[1] == 0x8000000000000000ull &&
 643               fstar.w[0] > 0x0ull)) { // f* > 1/2 and the result may be exact
 644             tmp64 = fstar.w[1] - 0x8000000000000000ull; // f* - 1/2
 645             if ((tmp64 > __bid_ten2mk128trunc[ind - 1].w[1]
 646                  || (tmp64 == __bid_ten2mk128trunc[ind - 1].w[1]
 647                  && fstar.w[0] >= __bid_ten2mk128trunc[ind - 1].w[0]))) {
 648               // set the inexact flag
 649               *pfpsf |= INEXACT_EXCEPTION;
 650             } // else the result is exact
 651           } else { // the result is inexact; f2* <= 1/2
 652             // set the inexact flag
 653             *pfpsf |= INEXACT_EXCEPTION;
 654           }
 655         } else if (ind - 1 <= 21) { // if 3 <= ind <= 21
 656           if (fstar.w[3] > 0x0 || (fstar.w[3] == 0x0
 657               && fstar.w[2] > __bid_one_half128[ind - 1]) || (fstar.w[3] == 0x0
 658               && fstar.w[2] == __bid_one_half128[ind - 1] && (fstar.w[1]
 659                                                        || fstar.w[0]))) {
 660             // f2* > 1/2 and the result may be exact
 661             // Calculate f2* - 1/2
 662             tmp64 = fstar.w[2] - __bid_one_half128[ind - 1];
 663             tmp64A = fstar.w[3];
 664             if (tmp64 > fstar.w[2])
 665               tmp64A--;
 666             if (tmp64A || tmp64
 667                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
 668                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
 669                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
 670               // set the inexact flag
 671               *pfpsf |= INEXACT_EXCEPTION;
 672             } // else the result is exact
 673           } else { // the result is inexact; f2* <= 1/2
 674             // set the inexact flag
 675             *pfpsf |= INEXACT_EXCEPTION;
 676           }
 677         } else { // if 22 <= ind <= 33
 678           if (fstar.w[3] > __bid_one_half128[ind - 1]
 679               || (fstar.w[3] == __bid_one_half128[ind - 1] && (fstar.w[2]
 680               || fstar.w[1] || fstar.w[0]))) {
 681             // f2* > 1/2 and the result may be exact
 682             // Calculate f2* - 1/2
 683             tmp64 = fstar.w[3] - __bid_one_half128[ind - 1];
 684             if (tmp64 || fstar.w[2]
 685                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
 686                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
 687                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
 688               // set the inexact flag
 689               *pfpsf |= INEXACT_EXCEPTION;
 690             } // else the result is exact
 691           } else { // the result is inexact; f2* <= 1/2
 692             // set the inexact flag
 693             *pfpsf |= INEXACT_EXCEPTION;
 694           }
 695         }
 696         // if the result was a midpoint it was rounded away from zero, so
 697         // it will need a correction
 698         // check for midpoints
 699         if ((fstar.w[3] == 0) && (fstar.w[2] == 0)
 700             && (fstar.w[1] || fstar.w[0])
 701             && (fstar.w[1] < __bid_ten2mk128trunc[ind - 1].w[1]
 702                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
 703                 && fstar.w[0] <= __bid_ten2mk128trunc[ind - 1].w[0]))) {
 704           // the result is a midpoint; round to nearest
 705           if (Cstar.w[0] & 0x01) { // Cstar.w[0] is odd; MP in [EVEN, ODD]
 706             // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
 707             Cstar.w[0]--; // Cstar.w[0] is now even
 708           } // else MP in [ODD, EVEN]
 709         }
 710         if (x_sign)
 711           res = -Cstar.w[0];
 712         else
 713           res = Cstar.w[0];
 714       } else if (exp == 0) {
 715         // 1 <= q <= 10
 716         // res = +/-C (exact)
 717         if (x_sign)
 718           res = -C1.w[0];
 719         else
 720           res = C1.w[0];
 721       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
 722         // res = +/-C * 10^exp (exact)
 723         if (x_sign)
 724           res = -C1.w[0] * __bid_ten2k64[exp];
 725         else
 726           res = C1.w[0] * __bid_ten2k64[exp];
 727       }
 728     }
 729   }
 730   BID_RETURN (res);
 731 }
 732
 733 /*****************************************************************************
 734  *  BID128_to_int32_floor
 735  ****************************************************************************/
 736
 737 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_floor, x)
 738
 739   int res;
 740   UINT64 x_sign;
 741   UINT64 x_exp;
 742   int exp;      // unbiased exponent
 743   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
 744   UINT64 tmp64, tmp64A;
 745   BID_UI64DOUBLE tmp1;
 746   unsigned int x_nr_bits;
 747   int q, ind, shift;
 748   UINT128 C1, C;
 749   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
 750   UINT256 fstar;
 751   UINT256 P256;
 752   int is_inexact_lt_midpoint = 0;
 753   int is_inexact_gt_midpoint = 0;
 754   int is_midpoint_lt_even = 0;
 755   int is_midpoint_gt_even = 0;
 756
 757   // unpack x
 758   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
 759   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
 760   C1.w[1] = x.w[1] & MASK_COEFF;
 761   C1.w[0] = x.w[0];
 762
 763   // check for NaN or Infinity
 764   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
 765     // x is special
 766     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
 767       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
 768         // set invalid flag
 769         *pfpsf |= INVALID_EXCEPTION;
 770         // return Integer Indefinite
 771         res = 0x80000000;
 772       } else { // x is QNaN
 773         // set invalid flag
 774         *pfpsf |= INVALID_EXCEPTION;
 775         // return Integer Indefinite
 776         res = 0x80000000;
 777       }
 778       BID_RETURN (res);
 779     } else { // x is not a NaN, so it must be infinity
 780       if (!x_sign) { // x is +inf
 781         // set invalid flag
 782         *pfpsf |= INVALID_EXCEPTION;
 783         // return Integer Indefinite
 784         res = 0x80000000;
 785       } else { // x is -inf
 786         // set invalid flag
 787         *pfpsf |= INVALID_EXCEPTION;
 788         // return Integer Indefinite
 789         res = 0x80000000;
 790       }
 791       BID_RETURN (res);
 792     }
 793   }
 794   // check for non-canonical values (after the check for special values)
 795   if ((C1.w[1] > 0x0001ed09bead87c0ull)
 796       || (C1.w[1] == 0x0001ed09bead87c0ull
 797           && (C1.w[0] > 0x378d8e63ffffffffull))
 798       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
 799     res = 0x00000000;
 800     BID_RETURN (res);
 801   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
 802     // x is 0
 803     res = 0x00000000;
 804     BID_RETURN (res);
 805   } else { // x is not special and is not zero
 806
 807     // q = nr. of decimal digits in x
 808     //  determine first the nr. of bits in x
 809     if (C1.w[1] == 0) {
 810       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
 811         // split the 64-bit value in two 32-bit halves to avoid rounding errors
 812         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
 813           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
 814           x_nr_bits =
 815             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 816         } else { // x < 2^32
 817           tmp1.d = (double) (C1.w[0]); // exact conversion
 818           x_nr_bits =
 819             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 820         }
 821       } else { // if x < 2^53
 822         tmp1.d = (double) C1.w[0]; // exact conversion
 823         x_nr_bits =
 824           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 825       }
 826     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
 827       tmp1.d = (double) C1.w[1]; // exact conversion
 828       x_nr_bits =
 829         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
 830     }
 831     q = __bid_nr_digits[x_nr_bits - 1].digits;
 832     if (q == 0) {
 833       q = __bid_nr_digits[x_nr_bits - 1].digits1;
 834       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
 835           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
 836           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
 837         q++;
 838     }
 839     exp = (x_exp >> 49) - 6176;
 840     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
 841       // set invalid flag
 842       *pfpsf |= INVALID_EXCEPTION;
 843       // return Integer Indefinite
 844       res = 0x80000000;
 845       BID_RETURN (res);
 846     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
 847       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
 848       // so x rounded to an integer may or may not fit in a signed 32-bit int
 849       // the cases that do not fit are identified here; the ones that fit
 850       // fall through and will be handled with other cases further,
 851       // under '1 <= q + exp <= 10'
 852       if (x_sign) { // if n < 0 and q + exp = 10
 853         // if n < -2^31 then n is too large
 854         // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31
 855         // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000000, 1<=q<=34
 856         if (q <= 11) {
 857           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
 858           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
 859           if (tmp64 > 0x500000000ull) {
 860             // set invalid flag
 861             *pfpsf |= INVALID_EXCEPTION;
 862             // return Integer Indefinite
 863             res = 0x80000000;
 864             BID_RETURN (res);
 865           }
 866           // else cases that can be rounded to a 32-bit int fall through
 867           // to '1 <= q + exp <= 10'
 868         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
 869           // 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000000 <=>
 870           // C > 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 23
 871           // (scale 2^31 up)
 872           tmp64 = 0x500000000ull;
 873           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
 874             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
 875           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
 876             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
 877           }
 878           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] > C.w[0])) {
 879             // set invalid flag
 880             *pfpsf |= INVALID_EXCEPTION;
 881             // return Integer Indefinite
 882             res = 0x80000000;
 883             BID_RETURN (res);
 884           }
 885           // else cases that can be rounded to a 32-bit int fall through
 886           // to '1 <= q + exp <= 10'
 887         }
 888       } else { // if n > 0 and q + exp = 10
 889         // if n >= 2^31 then n is too large
 890         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31
 891         // too large if 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000, 1<=q<=34
 892         if (q <= 11) {
 893           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
 894           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
 895           if (tmp64 >= 0x500000000ull) {
 896             // set invalid flag
 897             *pfpsf |= INVALID_EXCEPTION;
 898             // return Integer Indefinite
 899             res = 0x80000000;
 900             BID_RETURN (res);
 901           }
 902           // else cases that can be rounded to a 32-bit int fall through
 903           // to '1 <= q + exp <= 10'
 904         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
 905           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000 <=>
 906           // C >= 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 23
 907           // (scale 2^31 up)
 908           tmp64 = 0x500000000ull;
 909           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
 910             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
 911           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
 912             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
 913           }
 914           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
 915             // set invalid flag
 916             *pfpsf |= INVALID_EXCEPTION;
 917             // return Integer Indefinite
 918             res = 0x80000000;
 919             BID_RETURN (res);
 920           }
 921           // else cases that can be rounded to a 32-bit int fall through
 922           // to '1 <= q + exp <= 10'
 923         }
 924       }
 925     }
 926     // n is not too large to be converted to int32: -2^31 <= n < 2^31
 927     // Note: some of the cases tested for above fall through to this point
 928     if ((q + exp) <= 0) {
 929       // n = +/-0.0...c(0)c(1)...c(q-1) or n = +/-0.c(0)c(1)...c(q-1)
 930       // return 0
 931       if (x_sign)
 932         res = 0xffffffff;
 933       else
 934         res = 0x00000000;
 935       BID_RETURN (res);
 936     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
 937       // -2^31 <= x <= -1 or 1 <= x < 2^31 so x can be rounded
 938       // toward negative infinity to a 32-bit signed integer
 939       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
 940         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
 941         // chop off ind digits from the lower part of C1
 942         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
 943         tmp64 = C1.w[0];
 944         if (ind <= 19) {
 945           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
 946         } else {
 947           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
 948           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
 949         }
 950         if (C1.w[0] < tmp64)
 951           C1.w[1]++;
 952         // calculate C* and f*
 953         // C* is actually floor(C*) in this case
 954         // C* and f* need shifting and masking, as shown by
 955         // __bid_shiftright128[] and __bid_maskhigh128[]
 956         // 1 <= x <= 33
 957         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
 958         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
 959         // the approximation of 10^(-x) was rounded up to 118 bits
 960         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
 961         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
 962           Cstar.w[1] = P256.w[3];
 963           Cstar.w[0] = P256.w[2];
 964           fstar.w[3] = 0;
 965           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
 966           fstar.w[1] = P256.w[1];
 967           fstar.w[0] = P256.w[0];
 968         } else { // 22 <= ind - 1 <= 33
 969           Cstar.w[1] = 0;
 970           Cstar.w[0] = P256.w[3];
 971           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
 972           fstar.w[2] = P256.w[2];
 973           fstar.w[1] = P256.w[1];
 974           fstar.w[0] = P256.w[0];
 975         }
 976         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
 977         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
 978         // if (0 < f* < 10^(-x)) then the result is a midpoint
 979         //   if floor(C*) is even then C* = floor(C*) - logical right
 980         //       shift; C* has p decimal digits, correct by Prop. 1)
 981         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
 982         //       shift; C* has p decimal digits, correct by Pr. 1)
 983         // else
 984         //   C* = floor(C*) (logical right shift; C has p decimal digits,
 985         //       correct by Property 1)
 986         // n = C* * 10^(e+x)
 987
 988         // shift right C* by Ex-128 = __bid_shiftright128[ind]
 989         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
 990         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
 991           Cstar.w[0] =
 992             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
 993           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
 994         } else { // 22 <= ind - 1 <= 33
 995           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
 996         }
 997         // determine inexactness of the rounding of C*
 998         // if (0 < f* - 1/2 < 10^(-x)) then
 999         //   the result is exact
1000         // else // if (f* - 1/2 > T*) then
1001         //   the result is inexact
1002         if (ind - 1 <= 2) {
1003           if (fstar.w[1] > 0x8000000000000000ull ||
1004               (fstar.w[1] == 0x8000000000000000ull &&
1005               fstar.w[0] > 0x0ull)) { // f* > 1/2 and the result may be exact
1006             tmp64 = fstar.w[1] - 0x8000000000000000ull; // f* - 1/2
1007             if ((tmp64 > __bid_ten2mk128trunc[ind - 1].w[1]
1008                  || (tmp64 == __bid_ten2mk128trunc[ind - 1].w[1]
1009                  && fstar.w[0] >= __bid_ten2mk128trunc[ind - 1].w[0]))) {
1010               is_inexact_lt_midpoint = 1;
1011             } // else the result is exact
1012           } else { // the result is inexact; f2* <= 1/2
1013             is_inexact_gt_midpoint = 1;
1014           }
1015         } else if (ind - 1 <= 21) { // if 3 <= ind <= 21
1016           if (fstar.w[3] > 0x0 || (fstar.w[3] == 0x0
1017               && fstar.w[2] > __bid_one_half128[ind - 1]) || (fstar.w[3] == 0x0
1018               && fstar.w[2] == __bid_one_half128[ind - 1] && (fstar.w[1]
1019                                                        || fstar.w[0]))) {
1020             // f2* > 1/2 and the result may be exact
1021             // Calculate f2* - 1/2
1022             tmp64 = fstar.w[2] - __bid_one_half128[ind - 1];
1023             tmp64A = fstar.w[3];
1024             if (tmp64 > fstar.w[2])
1025               tmp64A--;
1026             if (tmp64A || tmp64
1027                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
1028                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1029                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
1030               is_inexact_lt_midpoint = 1;
1031             } // else the result is exact
1032           } else { // the result is inexact; f2* <= 1/2
1033             is_inexact_gt_midpoint = 1;
1034           }
1035         } else { // if 22 <= ind <= 33
1036           if (fstar.w[3] > __bid_one_half128[ind - 1]
1037               || (fstar.w[3] == __bid_one_half128[ind - 1] && (fstar.w[2]
1038               || fstar.w[1] || fstar.w[0]))) {
1039             // f2* > 1/2 and the result may be exact
1040             // Calculate f2* - 1/2
1041             tmp64 = fstar.w[3] - __bid_one_half128[ind - 1];
1042             if (tmp64 || fstar.w[2]
1043                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
1044                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1045                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
1046               is_inexact_lt_midpoint = 1;
1047             } // else the result is exact
1048           } else { // the result is inexact; f2* <= 1/2
1049             is_inexact_gt_midpoint = 1;
1050           }
1051         }
1052
1053         // if the result was a midpoint it was rounded away from zero, so
1054         // it will need a correction
1055         // check for midpoints
1056         if ((fstar.w[3] == 0) && (fstar.w[2] == 0)
1057             && (fstar.w[1] || fstar.w[0])
1058             && (fstar.w[1] < __bid_ten2mk128trunc[ind - 1].w[1]
1059                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1060                 && fstar.w[0] <= __bid_ten2mk128trunc[ind - 1].w[0]))) {
1061           // the result is a midpoint; round to nearest
1062           if (Cstar.w[0] & 0x01) { // Cstar.w[0] is odd; MP in [EVEN, ODD]
1063             // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
1064             Cstar.w[0]--; // Cstar.w[0] is now even
1065             is_midpoint_gt_even = 1;
1066             is_inexact_lt_midpoint = 0;
1067             is_inexact_gt_midpoint = 0;
1068           } else { // else MP in [ODD, EVEN]
1069             is_midpoint_lt_even = 1;
1070             is_inexact_lt_midpoint = 0;
1071             is_inexact_gt_midpoint = 0;
1072           }
1073         }
1074         // general correction for RM
1075         if (x_sign && (is_midpoint_gt_even || is_inexact_lt_midpoint)) {
1076           Cstar.w[0] = Cstar.w[0] + 1;
1077         } else if (!x_sign
1078                    && (is_midpoint_lt_even || is_inexact_gt_midpoint)) {
1079           Cstar.w[0] = Cstar.w[0] - 1;
1080         } else {
1081           ; // the result is already correct
1082         }
1083         if (x_sign)
1084           res = -Cstar.w[0];
1085         else
1086           res = Cstar.w[0];
1087       } else if (exp == 0) {
1088         // 1 <= q <= 10
1089         // res = +/-C (exact)
1090         if (x_sign)
1091           res = -C1.w[0];
1092         else
1093           res = C1.w[0];
1094       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
1095         // res = +/-C * 10^exp (exact)
1096         if (x_sign)
1097           res = -C1.w[0] * __bid_ten2k64[exp];
1098         else
1099           res = C1.w[0] * __bid_ten2k64[exp];
1100       }
1101     }
1102   }
1103   BID_RETURN (res);
1104 }
1105
1106
1107 /*****************************************************************************
1108  *  BID128_to_int32_xfloor
1109  ****************************************************************************/
1110
1111 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_xfloor, x)
1112
1113   int res;
1114   UINT64 x_sign;
1115   UINT64 x_exp;
1116   int exp;      // unbiased exponent
1117   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
1118   UINT64 tmp64, tmp64A;
1119   BID_UI64DOUBLE tmp1;
1120   unsigned int x_nr_bits;
1121   int q, ind, shift;
1122   UINT128 C1, C;
1123   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
1124   UINT256 fstar;
1125   UINT256 P256;
1126   int is_inexact_lt_midpoint = 0;
1127   int is_inexact_gt_midpoint = 0;
1128   int is_midpoint_lt_even = 0;
1129   int is_midpoint_gt_even = 0;
1130
1131   // unpack x
1132   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
1133   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
1134   C1.w[1] = x.w[1] & MASK_COEFF;
1135   C1.w[0] = x.w[0];
1136
1137   // check for NaN or Infinity
1138   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
1139     // x is special
1140     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
1141       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
1142         // set invalid flag
1143         *pfpsf |= INVALID_EXCEPTION;
1144         // return Integer Indefinite
1145         res = 0x80000000;
1146       } else { // x is QNaN
1147         // set invalid flag
1148         *pfpsf |= INVALID_EXCEPTION;
1149         // return Integer Indefinite
1150         res = 0x80000000;
1151       }
1152       BID_RETURN (res);
1153     } else { // x is not a NaN, so it must be infinity
1154       if (!x_sign) { // x is +inf
1155         // set invalid flag
1156         *pfpsf |= INVALID_EXCEPTION;
1157         // return Integer Indefinite
1158         res = 0x80000000;
1159       } else { // x is -inf
1160         // set invalid flag
1161         *pfpsf |= INVALID_EXCEPTION;
1162         // return Integer Indefinite
1163         res = 0x80000000;
1164       }
1165       BID_RETURN (res);
1166     }
1167   }
1168   // check for non-canonical values (after the check for special values)
1169   if ((C1.w[1] > 0x0001ed09bead87c0ull)
1170       || (C1.w[1] == 0x0001ed09bead87c0ull
1171           && (C1.w[0] > 0x378d8e63ffffffffull))
1172       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
1173     res = 0x00000000;
1174     BID_RETURN (res);
1175   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
1176     // x is 0
1177     res = 0x00000000;
1178     BID_RETURN (res);
1179   } else { // x is not special and is not zero
1180
1181     // q = nr. of decimal digits in x
1182     //  determine first the nr. of bits in x
1183     if (C1.w[1] == 0) {
1184       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
1185         // split the 64-bit value in two 32-bit halves to avoid rounding errors
1186         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
1187           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
1188           x_nr_bits =
1189             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1190         } else { // x < 2^32
1191           tmp1.d = (double) (C1.w[0]); // exact conversion
1192           x_nr_bits =
1193             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1194         }
1195       } else { // if x < 2^53
1196         tmp1.d = (double) C1.w[0]; // exact conversion
1197         x_nr_bits =
1198           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1199       }
1200     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
1201       tmp1.d = (double) C1.w[1]; // exact conversion
1202       x_nr_bits =
1203         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1204     }
1205     q = __bid_nr_digits[x_nr_bits - 1].digits;
1206     if (q == 0) {
1207       q = __bid_nr_digits[x_nr_bits - 1].digits1;
1208       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
1209           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
1210           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
1211         q++;
1212     }
1213     exp = (x_exp >> 49) - 6176;
1214     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
1215       // set invalid flag
1216       *pfpsf |= INVALID_EXCEPTION;
1217       // return Integer Indefinite
1218       res = 0x80000000;
1219       BID_RETURN (res);
1220     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
1221       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
1222       // so x rounded to an integer may or may not fit in a signed 32-bit int
1223       // the cases that do not fit are identified here; the ones that fit
1224       // fall through and will be handled with other cases further,
1225       // under '1 <= q + exp <= 10'
1226       if (x_sign) { // if n < 0 and q + exp = 10
1227         // if n < -2^31 then n is too large
1228         // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31
1229         // <=> 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000000, 1<=q<=34
1230         if (q <= 11) {
1231           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
1232           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
1233           if (tmp64 > 0x500000000ull) {
1234             // set invalid flag
1235             *pfpsf |= INVALID_EXCEPTION;
1236             // return Integer Indefinite
1237             res = 0x80000000;
1238             BID_RETURN (res);
1239           }
1240           // else cases that can be rounded to a 32-bit int fall through
1241           // to '1 <= q + exp <= 10'
1242         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
1243           // 0.c(0)c(1)...c(q-1) * 10^11 > 0x500000000 <=>
1244           // C > 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 23
1245           // (scale 2^31 up)
1246           tmp64 = 0x500000000ull;
1247           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
1248             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
1249           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
1250             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
1251           }
1252           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] > C.w[0])) {
1253             // set invalid flag
1254             *pfpsf |= INVALID_EXCEPTION;
1255             // return Integer Indefinite
1256             res = 0x80000000;
1257             BID_RETURN (res);
1258           }
1259           // else cases that can be rounded to a 32-bit int fall through
1260           // to '1 <= q + exp <= 10'
1261         }
1262       } else { // if n > 0 and q + exp = 10
1263         // if n >= 2^31 then n is too large
1264         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31
1265         // too large if 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000, 1<=q<=34
1266         if (q <= 11) {
1267           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
1268           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
1269           if (tmp64 >= 0x500000000ull) {
1270             // set invalid flag
1271             *pfpsf |= INVALID_EXCEPTION;
1272             // return Integer Indefinite
1273             res = 0x80000000;
1274             BID_RETURN (res);
1275           }
1276           // else cases that can be rounded to a 32-bit int fall through
1277           // to '1 <= q + exp <= 10'
1278         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
1279           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000 <=>
1280           // C >= 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 23
1281           // (scale 2^31 up)
1282           tmp64 = 0x500000000ull;
1283           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
1284             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
1285           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
1286             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
1287           }
1288           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
1289             // set invalid flag
1290             *pfpsf |= INVALID_EXCEPTION;
1291             // return Integer Indefinite
1292             res = 0x80000000;
1293             BID_RETURN (res);
1294           }
1295           // else cases that can be rounded to a 32-bit int fall through
1296           // to '1 <= q + exp <= 10'
1297         }
1298       }
1299     }
1300     // n is not too large to be converted to int32: -2^31 <= n < 2^31
1301     // Note: some of the cases tested for above fall through to this point
1302     if ((q + exp) <= 0) {
1303       // n = +/-0.0...c(0)c(1)...c(q-1) or n = +/-0.c(0)c(1)...c(q-1)
1304       // set inexact flag
1305       *pfpsf |= INEXACT_EXCEPTION;
1306       // return 0
1307       if (x_sign)
1308         res = 0xffffffff;
1309       else
1310         res = 0x00000000;
1311       BID_RETURN (res);
1312     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
1313       // -2^31 <= x <= -1 or 1 <= x < 2^31 so x can be rounded
1314       // toward negative infinity to a 32-bit signed integer
1315       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
1316         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
1317         // chop off ind digits from the lower part of C1
1318         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
1319         tmp64 = C1.w[0];
1320         if (ind <= 19) {
1321           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
1322         } else {
1323           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
1324           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
1325         }
1326         if (C1.w[0] < tmp64)
1327           C1.w[1]++;
1328         // calculate C* and f*
1329         // C* is actually floor(C*) in this case
1330         // C* and f* need shifting and masking, as shown by
1331         // __bid_shiftright128[] and __bid_maskhigh128[]
1332         // 1 <= x <= 33
1333         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
1334         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
1335         // the approximation of 10^(-x) was rounded up to 118 bits
1336         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
1337         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
1338           Cstar.w[1] = P256.w[3];
1339           Cstar.w[0] = P256.w[2];
1340           fstar.w[3] = 0;
1341           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
1342           fstar.w[1] = P256.w[1];
1343           fstar.w[0] = P256.w[0];
1344         } else { // 22 <= ind - 1 <= 33
1345           Cstar.w[1] = 0;
1346           Cstar.w[0] = P256.w[3];
1347           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
1348           fstar.w[2] = P256.w[2];
1349           fstar.w[1] = P256.w[1];
1350           fstar.w[0] = P256.w[0];
1351         }
1352         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
1353         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
1354         // if (0 < f* < 10^(-x)) then the result is a midpoint
1355         //   if floor(C*) is even then C* = floor(C*) - logical right
1356         //       shift; C* has p decimal digits, correct by Prop. 1)
1357         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
1358         //       shift; C* has p decimal digits, correct by Pr. 1)
1359         // else
1360         //   C* = floor(C*) (logical right shift; C has p decimal digits,
1361         //       correct by Property 1)
1362         // n = C* * 10^(e+x)
1363
1364         // shift right C* by Ex-128 = __bid_shiftright128[ind]
1365         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
1366         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
1367           Cstar.w[0] =
1368             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
1369           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
1370         } else { // 22 <= ind - 1 <= 33
1371           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
1372         }
1373         // determine inexactness of the rounding of C*
1374         // if (0 < f* - 1/2 < 10^(-x)) then
1375         //   the result is exact
1376         // else // if (f* - 1/2 > T*) then
1377         //   the result is inexact
1378         if (ind - 1 <= 2) {
1379           if (fstar.w[1] > 0x8000000000000000ull ||
1380               (fstar.w[1] == 0x8000000000000000ull &&
1381               fstar.w[0] > 0x0ull)) { // f* > 1/2 and the result may be exact
1382             tmp64 = fstar.w[1] - 0x8000000000000000ull; // f* - 1/2
1383             if ((tmp64 > __bid_ten2mk128trunc[ind - 1].w[1]
1384                  || (tmp64 == __bid_ten2mk128trunc[ind - 1].w[1]
1385                  && fstar.w[0] >= __bid_ten2mk128trunc[ind - 1].w[0]))) {
1386               // set the inexact flag
1387               *pfpsf |= INEXACT_EXCEPTION;
1388               is_inexact_lt_midpoint = 1;
1389             } // else the result is exact
1390           } else { // the result is inexact; f2* <= 1/2
1391             // set the inexact flag
1392             *pfpsf |= INEXACT_EXCEPTION;
1393             is_inexact_gt_midpoint = 1;
1394           }
1395         } else if (ind - 1 <= 21) { // if 3 <= ind <= 21
1396           if (fstar.w[3] > 0x0 || (fstar.w[3] == 0x0
1397               && fstar.w[2] > __bid_one_half128[ind - 1]) || (fstar.w[3] == 0x0
1398               && fstar.w[2] == __bid_one_half128[ind - 1] && (fstar.w[1]
1399                                                        || fstar.w[0]))) {
1400             // f2* > 1/2 and the result may be exact
1401             // Calculate f2* - 1/2
1402             tmp64 = fstar.w[2] - __bid_one_half128[ind - 1];
1403             tmp64A = fstar.w[3];
1404             if (tmp64 > fstar.w[2])
1405               tmp64A--;
1406             if (tmp64A || tmp64
1407                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
1408                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1409                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
1410               // set the inexact flag
1411               *pfpsf |= INEXACT_EXCEPTION;
1412               is_inexact_lt_midpoint = 1;
1413             } // else the result is exact
1414           } else { // the result is inexact; f2* <= 1/2
1415             // set the inexact flag
1416             *pfpsf |= INEXACT_EXCEPTION;
1417             is_inexact_gt_midpoint = 1;
1418           }
1419         } else { // if 22 <= ind <= 33
1420           if (fstar.w[3] > __bid_one_half128[ind - 1]
1421               || (fstar.w[3] == __bid_one_half128[ind - 1] && (fstar.w[2]
1422               || fstar.w[1] || fstar.w[0]))) {
1423             // f2* > 1/2 and the result may be exact
1424             // Calculate f2* - 1/2
1425             tmp64 = fstar.w[3] - __bid_one_half128[ind - 1];
1426             if (tmp64 || fstar.w[2]
1427                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
1428                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1429                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
1430               // set the inexact flag
1431               *pfpsf |= INEXACT_EXCEPTION;
1432               is_inexact_lt_midpoint = 1;
1433             } // else the result is exact
1434           } else { // the result is inexact; f2* <= 1/2
1435             // set the inexact flag
1436             *pfpsf |= INEXACT_EXCEPTION;
1437             is_inexact_gt_midpoint = 1;
1438           }
1439         }
1440
1441         // if the result was a midpoint it was rounded away from zero, so
1442         // it will need a correction
1443         // check for midpoints
1444         if ((fstar.w[3] == 0) && (fstar.w[2] == 0)
1445             && (fstar.w[1] || fstar.w[0])
1446             && (fstar.w[1] < __bid_ten2mk128trunc[ind - 1].w[1]
1447                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1448                 && fstar.w[0] <= __bid_ten2mk128trunc[ind - 1].w[0]))) {
1449           // the result is a midpoint; round to nearest
1450           if (Cstar.w[0] & 0x01) { // Cstar.w[0] is odd; MP in [EVEN, ODD]
1451             // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
1452             Cstar.w[0]--; // Cstar.w[0] is now even
1453             is_midpoint_gt_even = 1;
1454             is_inexact_lt_midpoint = 0;
1455             is_inexact_gt_midpoint = 0;
1456           } else { // else MP in [ODD, EVEN]
1457             is_midpoint_lt_even = 1;
1458             is_inexact_lt_midpoint = 0;
1459             is_inexact_gt_midpoint = 0;
1460           }
1461         }
1462         // general correction for RM
1463         if (x_sign && (is_midpoint_gt_even || is_inexact_lt_midpoint)) {
1464           Cstar.w[0] = Cstar.w[0] + 1;
1465         } else if (!x_sign
1466                    && (is_midpoint_lt_even || is_inexact_gt_midpoint)) {
1467           Cstar.w[0] = Cstar.w[0] - 1;
1468         } else {
1469           ; // the result is already correct
1470         }
1471         if (x_sign)
1472           res = -Cstar.w[0];
1473         else
1474           res = Cstar.w[0];
1475       } else if (exp == 0) {
1476         // 1 <= q <= 10
1477         // res = +/-C (exact)
1478         if (x_sign)
1479           res = -C1.w[0];
1480         else
1481           res = C1.w[0];
1482       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
1483         // res = +/-C * 10^exp (exact)
1484         if (x_sign)
1485           res = -C1.w[0] * __bid_ten2k64[exp];
1486         else
1487           res = C1.w[0] * __bid_ten2k64[exp];
1488       }
1489     }
1490   }
1491   BID_RETURN (res);
1492 }
1493
1494 /*****************************************************************************
1495  *  BID128_to_int32_ceil
1496  ****************************************************************************/
1497
1498 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_ceil, x)
1499
1500   int res;
1501   UINT64 x_sign;
1502   UINT64 x_exp;
1503   int exp;      // unbiased exponent
1504   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
1505   UINT64 tmp64, tmp64A;
1506   BID_UI64DOUBLE tmp1;
1507   unsigned int x_nr_bits;
1508   int q, ind, shift;
1509   UINT128 C1, C;
1510   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
1511   UINT256 fstar;
1512   UINT256 P256;
1513   int is_inexact_lt_midpoint = 0;
1514   int is_inexact_gt_midpoint = 0;
1515   int is_midpoint_lt_even = 0;
1516   int is_midpoint_gt_even = 0;
1517
1518   // unpack x
1519   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
1520   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
1521   C1.w[1] = x.w[1] & MASK_COEFF;
1522   C1.w[0] = x.w[0];
1523
1524   // check for NaN or Infinity
1525   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
1526     // x is special
1527     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
1528       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
1529         // set invalid flag
1530         *pfpsf |= INVALID_EXCEPTION;
1531         // return Integer Indefinite
1532         res = 0x80000000;
1533       } else { // x is QNaN
1534         // set invalid flag
1535         *pfpsf |= INVALID_EXCEPTION;
1536         // return Integer Indefinite
1537         res = 0x80000000;
1538       }
1539       BID_RETURN (res);
1540     } else { // x is not a NaN, so it must be infinity
1541       if (!x_sign) { // x is +inf
1542         // set invalid flag
1543         *pfpsf |= INVALID_EXCEPTION;
1544         // return Integer Indefinite
1545         res = 0x80000000;
1546       } else { // x is -inf
1547         // set invalid flag
1548         *pfpsf |= INVALID_EXCEPTION;
1549         // return Integer Indefinite
1550         res = 0x80000000;
1551       }
1552       BID_RETURN (res);
1553     }
1554   }
1555   // check for non-canonical values (after the check for special values)
1556   if ((C1.w[1] > 0x0001ed09bead87c0ull)
1557       || (C1.w[1] == 0x0001ed09bead87c0ull
1558           && (C1.w[0] > 0x378d8e63ffffffffull))
1559       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
1560     res = 0x00000000;
1561     BID_RETURN (res);
1562   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
1563     // x is 0
1564     res = 0x00000000;
1565     BID_RETURN (res);
1566   } else { // x is not special and is not zero
1567
1568     // q = nr. of decimal digits in x
1569     //  determine first the nr. of bits in x
1570     if (C1.w[1] == 0) {
1571       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
1572         // split the 64-bit value in two 32-bit halves to avoid rounding errors
1573         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
1574           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
1575           x_nr_bits =
1576             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1577         } else { // x < 2^32
1578           tmp1.d = (double) (C1.w[0]); // exact conversion
1579           x_nr_bits =
1580             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1581         }
1582       } else { // if x < 2^53
1583         tmp1.d = (double) C1.w[0]; // exact conversion
1584         x_nr_bits =
1585           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1586       }
1587     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
1588       tmp1.d = (double) C1.w[1]; // exact conversion
1589       x_nr_bits =
1590         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1591     }
1592     q = __bid_nr_digits[x_nr_bits - 1].digits;
1593     if (q == 0) {
1594       q = __bid_nr_digits[x_nr_bits - 1].digits1;
1595       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
1596           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
1597           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
1598         q++;
1599     }
1600     exp = (x_exp >> 49) - 6176;
1601     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
1602       // set invalid flag
1603       *pfpsf |= INVALID_EXCEPTION;
1604       // return Integer Indefinite
1605       res = 0x80000000;
1606       BID_RETURN (res);
1607     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
1608       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
1609       // so x rounded to an integer may or may not fit in a signed 32-bit int
1610       // the cases that do not fit are identified here; the ones that fit
1611       // fall through and will be handled with other cases further,
1612       // under '1 <= q + exp <= 10'
1613       if (x_sign) { // if n < 0 and q + exp = 10
1614         // if n <= -2^31-1 then n is too large
1615         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1
1616         // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x50000000a, 1<=q<=34
1617         if (q <= 11) {
1618           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
1619           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
1620           if (tmp64 >= 0x50000000aull) {
1621             // set invalid flag
1622             *pfpsf |= INVALID_EXCEPTION;
1623             // return Integer Indefinite
1624             res = 0x80000000;
1625             BID_RETURN (res);
1626           }
1627           // else cases that can be rounded to a 32-bit int fall through
1628           // to '1 <= q + exp <= 10'
1629         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
1630           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x50000000a <=>
1631           // C >= 0x50000000a * 10^(q-11) where 1 <= q - 11 <= 23
1632           // (scale 2^31+1 up)
1633           tmp64 = 0x50000000aull;
1634           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
1635             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
1636           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
1637             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
1638           }
1639           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
1640             // set invalid flag
1641             *pfpsf |= INVALID_EXCEPTION;
1642             // return Integer Indefinite
1643             res = 0x80000000;
1644             BID_RETURN (res);
1645           }
1646           // else cases that can be rounded to a 32-bit int fall through
1647           // to '1 <= q + exp <= 10'
1648         }
1649       } else { // if n > 0 and q + exp = 10
1650         // if n > 2^31 - 1 then n is too large
1651         // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31 - 1
1652         // too large if 0.c(0)c(1)...c(q-1) * 10^11 > 0x4fffffff6, 1<=q<=34
1653         if (q <= 11) {
1654           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
1655           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
1656           if (tmp64 > 0x4fffffff6ull) {
1657             // set invalid flag
1658             *pfpsf |= INVALID_EXCEPTION;
1659             // return Integer Indefinite
1660             res = 0x80000000;
1661             BID_RETURN (res);
1662           }
1663           // else cases that can be rounded to a 32-bit int fall through
1664           // to '1 <= q + exp <= 10'
1665         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
1666           // 0.c(0)c(1)...c(q-1) * 10^11 > 0x4fffffff6 <=>
1667           // C > 0x4fffffff6 * 10^(q-11) where 1 <= q - 11 <= 23
1668           // (scale 2^31 up)
1669           tmp64 = 0x4fffffff6ull;
1670           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
1671             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
1672           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
1673             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
1674           }
1675           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] > C.w[0])) {
1676             // set invalid flag
1677             *pfpsf |= INVALID_EXCEPTION;
1678             // return Integer Indefinite
1679             res = 0x80000000;
1680             BID_RETURN (res);
1681           }
1682           // else cases that can be rounded to a 32-bit int fall through
1683           // to '1 <= q + exp <= 10'
1684         }
1685       }
1686     }
1687     // n is not too large to be converted to int32: -2^31-1 < n <= 2^31-1
1688     // Note: some of the cases tested for above fall through to this point
1689     if ((q + exp) <= 0) {
1690       // n = +/-0.0...c(0)c(1)...c(q-1) or n = +/-0.c(0)c(1)...c(q-1)
1691       // return 0
1692       if (x_sign)
1693         res = 0x00000000;
1694       else
1695         res = 0x00000001;
1696       BID_RETURN (res);
1697     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
1698       // -2^31-1 < x <= -1 or 1 <= x <= 2^31-1 so x can be rounded
1699       // toward positive infinity to a 32-bit signed integer
1700       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
1701         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
1702         // chop off ind digits from the lower part of C1
1703         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
1704         tmp64 = C1.w[0];
1705         if (ind <= 19) {
1706           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
1707         } else {
1708           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
1709           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
1710         }
1711         if (C1.w[0] < tmp64)
1712           C1.w[1]++;
1713         // calculate C* and f*
1714         // C* is actually floor(C*) in this case
1715         // C* and f* need shifting and masking, as shown by
1716         // __bid_shiftright128[] and __bid_maskhigh128[]
1717         // 1 <= x <= 33
1718         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
1719         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
1720         // the approximation of 10^(-x) was rounded up to 118 bits
1721         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
1722         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
1723           Cstar.w[1] = P256.w[3];
1724           Cstar.w[0] = P256.w[2];
1725           fstar.w[3] = 0;
1726           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
1727           fstar.w[1] = P256.w[1];
1728           fstar.w[0] = P256.w[0];
1729         } else { // 22 <= ind - 1 <= 33
1730           Cstar.w[1] = 0;
1731           Cstar.w[0] = P256.w[3];
1732           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
1733           fstar.w[2] = P256.w[2];
1734           fstar.w[1] = P256.w[1];
1735           fstar.w[0] = P256.w[0];
1736         }
1737         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
1738         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
1739         // if (0 < f* < 10^(-x)) then the result is a midpoint
1740         //   if floor(C*) is even then C* = floor(C*) - logical right
1741         //       shift; C* has p decimal digits, correct by Prop. 1)
1742         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
1743         //       shift; C* has p decimal digits, correct by Pr. 1)
1744         // else
1745         //   C* = floor(C*) (logical right shift; C has p decimal digits,
1746         //       correct by Property 1)
1747         // n = C* * 10^(e+x)
1748
1749         // shift right C* by Ex-128 = __bid_shiftright128[ind]
1750         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
1751         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
1752           Cstar.w[0] =
1753             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
1754           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
1755         } else { // 22 <= ind - 1 <= 33
1756           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
1757         }
1758         // determine inexactness of the rounding of C*
1759         // if (0 < f* - 1/2 < 10^(-x)) then
1760         //   the result is exact
1761         // else // if (f* - 1/2 > T*) then
1762         //   the result is inexact
1763         if (ind - 1 <= 2) {
1764           if (fstar.w[1] > 0x8000000000000000ull ||
1765               (fstar.w[1] == 0x8000000000000000ull &&
1766               fstar.w[0] > 0x0ull)) { // f* > 1/2 and the result may be exact
1767             tmp64 = fstar.w[1] - 0x8000000000000000ull; // f* - 1/2
1768             if ((tmp64 > __bid_ten2mk128trunc[ind - 1].w[1]
1769                  || (tmp64 == __bid_ten2mk128trunc[ind - 1].w[1]
1770                  && fstar.w[0] >= __bid_ten2mk128trunc[ind - 1].w[0]))) {
1771               is_inexact_lt_midpoint = 1;
1772             } // else the result is exact
1773           } else { // the result is inexact; f2* <= 1/2
1774             is_inexact_gt_midpoint = 1;
1775           }
1776         } else if (ind - 1 <= 21) { // if 3 <= ind <= 21
1777           if (fstar.w[3] > 0x0 || (fstar.w[3] == 0x0
1778               && fstar.w[2] > __bid_one_half128[ind - 1]) || (fstar.w[3] == 0x0
1779               && fstar.w[2] == __bid_one_half128[ind - 1] && (fstar.w[1]
1780                                                        || fstar.w[0]))) {
1781             // f2* > 1/2 and the result may be exact
1782             // Calculate f2* - 1/2
1783             tmp64 = fstar.w[2] - __bid_one_half128[ind - 1];
1784             tmp64A = fstar.w[3];
1785             if (tmp64 > fstar.w[2])
1786               tmp64A--;
1787             if (tmp64A || tmp64
1788                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
1789                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1790                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
1791               is_inexact_lt_midpoint = 1;
1792             } // else the result is exact
1793           } else { // the result is inexact; f2* <= 1/2
1794             is_inexact_gt_midpoint = 1;
1795           }
1796         } else { // if 22 <= ind <= 33
1797           if (fstar.w[3] > __bid_one_half128[ind - 1]
1798               || (fstar.w[3] == __bid_one_half128[ind - 1] && (fstar.w[2]
1799                                                        || fstar.w[1]
1800                                                        || fstar.w[0]))) {
1801             // f2* > 1/2 and the result may be exact
1802             // Calculate f2* - 1/2
1803             tmp64 = fstar.w[3] - __bid_one_half128[ind - 1];
1804             if (tmp64 || fstar.w[2]
1805                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
1806                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1807                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
1808               is_inexact_lt_midpoint = 1;
1809             } // else the result is exact
1810           } else { // the result is inexact; f2* <= 1/2
1811             is_inexact_gt_midpoint = 1;
1812           }
1813         }
1814
1815         // if the result was a midpoint it was rounded away from zero, so
1816         // it will need a correction
1817         // check for midpoints
1818         if ((fstar.w[3] == 0) && (fstar.w[2] == 0)
1819             && (fstar.w[1] || fstar.w[0])
1820             && (fstar.w[1] < __bid_ten2mk128trunc[ind - 1].w[1]
1821                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
1822                 && fstar.w[0] <= __bid_ten2mk128trunc[ind - 1].w[0]))) {
1823           // the result is a midpoint; round to nearest
1824           if (Cstar.w[0] & 0x01) { // Cstar.w[0] is odd; MP in [EVEN, ODD]
1825             // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
1826             Cstar.w[0]--; // Cstar.w[0] is now even
1827             is_midpoint_gt_even = 1;
1828             is_inexact_lt_midpoint = 0;
1829             is_inexact_gt_midpoint = 0;
1830           } else { // else MP in [ODD, EVEN]
1831             is_midpoint_lt_even = 1;
1832             is_inexact_lt_midpoint = 0;
1833             is_inexact_gt_midpoint = 0;
1834           }
1835         }
1836         // general correction for RM
1837         if (x_sign && (is_midpoint_lt_even || is_inexact_gt_midpoint)) {
1838           Cstar.w[0] = Cstar.w[0] - 1;
1839         } else if (!x_sign
1840                    && (is_midpoint_gt_even || is_inexact_lt_midpoint)) {
1841           Cstar.w[0] = Cstar.w[0] + 1;
1842         } else {
1843           ; // the result is already correct
1844         }
1845         if (x_sign)
1846           res = -Cstar.w[0];
1847         else
1848           res = Cstar.w[0];
1849       } else if (exp == 0) {
1850         // 1 <= q <= 10
1851         // res = +/-C (exact)
1852         if (x_sign)
1853           res = -C1.w[0];
1854         else
1855           res = C1.w[0];
1856       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
1857         // res = +/-C * 10^exp (exact)
1858         if (x_sign)
1859           res = -C1.w[0] * __bid_ten2k64[exp];
1860         else
1861           res = C1.w[0] * __bid_ten2k64[exp];
1862       }
1863     }
1864   }
1865   BID_RETURN (res);
1866 }
1867
1868 /*****************************************************************************
1869  *  BID128_to_int32_xceil
1870  ****************************************************************************/
1871
1872 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_xceil, x)
1873
1874   int res;
1875   UINT64 x_sign;
1876   UINT64 x_exp;
1877   int exp;      // unbiased exponent
1878   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
1879   UINT64 tmp64, tmp64A;
1880   BID_UI64DOUBLE tmp1;
1881   unsigned int x_nr_bits;
1882   int q, ind, shift;
1883   UINT128 C1, C;
1884   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
1885   UINT256 fstar;
1886   UINT256 P256;
1887   int is_inexact_lt_midpoint = 0;
1888   int is_inexact_gt_midpoint = 0;
1889   int is_midpoint_lt_even = 0;
1890   int is_midpoint_gt_even = 0;
1891
1892   // unpack x
1893   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
1894   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
1895   C1.w[1] = x.w[1] & MASK_COEFF;
1896   C1.w[0] = x.w[0];
1897
1898   // check for NaN or Infinity
1899   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
1900     // x is special
1901     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
1902       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
1903         // set invalid flag
1904         *pfpsf |= INVALID_EXCEPTION;
1905         // return Integer Indefinite
1906         res = 0x80000000;
1907       } else { // x is QNaN
1908         // set invalid flag
1909         *pfpsf |= INVALID_EXCEPTION;
1910         // return Integer Indefinite
1911         res = 0x80000000;
1912       }
1913       BID_RETURN (res);
1914     } else { // x is not a NaN, so it must be infinity
1915       if (!x_sign) { // x is +inf
1916         // set invalid flag
1917         *pfpsf |= INVALID_EXCEPTION;
1918         // return Integer Indefinite
1919         res = 0x80000000;
1920       } else { // x is -inf
1921         // set invalid flag
1922         *pfpsf |= INVALID_EXCEPTION;
1923         // return Integer Indefinite
1924         res = 0x80000000;
1925       }
1926       BID_RETURN (res);
1927     }
1928   }
1929   // check for non-canonical values (after the check for special values)
1930   if ((C1.w[1] > 0x0001ed09bead87c0ull)
1931       || (C1.w[1] == 0x0001ed09bead87c0ull
1932           && (C1.w[0] > 0x378d8e63ffffffffull))
1933       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
1934     res = 0x00000000;
1935     BID_RETURN (res);
1936   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
1937     // x is 0
1938     res = 0x00000000;
1939     BID_RETURN (res);
1940   } else { // x is not special and is not zero
1941
1942     // q = nr. of decimal digits in x
1943     //  determine first the nr. of bits in x
1944     if (C1.w[1] == 0) {
1945       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
1946         // split the 64-bit value in two 32-bit halves to avoid rounding errors
1947         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
1948           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
1949           x_nr_bits =
1950             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1951         } else { // x < 2^32
1952           tmp1.d = (double) (C1.w[0]); // exact conversion
1953           x_nr_bits =
1954             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1955         }
1956       } else { // if x < 2^53
1957         tmp1.d = (double) C1.w[0]; // exact conversion
1958         x_nr_bits =
1959           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1960       }
1961     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
1962       tmp1.d = (double) C1.w[1]; // exact conversion
1963       x_nr_bits =
1964         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
1965     }
1966     q = __bid_nr_digits[x_nr_bits - 1].digits;
1967     if (q == 0) {
1968       q = __bid_nr_digits[x_nr_bits - 1].digits1;
1969       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
1970           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
1971           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
1972         q++;
1973     }
1974     exp = (x_exp >> 49) - 6176;
1975     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
1976       // set invalid flag
1977       *pfpsf |= INVALID_EXCEPTION;
1978       // return Integer Indefinite
1979       res = 0x80000000;
1980       BID_RETURN (res);
1981     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
1982       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
1983       // so x rounded to an integer may or may not fit in a signed 32-bit int
1984       // the cases that do not fit are identified here; the ones that fit
1985       // fall through and will be handled with other cases further,
1986       // under '1 <= q + exp <= 10'
1987       if (x_sign) { // if n < 0 and q + exp = 10
1988         // if n <= -2^31-1 then n is too large
1989         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1
1990         // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x50000000a, 1<=q<=34
1991         if (q <= 11) {
1992           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
1993           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
1994           if (tmp64 >= 0x50000000aull) {
1995             // set invalid flag
1996             *pfpsf |= INVALID_EXCEPTION;
1997             // return Integer Indefinite
1998             res = 0x80000000;
1999             BID_RETURN (res);
2000           }
2001           // else cases that can be rounded to a 32-bit int fall through
2002           // to '1 <= q + exp <= 10'
2003         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
2004           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x50000000a <=>
2005           // C >= 0x50000000a * 10^(q-11) where 1 <= q - 11 <= 23
2006           // (scale 2^31+1 up)
2007           tmp64 = 0x50000000aull;
2008           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
2009             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
2010           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
2011             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
2012           }
2013           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
2014             // set invalid flag
2015             *pfpsf |= INVALID_EXCEPTION;
2016             // return Integer Indefinite
2017             res = 0x80000000;
2018             BID_RETURN (res);
2019           }
2020           // else cases that can be rounded to a 32-bit int fall through
2021           // to '1 <= q + exp <= 10'
2022         }
2023       } else { // if n > 0 and q + exp = 10
2024         // if n > 2^31 - 1 then n is too large
2025         // too large if c(0)c(1)...c(9).c(10)...c(q-1) > 2^31 - 1
2026         // too large if 0.c(0)c(1)...c(q-1) * 10^11 > 0x4fffffff6, 1<=q<=34
2027         if (q <= 11) {
2028           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
2029           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
2030           if (tmp64 > 0x4fffffff6ull) {
2031             // set invalid flag
2032             *pfpsf |= INVALID_EXCEPTION;
2033             // return Integer Indefinite
2034             res = 0x80000000;
2035             BID_RETURN (res);
2036           }
2037           // else cases that can be rounded to a 32-bit int fall through
2038           // to '1 <= q + exp <= 10'
2039         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
2040           // 0.c(0)c(1)...c(q-1) * 10^11 > 0x4fffffff6 <=>
2041           // C > 0x4fffffff6 * 10^(q-11) where 1 <= q - 11 <= 23
2042           // (scale 2^31 up)
2043           tmp64 = 0x4fffffff6ull;
2044           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
2045             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
2046           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
2047             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
2048           }
2049           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] > C.w[0])) {
2050             // set invalid flag
2051             *pfpsf |= INVALID_EXCEPTION;
2052             // return Integer Indefinite
2053             res = 0x80000000;
2054             BID_RETURN (res);
2055           }
2056           // else cases that can be rounded to a 32-bit int fall through
2057           // to '1 <= q + exp <= 10'
2058         }
2059       }
2060     }
2061     // n is not too large to be converted to int32: -2^31-1 < n <= 2^31-1
2062     // Note: some of the cases tested for above fall through to this point
2063     if ((q + exp) <= 0) {
2064       // n = +/-0.0...c(0)c(1)...c(q-1) or n = +/-0.c(0)c(1)...c(q-1)
2065       // set inexact flag
2066       *pfpsf |= INEXACT_EXCEPTION;
2067       // return 0
2068       if (x_sign)
2069         res = 0x00000000;
2070       else
2071         res = 0x00000001;
2072       BID_RETURN (res);
2073     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
2074       // -2^31-1 < x <= -1 or 1 <= x <= 2^31-1 so x can be rounded
2075       // toward positive infinity to a 32-bit signed integer
2076       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
2077         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
2078         // chop off ind digits from the lower part of C1
2079         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
2080         tmp64 = C1.w[0];
2081         if (ind <= 19) {
2082           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
2083         } else {
2084           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
2085           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
2086         }
2087         if (C1.w[0] < tmp64)
2088           C1.w[1]++;
2089         // calculate C* and f*
2090         // C* is actually floor(C*) in this case
2091         // C* and f* need shifting and masking, as shown by
2092         // __bid_shiftright128[] and __bid_maskhigh128[]
2093         // 1 <= x <= 33
2094         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
2095         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
2096         // the approximation of 10^(-x) was rounded up to 118 bits
2097         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
2098         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
2099           Cstar.w[1] = P256.w[3];
2100           Cstar.w[0] = P256.w[2];
2101           fstar.w[3] = 0;
2102           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
2103           fstar.w[1] = P256.w[1];
2104           fstar.w[0] = P256.w[0];
2105         } else { // 22 <= ind - 1 <= 33
2106           Cstar.w[1] = 0;
2107           Cstar.w[0] = P256.w[3];
2108           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
2109           fstar.w[2] = P256.w[2];
2110           fstar.w[1] = P256.w[1];
2111           fstar.w[0] = P256.w[0];
2112         }
2113         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
2114         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
2115         // if (0 < f* < 10^(-x)) then the result is a midpoint
2116         //   if floor(C*) is even then C* = floor(C*) - logical right
2117         //       shift; C* has p decimal digits, correct by Prop. 1)
2118         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
2119         //       shift; C* has p decimal digits, correct by Pr. 1)
2120         // else
2121         //   C* = floor(C*) (logical right shift; C has p decimal digits,
2122         //       correct by Property 1)
2123         // n = C* * 10^(e+x)
2124
2125         // shift right C* by Ex-128 = __bid_shiftright128[ind]
2126         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
2127         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
2128           Cstar.w[0] =
2129             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
2130           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
2131         } else { // 22 <= ind - 1 <= 33
2132           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
2133         }
2134         // determine inexactness of the rounding of C*
2135         // if (0 < f* - 1/2 < 10^(-x)) then
2136         //   the result is exact
2137         // else // if (f* - 1/2 > T*) then
2138         //   the result is inexact
2139         if (ind - 1 <= 2) {
2140           if (fstar.w[1] > 0x8000000000000000ull ||
2141               (fstar.w[1] == 0x8000000000000000ull &&
2142               fstar.w[0] > 0x0ull)) { // f* > 1/2 and the result may be exact
2143             tmp64 = fstar.w[1] - 0x8000000000000000ull; // f* - 1/2
2144             if ((tmp64 > __bid_ten2mk128trunc[ind - 1].w[1]
2145                  || (tmp64 == __bid_ten2mk128trunc[ind - 1].w[1]
2146                  && fstar.w[0] >= __bid_ten2mk128trunc[ind - 1].w[0]))) {
2147               // set the inexact flag
2148               *pfpsf |= INEXACT_EXCEPTION;
2149               is_inexact_lt_midpoint = 1;
2150             } // else the result is exact
2151           } else { // the result is inexact; f2* <= 1/2
2152             // set the inexact flag
2153             *pfpsf |= INEXACT_EXCEPTION;
2154             is_inexact_gt_midpoint = 1;
2155           }
2156         } else if (ind - 1 <= 21) { // if 3 <= ind <= 21
2157           if (fstar.w[3] > 0x0 || (fstar.w[3] == 0x0
2158               && fstar.w[2] > __bid_one_half128[ind - 1]) || (fstar.w[3] == 0x0
2159               && fstar.w[2] == __bid_one_half128[ind - 1] && (fstar.w[1]
2160                                                        || fstar.w[0]))) {
2161             // f2* > 1/2 and the result may be exact
2162             // Calculate f2* - 1/2
2163             tmp64 = fstar.w[2] - __bid_one_half128[ind - 1];
2164             tmp64A = fstar.w[3];
2165             if (tmp64 > fstar.w[2])
2166               tmp64A--;
2167             if (tmp64A || tmp64
2168                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
2169                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2170                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
2171               // set the inexact flag
2172               *pfpsf |= INEXACT_EXCEPTION;
2173               is_inexact_lt_midpoint = 1;
2174             } // else the result is exact
2175           } else { // the result is inexact; f2* <= 1/2
2176             // set the inexact flag
2177             *pfpsf |= INEXACT_EXCEPTION;
2178             is_inexact_gt_midpoint = 1;
2179           }
2180         } else { // if 22 <= ind <= 33
2181           if (fstar.w[3] > __bid_one_half128[ind - 1]
2182               || (fstar.w[3] == __bid_one_half128[ind - 1] && (fstar.w[2]
2183                                                        || fstar.w[1]
2184                                                        || fstar.w[0]))) {
2185             // f2* > 1/2 and the result may be exact
2186             // Calculate f2* - 1/2
2187             tmp64 = fstar.w[3] - __bid_one_half128[ind - 1];
2188             if (tmp64 || fstar.w[2]
2189                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
2190                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2191                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
2192               // set the inexact flag
2193               *pfpsf |= INEXACT_EXCEPTION;
2194               is_inexact_lt_midpoint = 1;
2195             } // else the result is exact
2196           } else { // the result is inexact; f2* <= 1/2
2197             // set the inexact flag
2198             *pfpsf |= INEXACT_EXCEPTION;
2199             is_inexact_gt_midpoint = 1;
2200           }
2201         }
2202
2203         // if the result was a midpoint it was rounded away from zero, so
2204         // it will need a correction
2205         // check for midpoints
2206         if ((fstar.w[3] == 0) && (fstar.w[2] == 0)
2207             && (fstar.w[1] || fstar.w[0])
2208             && (fstar.w[1] < __bid_ten2mk128trunc[ind - 1].w[1]
2209                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2210                 && fstar.w[0] <= __bid_ten2mk128trunc[ind - 1].w[0]))) {
2211           // the result is a midpoint; round to nearest
2212           if (Cstar.w[0] & 0x01) { // Cstar.w[0] is odd; MP in [EVEN, ODD]
2213             // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
2214             Cstar.w[0]--; // Cstar.w[0] is now even
2215             is_midpoint_gt_even = 1;
2216             is_inexact_lt_midpoint = 0;
2217             is_inexact_gt_midpoint = 0;
2218           } else { // else MP in [ODD, EVEN]
2219             is_midpoint_lt_even = 1;
2220             is_inexact_lt_midpoint = 0;
2221             is_inexact_gt_midpoint = 0;
2222           }
2223         }
2224         // general correction for RM
2225         if (x_sign && (is_midpoint_lt_even || is_inexact_gt_midpoint)) {
2226           Cstar.w[0] = Cstar.w[0] - 1;
2227         } else if (!x_sign
2228                    && (is_midpoint_gt_even || is_inexact_lt_midpoint)) {
2229           Cstar.w[0] = Cstar.w[0] + 1;
2230         } else {
2231           ; // the result is already correct
2232         }
2233         if (x_sign)
2234           res = -Cstar.w[0];
2235         else
2236           res = Cstar.w[0];
2237       } else if (exp == 0) {
2238         // 1 <= q <= 10
2239         // res = +/-C (exact)
2240         if (x_sign)
2241           res = -C1.w[0];
2242         else
2243           res = C1.w[0];
2244       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
2245         // res = +/-C * 10^exp (exact)
2246         if (x_sign)
2247           res = -C1.w[0] * __bid_ten2k64[exp];
2248         else
2249           res = C1.w[0] * __bid_ten2k64[exp];
2250       }
2251     }
2252   }
2253   BID_RETURN (res);
2254 }
2255
2256 /*****************************************************************************
2257  *  BID128_to_int32_int
2258  ****************************************************************************/
2259
2260 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_int, x)
2261
2262   int res;
2263   UINT64 x_sign;
2264   UINT64 x_exp;
2265   int exp;      // unbiased exponent
2266   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
2267   UINT64 tmp64, tmp64A;
2268   BID_UI64DOUBLE tmp1;
2269   unsigned int x_nr_bits;
2270   int q, ind, shift;
2271   UINT128 C1, C;
2272   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
2273   UINT256 fstar;
2274   UINT256 P256;
2275   int is_inexact_gt_midpoint = 0;
2276   int is_midpoint_lt_even = 0;
2277
2278   // unpack x
2279   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
2280   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
2281   C1.w[1] = x.w[1] & MASK_COEFF;
2282   C1.w[0] = x.w[0];
2283
2284   // check for NaN or Infinity
2285   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
2286     // x is special
2287     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
2288       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
2289         // set invalid flag
2290         *pfpsf |= INVALID_EXCEPTION;
2291         // return Integer Indefinite
2292         res = 0x80000000;
2293       } else { // x is QNaN
2294         // set invalid flag
2295         *pfpsf |= INVALID_EXCEPTION;
2296         // return Integer Indefinite
2297         res = 0x80000000;
2298       }
2299       BID_RETURN (res);
2300     } else { // x is not a NaN, so it must be infinity
2301       if (!x_sign) { // x is +inf
2302         // set invalid flag
2303         *pfpsf |= INVALID_EXCEPTION;
2304         // return Integer Indefinite
2305         res = 0x80000000;
2306       } else { // x is -inf
2307         // set invalid flag
2308         *pfpsf |= INVALID_EXCEPTION;
2309         // return Integer Indefinite
2310         res = 0x80000000;
2311       }
2312       BID_RETURN (res);
2313     }
2314   }
2315   // check for non-canonical values (after the check for special values)
2316   if ((C1.w[1] > 0x0001ed09bead87c0ull)
2317       || (C1.w[1] == 0x0001ed09bead87c0ull
2318           && (C1.w[0] > 0x378d8e63ffffffffull))
2319       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
2320     res = 0x00000000;
2321     BID_RETURN (res);
2322   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
2323     // x is 0
2324     res = 0x00000000;
2325     BID_RETURN (res);
2326   } else { // x is not special and is not zero
2327
2328     // q = nr. of decimal digits in x
2329     //  determine first the nr. of bits in x
2330     if (C1.w[1] == 0) {
2331       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
2332         // split the 64-bit value in two 32-bit halves to avoid rounding errors
2333         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
2334           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
2335           x_nr_bits =
2336             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
2337         } else { // x < 2^32
2338           tmp1.d = (double) (C1.w[0]); // exact conversion
2339           x_nr_bits =
2340             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
2341         }
2342       } else { // if x < 2^53
2343         tmp1.d = (double) C1.w[0]; // exact conversion
2344         x_nr_bits =
2345           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
2346       }
2347     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
2348       tmp1.d = (double) C1.w[1]; // exact conversion
2349       x_nr_bits =
2350         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
2351     }
2352     q = __bid_nr_digits[x_nr_bits - 1].digits;
2353     if (q == 0) {
2354       q = __bid_nr_digits[x_nr_bits - 1].digits1;
2355       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
2356           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
2357           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
2358         q++;
2359     }
2360     exp = (x_exp >> 49) - 6176;
2361     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
2362       // set invalid flag
2363       *pfpsf |= INVALID_EXCEPTION;
2364       // return Integer Indefinite
2365       res = 0x80000000;
2366       BID_RETURN (res);
2367     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
2368       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
2369       // so x rounded to an integer may or may not fit in a signed 32-bit int
2370       // the cases that do not fit are identified here; the ones that fit
2371       // fall through and will be handled with other cases further,
2372       // under '1 <= q + exp <= 10'
2373       if (x_sign) { // if n < 0 and q + exp = 10
2374         // if n <= -2^31 - 1 then n is too large
2375         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1
2376         // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x50000000a, 1<=q<=34
2377         if (q <= 11) {
2378           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
2379           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
2380           if (tmp64 >= 0x50000000aull) {
2381             // set invalid flag
2382             *pfpsf |= INVALID_EXCEPTION;
2383             // return Integer Indefinite
2384             res = 0x80000000;
2385             BID_RETURN (res);
2386           }
2387           // else cases that can be rounded to a 32-bit int fall through
2388           // to '1 <= q + exp <= 10'
2389         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
2390           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x50000000a <=>
2391           // C >= 0x50000000a * 10^(q-11) where 1 <= q - 11 <= 23
2392           // (scale 2^31+1 up)
2393           tmp64 = 0x50000000aull;
2394           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
2395             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
2396           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
2397             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
2398           }
2399           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
2400             // set invalid flag
2401             *pfpsf |= INVALID_EXCEPTION;
2402             // return Integer Indefinite
2403             res = 0x80000000;
2404             BID_RETURN (res);
2405           }
2406           // else cases that can be rounded to a 32-bit int fall through
2407           // to '1 <= q + exp <= 10'
2408         }
2409       } else { // if n > 0 and q + exp = 10
2410         // if n >= 2^31 then n is too large
2411         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31
2412         // too large if 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000, 1<=q<=34
2413         if (q <= 11) {
2414           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
2415           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
2416           if (tmp64 >= 0x500000000ull) {
2417             // set invalid flag
2418             *pfpsf |= INVALID_EXCEPTION;
2419             // return Integer Indefinite
2420             res = 0x80000000;
2421             BID_RETURN (res);
2422           }
2423           // else cases that can be rounded to a 32-bit int fall through
2424           // to '1 <= q + exp <= 10'
2425         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
2426           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000 <=>
2427           // C >= 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 23
2428           // (scale 2^31-1/2 up)
2429           tmp64 = 0x500000000ull;
2430           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
2431             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
2432           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
2433             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
2434           }
2435           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
2436             // set invalid flag
2437             *pfpsf |= INVALID_EXCEPTION;
2438             // return Integer Indefinite
2439             res = 0x80000000;
2440             BID_RETURN (res);
2441           }
2442           // else cases that can be rounded to a 32-bit int fall through
2443           // to '1 <= q + exp <= 10'
2444         }
2445       }
2446     }
2447     // n is not too large to be converted to int32: -2^31 - 1 < n < 2^31
2448     // Note: some of the cases tested for above fall through to this point
2449     if ((q + exp) <= 0) {
2450       // n = +/-0.0...c(0)c(1)...c(q-1) or n = +/-0.c(0)c(1)...c(q-1)
2451       // return 0
2452       res = 0x00000000;
2453       BID_RETURN (res);
2454     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
2455       // -2^31-1 < x <= -1 or 1 <= x < 2^31 so x can be rounded
2456       // toward zero to a 32-bit signed integer
2457       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
2458         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
2459         // chop off ind digits from the lower part of C1
2460         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
2461         tmp64 = C1.w[0];
2462         if (ind <= 19) {
2463           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
2464         } else {
2465           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
2466           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
2467         }
2468         if (C1.w[0] < tmp64)
2469           C1.w[1]++;
2470         // calculate C* and f*
2471         // C* is actually floor(C*) in this case
2472         // C* and f* need shifting and masking, as shown by
2473         // __bid_shiftright128[] and __bid_maskhigh128[]
2474         // 1 <= x <= 33
2475         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
2476         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
2477         // the approximation of 10^(-x) was rounded up to 118 bits
2478         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
2479         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
2480           Cstar.w[1] = P256.w[3];
2481           Cstar.w[0] = P256.w[2];
2482           fstar.w[3] = 0;
2483           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
2484           fstar.w[1] = P256.w[1];
2485           fstar.w[0] = P256.w[0];
2486         } else { // 22 <= ind - 1 <= 33
2487           Cstar.w[1] = 0;
2488           Cstar.w[0] = P256.w[3];
2489           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
2490           fstar.w[2] = P256.w[2];
2491           fstar.w[1] = P256.w[1];
2492           fstar.w[0] = P256.w[0];
2493         }
2494         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
2495         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
2496         // if (0 < f* < 10^(-x)) then the result is a midpoint
2497         //   if floor(C*) is even then C* = floor(C*) - logical right
2498         //       shift; C* has p decimal digits, correct by Prop. 1)
2499         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
2500         //       shift; C* has p decimal digits, correct by Pr. 1)
2501         // else
2502         //   C* = floor(C*) (logical right shift; C has p decimal digits,
2503         //       correct by Property 1)
2504         // n = C* * 10^(e+x)
2505
2506         // shift right C* by Ex-128 = __bid_shiftright128[ind]
2507         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
2508         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
2509           Cstar.w[0] =
2510             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
2511           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
2512         } else { // 22 <= ind - 1 <= 33
2513           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
2514         }
2515         // determine inexactness of the rounding of C*
2516         // if (0 < f* - 1/2 < 10^(-x)) then
2517         //   the result is exact
2518         // else // if (f* - 1/2 > T*) then
2519         //   the result is inexact
2520         if (ind - 1 <= 2) {
2521           if (fstar.w[1] > 0x8000000000000000ull ||
2522               (fstar.w[1] == 0x8000000000000000ull &&
2523               fstar.w[0] > 0x0ull)) { // f* > 1/2 and the result may be exact
2524             tmp64 = fstar.w[1] - 0x8000000000000000ull; // f* - 1/2
2525             if ((tmp64 > __bid_ten2mk128trunc[ind - 1].w[1]
2526                  || (tmp64 == __bid_ten2mk128trunc[ind - 1].w[1]
2527                  && fstar.w[0] >= __bid_ten2mk128trunc[ind - 1].w[0]))) {
2528             } // else the result is exact
2529           } else { // the result is inexact; f2* <= 1/2
2530             is_inexact_gt_midpoint = 1;
2531           }
2532         } else if (ind - 1 <= 21) { // if 3 <= ind <= 21
2533           if (fstar.w[3] > 0x0 || (fstar.w[3] == 0x0
2534               && fstar.w[2] > __bid_one_half128[ind - 1]) || (fstar.w[3] == 0x0
2535               && fstar.w[2] == __bid_one_half128[ind - 1] && (fstar.w[1]
2536                                                        || fstar.w[0]))) {
2537             // f2* > 1/2 and the result may be exact
2538             // Calculate f2* - 1/2
2539             tmp64 = fstar.w[2] - __bid_one_half128[ind - 1];
2540             tmp64A = fstar.w[3];
2541             if (tmp64 > fstar.w[2])
2542               tmp64A--;
2543             if (tmp64A || tmp64
2544                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
2545                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2546                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
2547             } // else the result is exact
2548           } else { // the result is inexact; f2* <= 1/2
2549             is_inexact_gt_midpoint = 1;
2550           }
2551         } else { // if 22 <= ind <= 33
2552           if (fstar.w[3] > __bid_one_half128[ind - 1]
2553               || (fstar.w[3] == __bid_one_half128[ind - 1] && (fstar.w[2]
2554                                                        || fstar.w[1]
2555                                                        || fstar.w[0]))) {
2556             // f2* > 1/2 and the result may be exact
2557             // Calculate f2* - 1/2
2558             tmp64 = fstar.w[3] - __bid_one_half128[ind - 1];
2559             if (tmp64 || fstar.w[2]
2560                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
2561                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2562                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
2563             } // else the result is exact
2564           } else { // the result is inexact; f2* <= 1/2
2565             is_inexact_gt_midpoint = 1;
2566           }
2567         }
2568
2569         // if the result was a midpoint it was rounded away from zero, so
2570         // it will need a correction
2571         // check for midpoints
2572         if ((fstar.w[3] == 0) && (fstar.w[2] == 0)
2573             && (fstar.w[1] || fstar.w[0])
2574             && (fstar.w[1] < __bid_ten2mk128trunc[ind - 1].w[1]
2575                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2576                 && fstar.w[0] <= __bid_ten2mk128trunc[ind - 1].w[0]))) {
2577           // the result is a midpoint; round to nearest
2578           if (Cstar.w[0] & 0x01) { // Cstar.w[0] is odd; MP in [EVEN, ODD]
2579             // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
2580             Cstar.w[0]--; // Cstar.w[0] is now even
2581             is_inexact_gt_midpoint = 0;
2582           } else { // else MP in [ODD, EVEN]
2583             is_midpoint_lt_even = 1;
2584             is_inexact_gt_midpoint = 0;
2585           }
2586         }
2587         // general correction for RZ
2588         if (is_midpoint_lt_even || is_inexact_gt_midpoint) {
2589           Cstar.w[0] = Cstar.w[0] - 1;
2590         } else {
2591           ; // exact, the result is already correct
2592         }
2593         if (x_sign)
2594           res = -Cstar.w[0];
2595         else
2596           res = Cstar.w[0];
2597       } else if (exp == 0) {
2598         // 1 <= q <= 10
2599         // res = +/-C (exact)
2600         if (x_sign)
2601           res = -C1.w[0];
2602         else
2603           res = C1.w[0];
2604       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
2605         // res = +/-C * 10^exp (exact)
2606         if (x_sign)
2607           res = -C1.w[0] * __bid_ten2k64[exp];
2608         else
2609           res = C1.w[0] * __bid_ten2k64[exp];
2610       }
2611     }
2612   }
2613   BID_RETURN (res);
2614 }
2615
2616 /*****************************************************************************
2617  *  BID128_to_int32_xint
2618  ****************************************************************************/
2619
2620 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_xint, x)
2621
2622   int res;
2623   UINT64 x_sign;
2624   UINT64 x_exp;
2625   int exp;      // unbiased exponent
2626   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
2627   UINT64 tmp64, tmp64A;
2628   BID_UI64DOUBLE tmp1;
2629   unsigned int x_nr_bits;
2630   int q, ind, shift;
2631   UINT128 C1, C;
2632   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
2633   UINT256 fstar;
2634   UINT256 P256;
2635   int is_inexact_gt_midpoint = 0;
2636   int is_midpoint_lt_even = 0;
2637
2638   // unpack x
2639   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
2640   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
2641   C1.w[1] = x.w[1] & MASK_COEFF;
2642   C1.w[0] = x.w[0];
2643
2644   // check for NaN or Infinity
2645   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
2646     // x is special
2647     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
2648       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
2649         // set invalid flag
2650         *pfpsf |= INVALID_EXCEPTION;
2651         // return Integer Indefinite
2652         res = 0x80000000;
2653       } else { // x is QNaN
2654         // set invalid flag
2655         *pfpsf |= INVALID_EXCEPTION;
2656         // return Integer Indefinite
2657         res = 0x80000000;
2658       }
2659       BID_RETURN (res);
2660     } else { // x is not a NaN, so it must be infinity
2661       if (!x_sign) { // x is +inf
2662         // set invalid flag
2663         *pfpsf |= INVALID_EXCEPTION;
2664         // return Integer Indefinite
2665         res = 0x80000000;
2666       } else { // x is -inf
2667         // set invalid flag
2668         *pfpsf |= INVALID_EXCEPTION;
2669         // return Integer Indefinite
2670         res = 0x80000000;
2671       }
2672       BID_RETURN (res);
2673     }
2674   }
2675   // check for non-canonical values (after the check for special values)
2676   if ((C1.w[1] > 0x0001ed09bead87c0ull)
2677       || (C1.w[1] == 0x0001ed09bead87c0ull
2678           && (C1.w[0] > 0x378d8e63ffffffffull))
2679       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
2680     res = 0x00000000;
2681     BID_RETURN (res);
2682   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
2683     // x is 0
2684     res = 0x00000000;
2685     BID_RETURN (res);
2686   } else { // x is not special and is not zero
2687
2688     // q = nr. of decimal digits in x
2689     //  determine first the nr. of bits in x
2690     if (C1.w[1] == 0) {
2691       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
2692         // split the 64-bit value in two 32-bit halves to avoid rounding errors
2693         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
2694           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
2695           x_nr_bits =
2696             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
2697         } else { // x < 2^32
2698           tmp1.d = (double) (C1.w[0]); // exact conversion
2699           x_nr_bits =
2700             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
2701         }
2702       } else { // if x < 2^53
2703         tmp1.d = (double) C1.w[0]; // exact conversion
2704         x_nr_bits =
2705           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
2706       }
2707     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
2708       tmp1.d = (double) C1.w[1]; // exact conversion
2709       x_nr_bits =
2710         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
2711     }
2712     q = __bid_nr_digits[x_nr_bits - 1].digits;
2713     if (q == 0) {
2714       q = __bid_nr_digits[x_nr_bits - 1].digits1;
2715       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
2716           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
2717           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
2718         q++;
2719     }
2720     exp = (x_exp >> 49) - 6176;
2721     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
2722       // set invalid flag
2723       *pfpsf |= INVALID_EXCEPTION;
2724       // return Integer Indefinite
2725       res = 0x80000000;
2726       BID_RETURN (res);
2727     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
2728       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
2729       // so x rounded to an integer may or may not fit in a signed 32-bit int
2730       // the cases that do not fit are identified here; the ones that fit
2731       // fall through and will be handled with other cases further,
2732       // under '1 <= q + exp <= 10'
2733       if (x_sign) { // if n < 0 and q + exp = 10
2734         // if n <= -2^31 - 1 then n is too large
2735         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1
2736         // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x50000000a, 1<=q<=34
2737         if (q <= 11) {
2738           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
2739           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
2740           if (tmp64 >= 0x50000000aull) {
2741             // set invalid flag
2742             *pfpsf |= INVALID_EXCEPTION;
2743             // return Integer Indefinite
2744             res = 0x80000000;
2745             BID_RETURN (res);
2746           }
2747           // else cases that can be rounded to a 32-bit int fall through
2748           // to '1 <= q + exp <= 10'
2749         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
2750           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x50000000a <=>
2751           // C >= 0x50000000a * 10^(q-11) where 1 <= q - 11 <= 23
2752           // (scale 2^31+1 up)
2753           tmp64 = 0x50000000aull;
2754           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
2755             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
2756           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
2757             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
2758           }
2759           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
2760             // set invalid flag
2761             *pfpsf |= INVALID_EXCEPTION;
2762             // return Integer Indefinite
2763             res = 0x80000000;
2764             BID_RETURN (res);
2765           }
2766           // else cases that can be rounded to a 32-bit int fall through
2767           // to '1 <= q + exp <= 10'
2768         }
2769       } else { // if n > 0 and q + exp = 10
2770         // if n >= 2^31 then n is too large
2771         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31
2772         // too large if 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000, 1<=q<=34
2773         if (q <= 11) {
2774           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
2775           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
2776           if (tmp64 >= 0x500000000ull) {
2777             // set invalid flag
2778             *pfpsf |= INVALID_EXCEPTION;
2779             // return Integer Indefinite
2780             res = 0x80000000;
2781             BID_RETURN (res);
2782           }
2783           // else cases that can be rounded to a 32-bit int fall through
2784           // to '1 <= q + exp <= 10'
2785         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
2786           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000000 <=>
2787           // C >= 0x500000000 * 10^(q-11) where 1 <= q - 11 <= 23
2788           // (scale 2^31-1/2 up)
2789           tmp64 = 0x500000000ull;
2790           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
2791             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
2792           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
2793             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
2794           }
2795           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
2796             // set invalid flag
2797             *pfpsf |= INVALID_EXCEPTION;
2798             // return Integer Indefinite
2799             res = 0x80000000;
2800             BID_RETURN (res);
2801           }
2802           // else cases that can be rounded to a 32-bit int fall through
2803           // to '1 <= q + exp <= 10'
2804         }
2805       }
2806     }
2807     // n is not too large to be converted to int32: -2^31 - 1 < n < 2^31
2808     // Note: some of the cases tested for above fall through to this point
2809     if ((q + exp) <= 0) {
2810       // n = +/-0.0...c(0)c(1)...c(q-1) or n = +/-0.c(0)c(1)...c(q-1)
2811       // set inexact flag
2812       *pfpsf |= INEXACT_EXCEPTION;
2813       // return 0
2814       res = 0x00000000;
2815       BID_RETURN (res);
2816     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
2817       // -2^31-1 < x <= -1 or 1 <= x < 2^31 so x can be rounded
2818       // toward zero to a 32-bit signed integer
2819       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
2820         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
2821         // chop off ind digits from the lower part of C1
2822         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
2823         tmp64 = C1.w[0];
2824         if (ind <= 19) {
2825           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
2826         } else {
2827           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
2828           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
2829         }
2830         if (C1.w[0] < tmp64)
2831           C1.w[1]++;
2832         // calculate C* and f*
2833         // C* is actually floor(C*) in this case
2834         // C* and f* need shifting and masking, as shown by
2835         // __bid_shiftright128[] and __bid_maskhigh128[]
2836         // 1 <= x <= 33
2837         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
2838         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
2839         // the approximation of 10^(-x) was rounded up to 118 bits
2840         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
2841         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
2842           Cstar.w[1] = P256.w[3];
2843           Cstar.w[0] = P256.w[2];
2844           fstar.w[3] = 0;
2845           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
2846           fstar.w[1] = P256.w[1];
2847           fstar.w[0] = P256.w[0];
2848         } else { // 22 <= ind - 1 <= 33
2849           Cstar.w[1] = 0;
2850           Cstar.w[0] = P256.w[3];
2851           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
2852           fstar.w[2] = P256.w[2];
2853           fstar.w[1] = P256.w[1];
2854           fstar.w[0] = P256.w[0];
2855         }
2856         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
2857         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
2858         // if (0 < f* < 10^(-x)) then the result is a midpoint
2859         //   if floor(C*) is even then C* = floor(C*) - logical right
2860         //       shift; C* has p decimal digits, correct by Prop. 1)
2861         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
2862         //       shift; C* has p decimal digits, correct by Pr. 1)
2863         // else
2864         //   C* = floor(C*) (logical right shift; C has p decimal digits,
2865         //       correct by Property 1)
2866         // n = C* * 10^(e+x)
2867
2868         // shift right C* by Ex-128 = __bid_shiftright128[ind]
2869         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
2870         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
2871           Cstar.w[0] =
2872             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
2873           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
2874         } else { // 22 <= ind - 1 <= 33
2875           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
2876         }
2877         // determine inexactness of the rounding of C*
2878         // if (0 < f* - 1/2 < 10^(-x)) then
2879         //   the result is exact
2880         // else // if (f* - 1/2 > T*) then
2881         //   the result is inexact
2882         if (ind - 1 <= 2) {
2883           if (fstar.w[1] > 0x8000000000000000ull ||
2884               (fstar.w[1] == 0x8000000000000000ull &&
2885               fstar.w[0] > 0x0ull)) { // f* > 1/2 and the result may be exact
2886             tmp64 = fstar.w[1] - 0x8000000000000000ull; // f* - 1/2
2887             if ((tmp64 > __bid_ten2mk128trunc[ind - 1].w[1]
2888                  || (tmp64 == __bid_ten2mk128trunc[ind - 1].w[1]
2889                  && fstar.w[0] >= __bid_ten2mk128trunc[ind - 1].w[0]))) {
2890               // set the inexact flag
2891               *pfpsf |= INEXACT_EXCEPTION;
2892             } // else the result is exact
2893           } else { // the result is inexact; f2* <= 1/2
2894             // set the inexact flag
2895             *pfpsf |= INEXACT_EXCEPTION;
2896             is_inexact_gt_midpoint = 1;
2897           }
2898         } else if (ind - 1 <= 21) { // if 3 <= ind <= 21
2899           if (fstar.w[3] > 0x0 || (fstar.w[3] == 0x0
2900               && fstar.w[2] > __bid_one_half128[ind - 1]) || (fstar.w[3] == 0x0
2901               && fstar.w[2] == __bid_one_half128[ind - 1] && (fstar.w[1]
2902                                                        || fstar.w[0]))) {
2903             // f2* > 1/2 and the result may be exact
2904             // Calculate f2* - 1/2
2905             tmp64 = fstar.w[2] - __bid_one_half128[ind - 1];
2906             tmp64A = fstar.w[3];
2907             if (tmp64 > fstar.w[2])
2908               tmp64A--;
2909             if (tmp64A || tmp64
2910                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
2911                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2912                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
2913               // set the inexact flag
2914               *pfpsf |= INEXACT_EXCEPTION;
2915             } // else the result is exact
2916           } else { // the result is inexact; f2* <= 1/2
2917             // set the inexact flag
2918             *pfpsf |= INEXACT_EXCEPTION;
2919             is_inexact_gt_midpoint = 1;
2920           }
2921         } else { // if 22 <= ind <= 33
2922           if (fstar.w[3] > __bid_one_half128[ind - 1]
2923               || (fstar.w[3] == __bid_one_half128[ind - 1] && (fstar.w[2]
2924                                                        || fstar.w[1]
2925                                                        || fstar.w[0]))) {
2926             // f2* > 1/2 and the result may be exact
2927             // Calculate f2* - 1/2
2928             tmp64 = fstar.w[3] - __bid_one_half128[ind - 1];
2929             if (tmp64 || fstar.w[2]
2930                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
2931                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2932                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
2933               // set the inexact flag
2934               *pfpsf |= INEXACT_EXCEPTION;
2935             } // else the result is exact
2936           } else { // the result is inexact; f2* <= 1/2
2937             // set the inexact flag
2938             *pfpsf |= INEXACT_EXCEPTION;
2939             is_inexact_gt_midpoint = 1;
2940           }
2941         }
2942
2943         // if the result was a midpoint it was rounded away from zero, so
2944         // it will need a correction
2945         // check for midpoints
2946         if ((fstar.w[3] == 0) && (fstar.w[2] == 0)
2947             && (fstar.w[1] || fstar.w[0])
2948             && (fstar.w[1] < __bid_ten2mk128trunc[ind - 1].w[1]
2949                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
2950                 && fstar.w[0] <= __bid_ten2mk128trunc[ind - 1].w[0]))) {
2951           // the result is a midpoint; round to nearest
2952           if (Cstar.w[0] & 0x01) { // Cstar.w[0] is odd; MP in [EVEN, ODD]
2953             // if floor(C*) is odd C = floor(C*) - 1; the result >= 1
2954             Cstar.w[0]--; // Cstar.w[0] is now even
2955             is_inexact_gt_midpoint = 0;
2956           } else { // else MP in [ODD, EVEN]
2957             is_midpoint_lt_even = 1;
2958             is_inexact_gt_midpoint = 0;
2959           }
2960         }
2961         // general correction for RZ
2962         if (is_midpoint_lt_even || is_inexact_gt_midpoint) {
2963           Cstar.w[0] = Cstar.w[0] - 1;
2964         } else {
2965           ; // exact, the result is already correct
2966         }
2967         if (x_sign)
2968           res = -Cstar.w[0];
2969         else
2970           res = Cstar.w[0];
2971       } else if (exp == 0) {
2972         // 1 <= q <= 10
2973         // res = +/-C (exact)
2974         if (x_sign)
2975           res = -C1.w[0];
2976         else
2977           res = C1.w[0];
2978       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
2979         // res = +/-C * 10^exp (exact)
2980         if (x_sign)
2981           res = -C1.w[0] * __bid_ten2k64[exp];
2982         else
2983           res = C1.w[0] * __bid_ten2k64[exp];
2984       }
2985     }
2986   }
2987   BID_RETURN (res);
2988 }
2989
2990 /*****************************************************************************
2991  *  BID128_to_int32_rninta
2992  ****************************************************************************/
2993
2994 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_rninta, x)
2995
2996   int res;
2997   UINT64 x_sign;
2998   UINT64 x_exp;
2999   int exp;      // unbiased exponent
3000   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
3001   UINT64 tmp64;
3002   BID_UI64DOUBLE tmp1;
3003   unsigned int x_nr_bits;
3004   int q, ind, shift;
3005   UINT128 C1, C;
3006   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
3007   UINT256 P256;
3008
3009   // unpack x
3010   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
3011   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
3012   C1.w[1] = x.w[1] & MASK_COEFF;
3013   C1.w[0] = x.w[0];
3014
3015   // check for NaN or Infinity
3016   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
3017     // x is special
3018     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
3019       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
3020         // set invalid flag
3021         *pfpsf |= INVALID_EXCEPTION;
3022         // return Integer Indefinite
3023         res = 0x80000000;
3024       } else { // x is QNaN
3025         // set invalid flag
3026         *pfpsf |= INVALID_EXCEPTION;
3027         // return Integer Indefinite
3028         res = 0x80000000;
3029       }
3030       BID_RETURN (res);
3031     } else { // x is not a NaN, so it must be infinity
3032       if (!x_sign) { // x is +inf
3033         // set invalid flag
3034         *pfpsf |= INVALID_EXCEPTION;
3035         // return Integer Indefinite
3036         res = 0x80000000;
3037       } else { // x is -inf
3038         // set invalid flag
3039         *pfpsf |= INVALID_EXCEPTION;
3040         // return Integer Indefinite
3041         res = 0x80000000;
3042       }
3043       BID_RETURN (res);
3044     }
3045   }
3046   // check for non-canonical values (after the check for special values)
3047   if ((C1.w[1] > 0x0001ed09bead87c0ull)
3048       || (C1.w[1] == 0x0001ed09bead87c0ull
3049           && (C1.w[0] > 0x378d8e63ffffffffull))
3050       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
3051     res = 0x00000000;
3052     BID_RETURN (res);
3053   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
3054     // x is 0
3055     res = 0x00000000;
3056     BID_RETURN (res);
3057   } else { // x is not special and is not zero
3058
3059     // q = nr. of decimal digits in x
3060     //  determine first the nr. of bits in x
3061     if (C1.w[1] == 0) {
3062       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
3063         // split the 64-bit value in two 32-bit halves to avoid rounding errors
3064         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
3065           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
3066           x_nr_bits =
3067             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
3068         } else { // x < 2^32
3069           tmp1.d = (double) (C1.w[0]); // exact conversion
3070           x_nr_bits =
3071             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
3072         }
3073       } else { // if x < 2^53
3074         tmp1.d = (double) C1.w[0]; // exact conversion
3075         x_nr_bits =
3076           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
3077       }
3078     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
3079       tmp1.d = (double) C1.w[1]; // exact conversion
3080       x_nr_bits =
3081         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
3082     }
3083     q = __bid_nr_digits[x_nr_bits - 1].digits;
3084     if (q == 0) {
3085       q = __bid_nr_digits[x_nr_bits - 1].digits1;
3086       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
3087           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
3088           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
3089         q++;
3090     }
3091     exp = (x_exp >> 49) - 6176;
3092     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
3093       // set invalid flag
3094       *pfpsf |= INVALID_EXCEPTION;
3095       // return Integer Indefinite
3096       res = 0x80000000;
3097       BID_RETURN (res);
3098     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
3099       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
3100       // so x rounded to an integer may or may not fit in a signed 32-bit int
3101       // the cases that do not fit are identified here; the ones that fit
3102       // fall through and will be handled with other cases further,
3103       // under '1 <= q + exp <= 10'
3104       if (x_sign) { // if n < 0 and q + exp = 10
3105         // if n <= -2^31 - 1/2 then n is too large
3106         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1/2
3107         // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000005, 1<=q<=34
3108         if (q <= 11) {
3109           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
3110           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
3111           if (tmp64 >= 0x500000005ull) {
3112             // set invalid flag
3113             *pfpsf |= INVALID_EXCEPTION;
3114             // return Integer Indefinite
3115             res = 0x80000000;
3116             BID_RETURN (res);
3117           }
3118           // else cases that can be rounded to a 32-bit int fall through
3119           // to '1 <= q + exp <= 10'
3120         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
3121           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000005 <=>
3122           // C >= 0x500000005 * 10^(q-11) where 1 <= q - 11 <= 23
3123           // (scale 2^31+1/2 up)
3124           tmp64 = 0x500000005ull;
3125           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
3126             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
3127           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
3128             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
3129           }
3130           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
3131             // set invalid flag
3132             *pfpsf |= INVALID_EXCEPTION;
3133             // return Integer Indefinite
3134             res = 0x80000000;
3135             BID_RETURN (res);
3136           }
3137           // else cases that can be rounded to a 32-bit int fall through
3138           // to '1 <= q + exp <= 10'
3139         }
3140       } else { // if n > 0 and q + exp = 10
3141         // if n >= 2^31 - 1/2 then n is too large
3142         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31-1/2
3143         // too large if 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb, 1<=q<=34
3144         if (q <= 11) {
3145           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
3146           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
3147           if (tmp64 >= 0x4fffffffbull) {
3148             // set invalid flag
3149             *pfpsf |= INVALID_EXCEPTION;
3150             // return Integer Indefinite
3151             res = 0x80000000;
3152             BID_RETURN (res);
3153           }
3154           // else cases that can be rounded to a 32-bit int fall through
3155           // to '1 <= q + exp <= 10'
3156         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
3157           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb <=>
3158           // C >= 0x4fffffffb * 10^(q-11) where 1 <= q - 11 <= 23
3159           // (scale 2^31-1/2 up)
3160           tmp64 = 0x4fffffffbull;
3161           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
3162             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
3163           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
3164             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
3165           }
3166           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
3167             // set invalid flag
3168             *pfpsf |= INVALID_EXCEPTION;
3169             // return Integer Indefinite
3170             res = 0x80000000;
3171             BID_RETURN (res);
3172           }
3173           // else cases that can be rounded to a 32-bit int fall through
3174           // to '1 <= q + exp <= 10'
3175         }
3176       }
3177     }
3178     // n is not too large to be converted to int32: -2^31 - 1/2 < n < 2^31 - 1/2
3179     // Note: some of the cases tested for above fall through to this point
3180     if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1)
3181       // return 0
3182       res = 0x00000000;
3183       BID_RETURN (res);
3184     } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1)
3185       // if 0.c(0)c(1)...c(q-1) < 0.5 <=> c(0)c(1)...c(q-1) < 5 * 10^(q-1)
3186       //   res = 0
3187       // else
3188       //   res = +/-1
3189       ind = q - 1;
3190       if (ind <= 18) { // 0 <= ind <= 18
3191         if ((C1.w[1] == 0) && (C1.w[0] < __bid_midpoint64[ind])) {
3192           res = 0x00000000; // return 0
3193         } else if (x_sign) { // n < 0
3194           res = 0xffffffff; // return -1
3195         } else { // n > 0
3196           res = 0x00000001; // return +1
3197         }
3198       } else { // 19 <= ind <= 33
3199         if ((C1.w[1] < __bid_midpoint128[ind - 19].w[1])
3200             || ((C1.w[1] == __bid_midpoint128[ind - 19].w[1])
3201             && (C1.w[0] < __bid_midpoint128[ind - 19].w[0]))) {
3202           res = 0x00000000; // return 0
3203         } else if (x_sign) { // n < 0
3204           res = 0xffffffff; // return -1
3205         } else { // n > 0
3206           res = 0x00000001; // return +1
3207         }
3208       }
3209     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
3210       // -2^31-1/2 < x <= -1 or 1 <= x < 2^31-1/2 so x can be rounded
3211       // to nearest-away to a 32-bit signed integer
3212       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
3213         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
3214         // chop off ind digits from the lower part of C1
3215         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
3216         tmp64 = C1.w[0];
3217         if (ind <= 19) {
3218           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
3219         } else {
3220           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
3221           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
3222         }
3223         if (C1.w[0] < tmp64)
3224           C1.w[1]++;
3225         // calculate C* and f*
3226         // C* is actually floor(C*) in this case
3227         // C* and f* need shifting and masking, as shown by
3228         // __bid_shiftright128[] and __bid_maskhigh128[]
3229         // 1 <= x <= 33
3230         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
3231         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
3232         // the approximation of 10^(-x) was rounded up to 118 bits
3233         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
3234         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
3235           Cstar.w[1] = P256.w[3];
3236           Cstar.w[0] = P256.w[2];
3237         } else { // 22 <= ind - 1 <= 33
3238           Cstar.w[1] = 0;
3239           Cstar.w[0] = P256.w[3];
3240         }
3241         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
3242         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
3243         // if (0 < f* < 10^(-x)) then the result is a midpoint
3244         //   if floor(C*) is even then C* = floor(C*) - logical right
3245         //       shift; C* has p decimal digits, correct by Prop. 1)
3246         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
3247         //       shift; C* has p decimal digits, correct by Pr. 1)
3248         // else
3249         //   C* = floor(C*) (logical right shift; C has p decimal digits,
3250         //       correct by Property 1)
3251         // n = C* * 10^(e+x)
3252
3253         // shift right C* by Ex-128 = __bid_shiftright128[ind]
3254         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
3255         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
3256           Cstar.w[0] =
3257             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
3258           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
3259         } else { // 22 <= ind - 1 <= 33
3260           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
3261         }
3262         // if the result was a midpoint, it was already rounded away from zero
3263         if (x_sign)
3264           res = -Cstar.w[0];
3265         else
3266           res = Cstar.w[0];
3267         // no need to check for midpoints - already rounded away from zero!
3268       } else if (exp == 0) {
3269         // 1 <= q <= 10
3270         // res = +/-C (exact)
3271         if (x_sign)
3272           res = -C1.w[0];
3273         else
3274           res = C1.w[0];
3275       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
3276         // res = +/-C * 10^exp (exact)
3277         if (x_sign)
3278           res = -C1.w[0] * __bid_ten2k64[exp];
3279         else
3280           res = C1.w[0] * __bid_ten2k64[exp];
3281       }
3282     }
3283   }
3284   BID_RETURN (res);
3285 }
3286
3287 /*****************************************************************************
3288  *  BID128_to_int32_xrninta
3289  ****************************************************************************/
3290
3291 BID128_FUNCTION_ARG1_NORND_CUSTOMRESTYPE(int, __bid128_to_int32_xrninta, x)
3292
3293   int res;
3294   UINT64 x_sign;
3295   UINT64 x_exp;
3296   int exp;      // unbiased exponent
3297   // Note: C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (all are UINT64)
3298   UINT64 tmp64, tmp64A;
3299   BID_UI64DOUBLE tmp1;
3300   unsigned int x_nr_bits;
3301   int q, ind, shift;
3302   UINT128 C1, C;
3303   UINT128 Cstar; // C* represents up to 34 decimal digits ~ 113 bits
3304   UINT256 fstar;
3305   UINT256 P256;
3306
3307   // unpack x
3308   x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative
3309   x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bit positions
3310   C1.w[1] = x.w[1] & MASK_COEFF;
3311   C1.w[0] = x.w[0];
3312
3313   // check for NaN or Infinity
3314   if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) {
3315     // x is special
3316     if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN
3317       if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN
3318         // set invalid flag
3319         *pfpsf |= INVALID_EXCEPTION;
3320         // return Integer Indefinite
3321         res = 0x80000000;
3322       } else { // x is QNaN
3323         // set invalid flag
3324         *pfpsf |= INVALID_EXCEPTION;
3325         // return Integer Indefinite
3326         res = 0x80000000;
3327       }
3328       BID_RETURN (res);
3329     } else { // x is not a NaN, so it must be infinity
3330       if (!x_sign) { // x is +inf
3331         // set invalid flag
3332         *pfpsf |= INVALID_EXCEPTION;
3333         // return Integer Indefinite
3334         res = 0x80000000;
3335       } else { // x is -inf
3336         // set invalid flag
3337         *pfpsf |= INVALID_EXCEPTION;
3338         // return Integer Indefinite
3339         res = 0x80000000;
3340       }
3341       BID_RETURN (res);
3342     }
3343   }
3344   // check for non-canonical values (after the check for special values)
3345   if ((C1.w[1] > 0x0001ed09bead87c0ull)
3346       || (C1.w[1] == 0x0001ed09bead87c0ull
3347           && (C1.w[0] > 0x378d8e63ffffffffull))
3348       || ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull)) {
3349     res = 0x00000000;
3350     BID_RETURN (res);
3351   } else if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) {
3352     // x is 0
3353     res = 0x00000000;
3354     BID_RETURN (res);
3355   } else { // x is not special and is not zero
3356
3357     // q = nr. of decimal digits in x
3358     //  determine first the nr. of bits in x
3359     if (C1.w[1] == 0) {
3360       if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
3361         // split the 64-bit value in two 32-bit halves to avoid rounding errors
3362         if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
3363           tmp1.d = (double) (C1.w[0] >> 32); // exact conversion
3364           x_nr_bits =
3365             33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
3366         } else { // x < 2^32
3367           tmp1.d = (double) (C1.w[0]); // exact conversion
3368           x_nr_bits =
3369             1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
3370         }
3371       } else { // if x < 2^53
3372         tmp1.d = (double) C1.w[0]; // exact conversion
3373         x_nr_bits =
3374           1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
3375       }
3376     } else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1])
3377       tmp1.d = (double) C1.w[1]; // exact conversion
3378       x_nr_bits =
3379         65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff);
3380     }
3381     q = __bid_nr_digits[x_nr_bits - 1].digits;
3382     if (q == 0) {
3383       q = __bid_nr_digits[x_nr_bits - 1].digits1;
3384       if (C1.w[1] > __bid_nr_digits[x_nr_bits - 1].threshold_hi
3385           || (C1.w[1] == __bid_nr_digits[x_nr_bits - 1].threshold_hi
3386           && C1.w[0] >= __bid_nr_digits[x_nr_bits - 1].threshold_lo))
3387         q++;
3388     }
3389     exp = (x_exp >> 49) - 6176;
3390     if ((q + exp) > 10) { // x >= 10^10 ~= 2^33.2... (cannot fit in 32 bits)
3391       // set invalid flag
3392       *pfpsf |= INVALID_EXCEPTION;
3393       // return Integer Indefinite
3394       res = 0x80000000;
3395       BID_RETURN (res);
3396     } else if ((q + exp) == 10) { // x = c(0)c(1)...c(9).c(10)...c(q-1)
3397       // in this case 2^29.89... ~= 10^9 <= x < 10^10 ~= 2^33.2...
3398       // so x rounded to an integer may or may not fit in a signed 32-bit int
3399       // the cases that do not fit are identified here; the ones that fit
3400       // fall through and will be handled with other cases further,
3401       // under '1 <= q + exp <= 10'
3402       if (x_sign) { // if n < 0 and q + exp = 10
3403         // if n <= -2^31 - 1/2 then n is too large
3404         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31+1/2
3405         // <=> 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000005, 1<=q<=34
3406         if (q <= 11) {
3407           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
3408           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
3409           if (tmp64 >= 0x500000005ull) {
3410             // set invalid flag
3411             *pfpsf |= INVALID_EXCEPTION;
3412             // return Integer Indefinite
3413             res = 0x80000000;
3414             BID_RETURN (res);
3415           }
3416           // else cases that can be rounded to a 32-bit int fall through
3417           // to '1 <= q + exp <= 10'
3418         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
3419           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x500000005 <=>
3420           // C >= 0x500000005 * 10^(q-11) where 1 <= q - 11 <= 23
3421           // (scale 2^31+1/2 up)
3422           tmp64 = 0x500000005ull;
3423           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
3424             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
3425           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
3426             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
3427           }
3428           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
3429             // set invalid flag
3430             *pfpsf |= INVALID_EXCEPTION;
3431             // return Integer Indefinite
3432             res = 0x80000000;
3433             BID_RETURN (res);
3434           }
3435           // else cases that can be rounded to a 32-bit int fall through
3436           // to '1 <= q + exp <= 10'
3437         }
3438       } else { // if n > 0 and q + exp = 10
3439         // if n >= 2^31 - 1/2 then n is too large
3440         // too large if c(0)c(1)...c(9).c(10)...c(q-1) >= 2^31-1/2
3441         // too large if 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb, 1<=q<=34
3442         if (q <= 11) {
3443           tmp64 = C1.w[0] * __bid_ten2k64[11 - q]; // C scaled up to 11-digit int
3444           // c(0)c(1)...c(9)c(10) or c(0)c(1)...c(q-1)0...0 (11 digits)
3445           if (tmp64 >= 0x4fffffffbull) {
3446             // set invalid flag
3447             *pfpsf |= INVALID_EXCEPTION;
3448             // return Integer Indefinite
3449             res = 0x80000000;
3450             BID_RETURN (res);
3451           }
3452           // else cases that can be rounded to a 32-bit int fall through
3453           // to '1 <= q + exp <= 10'
3454         } else { // if (q > 11), i.e. 12 <= q <= 34 and so -24 <= exp <= -2
3455           // 0.c(0)c(1)...c(q-1) * 10^11 >= 0x4fffffffb <=>
3456           // C >= 0x4fffffffb * 10^(q-11) where 1 <= q - 11 <= 23
3457           // (scale 2^31-1/2 up)
3458           tmp64 = 0x4fffffffbull;
3459           if (q - 11 <= 19) { // 1 <= q - 11 <= 19; 10^(q-11) requires 64 bits
3460             __mul_64x64_to_128MACH (C, tmp64, __bid_ten2k64[q - 11]);
3461           } else { // 20 <= q - 11 <= 23, and 10^(q-11) requires 128 bits
3462             __mul_128x64_to_128 (C, tmp64, __bid_ten2k128[q - 31]);
3463           }
3464           if (C1.w[1] > C.w[1] || (C1.w[1] == C.w[1] && C1.w[0] >= C.w[0])) {
3465             // set invalid flag
3466             *pfpsf |= INVALID_EXCEPTION;
3467             // return Integer Indefinite
3468             res = 0x80000000;
3469             BID_RETURN (res);
3470           }
3471           // else cases that can be rounded to a 32-bit int fall through
3472           // to '1 <= q + exp <= 10'
3473         }
3474       }
3475     }
3476     // n is not too large to be converted to int32: -2^31 - 1/2 < n < 2^31 - 1/2
3477     // Note: some of the cases tested for above fall through to this point
3478     if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1)
3479       // set inexact flag
3480       *pfpsf |= INEXACT_EXCEPTION;
3481       // return 0
3482       res = 0x00000000;
3483       BID_RETURN (res);
3484     } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1)
3485       // if 0.c(0)c(1)...c(q-1) < 0.5 <=> c(0)c(1)...c(q-1) < 5 * 10^(q-1)
3486       //   res = 0
3487       // else
3488       //   res = +/-1
3489       ind = q - 1;
3490       if (ind <= 18) { // 0 <= ind <= 18
3491         if ((C1.w[1] == 0) && (C1.w[0] < __bid_midpoint64[ind])) {
3492           res = 0x00000000; // return 0
3493         } else if (x_sign) { // n < 0
3494           res = 0xffffffff; // return -1
3495         } else { // n > 0
3496           res = 0x00000001; // return +1
3497         }
3498       } else { // 19 <= ind <= 33
3499         if ((C1.w[1] < __bid_midpoint128[ind - 19].w[1])
3500             || ((C1.w[1] == __bid_midpoint128[ind - 19].w[1])
3501             && (C1.w[0] < __bid_midpoint128[ind - 19].w[0]))) {
3502           res = 0x00000000; // return 0
3503         } else if (x_sign) { // n < 0
3504           res = 0xffffffff; // return -1
3505         } else { // n > 0
3506           res = 0x00000001; // return +1
3507         }
3508       }
3509       // set inexact flag
3510       *pfpsf |= INEXACT_EXCEPTION;
3511     } else { // if (1 <= q + exp <= 10, 1 <= q <= 34, -33 <= exp <= 9)
3512       // -2^31-1/2 < x <= -1 or 1 <= x < 2^31-1/2 so x can be rounded
3513       // to nearest-away to a 32-bit signed integer
3514       if (exp < 0) { // 2 <= q <= 34, -33 <= exp <= -1, 1 <= q + exp <= 10
3515         ind = -exp; // 1 <= ind <= 33; ind is a synonym for 'x'
3516         // chop off ind digits from the lower part of C1
3517         // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 127 bits
3518         tmp64 = C1.w[0];
3519         if (ind <= 19) {
3520           C1.w[0] = C1.w[0] + __bid_midpoint64[ind - 1];
3521         } else {
3522           C1.w[0] = C1.w[0] + __bid_midpoint128[ind - 20].w[0];
3523           C1.w[1] = C1.w[1] + __bid_midpoint128[ind - 20].w[1];
3524         }
3525         if (C1.w[0] < tmp64)
3526           C1.w[1]++;
3527         // calculate C* and f*
3528         // C* is actually floor(C*) in this case
3529         // C* and f* need shifting and masking, as shown by
3530         // __bid_shiftright128[] and __bid_maskhigh128[]
3531         // 1 <= x <= 33
3532         // kx = 10^(-x) = __bid_ten2mk128[ind - 1]
3533         // C* = (C1 + 1/2 * 10^x) * 10^(-x)
3534         // the approximation of 10^(-x) was rounded up to 118 bits
3535         __mul_128x128_to_256 (P256, C1, __bid_ten2mk128[ind - 1]);
3536         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
3537           Cstar.w[1] = P256.w[3];
3538           Cstar.w[0] = P256.w[2];
3539           fstar.w[3] = 0;
3540           fstar.w[2] = P256.w[2] & __bid_maskhigh128[ind - 1];
3541           fstar.w[1] = P256.w[1];
3542           fstar.w[0] = P256.w[0];
3543         } else { // 22 <= ind - 1 <= 33
3544           Cstar.w[1] = 0;
3545           Cstar.w[0] = P256.w[3];
3546           fstar.w[3] = P256.w[3] & __bid_maskhigh128[ind - 1];
3547           fstar.w[2] = P256.w[2];
3548           fstar.w[1] = P256.w[1];
3549           fstar.w[0] = P256.w[0];
3550         }
3551         // the top Ex bits of 10^(-x) are T* = __bid_ten2mk128trunc[ind], e.g.
3552         // if x=1, T*=__bid_ten2mk128trunc[0]=0x19999999999999999999999999999999
3553         // if (0 < f* < 10^(-x)) then the result is a midpoint
3554         //   if floor(C*) is even then C* = floor(C*) - logical right
3555         //       shift; C* has p decimal digits, correct by Prop. 1)
3556         //   else if floor(C*) is odd C* = floor(C*)-1 (logical right
3557         //       shift; C* has p decimal digits, correct by Pr. 1)
3558         // else
3559         //   C* = floor(C*) (logical right shift; C has p decimal digits,
3560         //       correct by Property 1)
3561         // n = C* * 10^(e+x)
3562
3563         // shift right C* by Ex-128 = __bid_shiftright128[ind]
3564         shift = __bid_shiftright128[ind - 1]; // 0 <= shift <= 102
3565         if (ind - 1 <= 21) { // 0 <= ind - 1 <= 21
3566           Cstar.w[0] =
3567             (Cstar.w[0] >> shift) | (Cstar.w[1] << (64 - shift));
3568           // redundant, it will be 0! Cstar.w[1] = (Cstar.w[1] >> shift);
3569         } else { // 22 <= ind - 1 <= 33
3570           Cstar.w[0] = (Cstar.w[0] >> (shift - 64)); // 2 <= shift - 64 <= 38
3571         }
3572         // if the result was a midpoint, it was already rounded away from zero
3573         if (x_sign)
3574           res = -Cstar.w[0];
3575         else
3576           res = Cstar.w[0];
3577         // determine inexactness of the rounding of C*
3578         // if (0 < f* - 1/2 < 10^(-x)) then
3579         //   the result is exact
3580         // else // if (f* - 1/2 > T*) then
3581         //   the result is inexact
3582         if (ind - 1 <= 2) {
3583           if (fstar.w[1] > 0x8000000000000000ull ||
3584               (fstar.w[1] == 0x8000000000000000ull &&
3585               fstar.w[0] > 0x0ull)) { // f* > 1/2 and the result may be exact
3586             tmp64 = fstar.w[1] - 0x8000000000000000ull; // f* - 1/2
3587             if ((tmp64 > __bid_ten2mk128trunc[ind - 1].w[1]
3588                  || (tmp64 == __bid_ten2mk128trunc[ind - 1].w[1]
3589                  && fstar.w[0] >= __bid_ten2mk128trunc[ind - 1].w[0]))) {
3590               // set the inexact flag
3591               *pfpsf |= INEXACT_EXCEPTION;
3592             } // else the result is exact
3593           } else { // the result is inexact; f2* <= 1/2
3594             // set the inexact flag
3595             *pfpsf |= INEXACT_EXCEPTION;
3596           }
3597         } else if (ind - 1 <= 21) { // if 3 <= ind <= 21
3598           if (fstar.w[3] > 0x0 || (fstar.w[3] == 0x0
3599               && fstar.w[2] > __bid_one_half128[ind - 1]) || (fstar.w[3] == 0x0
3600               && fstar.w[2] == __bid_one_half128[ind - 1] && (fstar.w[1]
3601                                                        || fstar.w[0]))) {
3602             // f2* > 1/2 and the result may be exact
3603             // Calculate f2* - 1/2
3604             tmp64 = fstar.w[2] - __bid_one_half128[ind - 1];
3605             tmp64A = fstar.w[3];
3606             if (tmp64 > fstar.w[2])
3607               tmp64A--;
3608             if (tmp64A || tmp64
3609                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
3610                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
3611                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
3612               // set the inexact flag
3613               *pfpsf |= INEXACT_EXCEPTION;
3614             } // else the result is exact
3615           } else { // the result is inexact; f2* <= 1/2
3616             // set the inexact flag
3617             *pfpsf |= INEXACT_EXCEPTION;
3618           }
3619         } else { // if 22 <= ind <= 33
3620           if (fstar.w[3] > __bid_one_half128[ind - 1]
3621               || (fstar.w[3] == __bid_one_half128[ind - 1] && (fstar.w[2]
3622                                                        || fstar.w[1]
3623                                                        || fstar.w[0]))) {
3624             // f2* > 1/2 and the result may be exact
3625             // Calculate f2* - 1/2
3626             tmp64 = fstar.w[3] - __bid_one_half128[ind - 1];
3627             if (tmp64 || fstar.w[2]
3628                 || fstar.w[1] > __bid_ten2mk128trunc[ind - 1].w[1]
3629                 || (fstar.w[1] == __bid_ten2mk128trunc[ind - 1].w[1]
3630                 && fstar.w[0] > __bid_ten2mk128trunc[ind - 1].w[0])) {
3631               // set the inexact flag
3632               *pfpsf |= INEXACT_EXCEPTION;
3633             } // else the result is exact
3634           } else { // the result is inexact; f2* <= 1/2
3635             // set the inexact flag
3636             *pfpsf |= INEXACT_EXCEPTION;
3637           }
3638         }
3639         // no need to check for midpoints - already rounded away from zero!
3640       } else if (exp == 0) {
3641         // 1 <= q <= 10
3642         // res = +/-C (exact)
3643         if (x_sign)
3644           res = -C1.w[0];
3645         else
3646           res = C1.w[0];
3647       } else { // if (exp > 0) => 1 <= exp <= 9, 1 <= q < 9, 2 <= q + exp <= 10
3648         // res = +/-C * 10^exp (exact)
3649         if (x_sign)
3650           res = -C1.w[0] * __bid_ten2k64[exp];
3651         else
3652           res = C1.w[0] * __bid_ten2k64[exp];
3653       }
3654     }
3655   }
3656   BID_RETURN (res);
3657 }