Subversion Repositories Kolibri OS

Rev

Blame | Last modification | View Log | RSS feed

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