OpenCores
URL https://opencores.org/ocsvn/openrisc/openrisc/trunk

Subversion Repositories openrisc

[/] [openrisc/] [trunk/] [gnu-stable/] [gcc-4.5.1/] [libgcc/] [config/] [libbid/] [bid128_to_int64.c] - Blame information for rev 853

Go to most recent revision | Details | Compare with Previous | View Log

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

powered by: WebSVN 2.1.0

© copyright 1999-2024 OpenCores.org, equivalent to Oliscience, all rights reserved. OpenCores®, registered trademark.