]>
Commit | Line | Data |
---|---|---|
748086b7 | 1 | /* Copyright (C) 2007, 2009 Free Software Foundation, Inc. |
200359e8 L |
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 | |
748086b7 | 7 | Software Foundation; either version 3, or (at your option) any later |
200359e8 L |
8 | version. |
9 | ||
200359e8 L |
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 | ||
748086b7 JJ |
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/>. */ | |
200359e8 L |
23 | |
24 | #include "bid_internal.h" | |
25 | ||
26 | /***************************************************************************** | |
27 | * BID64_to_int64_rnint | |
28 | ****************************************************************************/ | |
29 | ||
30 | #if DECIMAL_CALL_BY_REFERENCE | |
31 | void | |
b2a00c89 | 32 | bid64_to_int64_rnint (SINT64 * pres, UINT64 * px |
200359e8 L |
33 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
34 | _EXC_INFO_PARAM) { | |
35 | UINT64 x = *px; | |
36 | #else | |
37 | SINT64 | |
b2a00c89 | 38 | bid64_to_int64_rnint (UINT64 x |
200359e8 L |
39 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
40 | _EXC_INFO_PARAM) { | |
41 | #endif | |
42 | SINT64 res; | |
43 | UINT64 x_sign; | |
44 | UINT64 x_exp; | |
b2a00c89 | 45 | int exp; // unbiased exponent |
200359e8 L |
46 | // Note: C1 represents x_significand (UINT64) |
47 | BID_UI64DOUBLE tmp1; | |
48 | unsigned int x_nr_bits; | |
49 | int q, ind, shift; | |
50 | UINT64 C1; | |
51 | UINT128 C; | |
b2a00c89 | 52 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
53 | UINT128 fstar; |
54 | UINT128 P128; | |
55 | ||
56 | // check for NaN or Infinity | |
57 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
58 | // set invalid flag | |
59 | *pfpsf |= INVALID_EXCEPTION; | |
60 | // return Integer Indefinite | |
61 | res = 0x8000000000000000ull; | |
62 | BID_RETURN (res); | |
63 | } | |
64 | // unpack x | |
b2a00c89 | 65 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
66 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
67 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 68 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 69 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 70 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
71 | x_exp = 0; |
72 | C1 = 0; | |
73 | } | |
74 | } else { | |
b2a00c89 | 75 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
76 | C1 = x & MASK_BINARY_SIG1; |
77 | } | |
78 | ||
79 | // check for zeros (possibly from non-canonical values) | |
80 | if (C1 == 0x0ull) { | |
81 | // x is 0 | |
82 | res = 0x00000000; | |
83 | BID_RETURN (res); | |
84 | } | |
85 | // x is not special and is not zero | |
86 | ||
87 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
88 | // determine first the nr. of bits in x | |
b2a00c89 | 89 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 90 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
91 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
92 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
93 | x_nr_bits = |
94 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
95 | } else { // x < 2^32 |
96 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
97 | x_nr_bits = |
98 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
99 | } | |
b2a00c89 L |
100 | } else { // if x < 2^53 |
101 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
102 | x_nr_bits = |
103 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
104 | } | |
b2a00c89 | 105 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 106 | if (q == 0) { |
b2a00c89 L |
107 | q = nr_digits[x_nr_bits - 1].digits1; |
108 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
109 | q++; |
110 | } | |
b2a00c89 | 111 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 112 | |
b2a00c89 | 113 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
114 | // set invalid flag |
115 | *pfpsf |= INVALID_EXCEPTION; | |
116 | // return Integer Indefinite | |
117 | res = 0x8000000000000000ull; | |
118 | BID_RETURN (res); | |
b2a00c89 | 119 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
120 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
121 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
122 | // the cases that do not fit are identified here; the ones that fit | |
123 | // fall through and will be handled with other cases further, | |
124 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 125 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
126 | // if n < -2^63 - 1/2 then n is too large |
127 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) > 2^63+1/2 | |
128 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 > 5*(2^64+1), 1<=q<=16 | |
129 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 > 0x50000000000000005, 1<=q<=16 | |
130 | // <=> C * 10^(20-q) > 0x50000000000000005, 1<=q<=16 | |
131 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 132 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
133 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x50000000000000005, has 20 |
134 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] > 0x05ull)) { | |
135 | // set invalid flag | |
136 | *pfpsf |= INVALID_EXCEPTION; | |
137 | // return Integer Indefinite | |
138 | res = 0x8000000000000000ull; | |
139 | BID_RETURN (res); | |
140 | } | |
141 | // else cases that can be rounded to a 64-bit int fall through | |
142 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 143 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
144 | // if n >= 2^63 - 1/2 then n is too large |
145 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63-1/2 | |
146 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64-1), 1<=q<=16 | |
147 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 0x4fffffffffffffffb, 1<=q<=16 | |
148 | // <=> if C * 10^(20-q) >= 0x4fffffffffffffffb, 1<=q<=16 | |
149 | C.w[1] = 0x0000000000000004ull; | |
150 | C.w[0] = 0xfffffffffffffffbull; | |
151 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 152 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
153 | if (C.w[1] > 0x04ull || |
154 | (C.w[1] == 0x04ull && C.w[0] >= 0xfffffffffffffffbull)) { | |
155 | // set invalid flag | |
156 | *pfpsf |= INVALID_EXCEPTION; | |
157 | // return Integer Indefinite | |
158 | res = 0x8000000000000000ull; | |
159 | BID_RETURN (res); | |
160 | } | |
161 | // else cases that can be rounded to a 64-bit int fall through | |
162 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
163 | } // end else if n > 0 and q + exp = 19 |
164 | } // end else if ((q + exp) == 19) | |
200359e8 L |
165 | |
166 | // n is not too large to be converted to int64: -2^63-1/2 <= n < 2^63-1/2 | |
167 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 168 | if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
169 | // return 0 |
170 | res = 0x0000000000000000ull; | |
171 | BID_RETURN (res); | |
b2a00c89 | 172 | } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1) |
200359e8 L |
173 | // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1) |
174 | // res = 0 | |
175 | // else | |
176 | // res = +/-1 | |
b2a00c89 L |
177 | ind = q - 1; // 0 <= ind <= 15 |
178 | if (C1 <= midpoint64[ind]) { | |
179 | res = 0x0000000000000000ull; // return 0 | |
180 | } else if (x_sign) { // n < 0 | |
181 | res = 0xffffffffffffffffull; // return -1 | |
182 | } else { // n > 0 | |
183 | res = 0x0000000000000001ull; // return +1 | |
200359e8 | 184 | } |
b2a00c89 | 185 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
186 | // -2^63-1/2 <= x <= -1 or 1 <= x < 2^63-1/2 so x can be rounded |
187 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
188 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
189 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
190 | // chop off ind digits from the lower part of C1 |
191 | // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 64 bits | |
b2a00c89 | 192 | C1 = C1 + midpoint64[ind - 1]; |
200359e8 L |
193 | // calculate C* and f* |
194 | // C* is actually floor(C*) in this case | |
195 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 196 | // shiftright128[] and maskhigh128[] |
200359e8 | 197 | // 1 <= x <= 15 |
b2a00c89 | 198 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
199 | // C* = (C1 + 1/2 * 10^x) * 10^(-x) |
200 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 201 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 202 | Cstar = P128.w[1]; |
b2a00c89 | 203 | fstar.w[1] = P128.w[1] & maskhigh128[ind - 1]; |
200359e8 | 204 | fstar.w[0] = P128.w[0]; |
b2a00c89 L |
205 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
206 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
207 | // if (0 < f* < 10^(-x)) then the result is a midpoint |
208 | // if floor(C*) is even then C* = floor(C*) - logical right | |
209 | // shift; C* has p decimal digits, correct by Prop. 1) | |
210 | // else if floor(C*) is odd C* = floor(C*)-1 (logical right | |
211 | // shift; C* has p decimal digits, correct by Pr. 1) | |
212 | // else | |
213 | // C* = floor(C*) (logical right shift; C has p decimal digits, | |
214 | // correct by Property 1) | |
215 | // n = C* * 10^(e+x) | |
216 | ||
b2a00c89 L |
217 | // shift right C* by Ex-64 = shiftright128[ind] |
218 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
219 | Cstar = Cstar >> shift; |
220 | ||
221 | // if the result was a midpoint it was rounded away from zero, so | |
222 | // it will need a correction | |
223 | // check for midpoints | |
224 | if ((fstar.w[1] == 0) && fstar.w[0] && | |
b2a00c89 L |
225 | (fstar.w[0] <= ten2mk128trunc[ind - 1].w[1])) { |
226 | // ten2mk128trunc[ind -1].w[1] is identical to | |
227 | // ten2mk128[ind -1].w[1] | |
200359e8 | 228 | // the result is a midpoint; round to nearest |
b2a00c89 | 229 | if (Cstar & 0x01) { // Cstar is odd; MP in [EVEN, ODD] |
200359e8 | 230 | // if floor(C*) is odd C = floor(C*) - 1; the result >= 1 |
b2a00c89 L |
231 | Cstar--; // Cstar is now even |
232 | } // else MP in [ODD, EVEN] | |
200359e8 L |
233 | } |
234 | if (x_sign) | |
235 | res = -Cstar; | |
236 | else | |
237 | res = Cstar; | |
238 | } else if (exp == 0) { | |
239 | // 1 <= q <= 16 | |
240 | // res = +/-C (exact) | |
241 | if (x_sign) | |
242 | res = -C1; | |
243 | else | |
244 | res = C1; | |
b2a00c89 | 245 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
246 | // (the upper limit of 20 on q + exp is due to the fact that |
247 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
248 | // res = +/-C * 10^exp (exact) | |
249 | if (x_sign) | |
b2a00c89 | 250 | res = -C1 * ten2k64[exp]; |
200359e8 | 251 | else |
b2a00c89 | 252 | res = C1 * ten2k64[exp]; |
200359e8 L |
253 | } |
254 | } | |
255 | BID_RETURN (res); | |
256 | } | |
257 | ||
258 | /***************************************************************************** | |
259 | * BID64_to_int64_xrnint | |
260 | ****************************************************************************/ | |
261 | ||
262 | #if DECIMAL_CALL_BY_REFERENCE | |
263 | void | |
b2a00c89 | 264 | bid64_to_int64_xrnint (SINT64 * pres, UINT64 * px |
200359e8 L |
265 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
266 | _EXC_INFO_PARAM) { | |
267 | UINT64 x = *px; | |
268 | #else | |
269 | SINT64 | |
b2a00c89 | 270 | bid64_to_int64_xrnint (UINT64 x |
200359e8 L |
271 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
272 | _EXC_INFO_PARAM) { | |
273 | #endif | |
274 | SINT64 res; | |
275 | UINT64 x_sign; | |
276 | UINT64 x_exp; | |
b2a00c89 | 277 | int exp; // unbiased exponent |
200359e8 L |
278 | // Note: C1 represents x_significand (UINT64) |
279 | UINT64 tmp64; | |
280 | BID_UI64DOUBLE tmp1; | |
281 | unsigned int x_nr_bits; | |
282 | int q, ind, shift; | |
283 | UINT64 C1; | |
284 | UINT128 C; | |
b2a00c89 | 285 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
286 | UINT128 fstar; |
287 | UINT128 P128; | |
288 | ||
289 | // check for NaN or Infinity | |
290 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
291 | // set invalid flag | |
292 | *pfpsf |= INVALID_EXCEPTION; | |
293 | // return Integer Indefinite | |
294 | res = 0x8000000000000000ull; | |
295 | BID_RETURN (res); | |
296 | } | |
297 | // unpack x | |
b2a00c89 | 298 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
299 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
300 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 301 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 302 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 303 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
304 | x_exp = 0; |
305 | C1 = 0; | |
306 | } | |
307 | } else { | |
b2a00c89 | 308 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
309 | C1 = x & MASK_BINARY_SIG1; |
310 | } | |
311 | ||
312 | // check for zeros (possibly from non-canonical values) | |
313 | if (C1 == 0x0ull) { | |
314 | // x is 0 | |
315 | res = 0x00000000; | |
316 | BID_RETURN (res); | |
317 | } | |
318 | // x is not special and is not zero | |
319 | ||
320 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
321 | // determine first the nr. of bits in x | |
b2a00c89 | 322 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 323 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
324 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
325 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
326 | x_nr_bits = |
327 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
328 | } else { // x < 2^32 |
329 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
330 | x_nr_bits = |
331 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
332 | } | |
b2a00c89 L |
333 | } else { // if x < 2^53 |
334 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
335 | x_nr_bits = |
336 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
337 | } | |
b2a00c89 | 338 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 339 | if (q == 0) { |
b2a00c89 L |
340 | q = nr_digits[x_nr_bits - 1].digits1; |
341 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
342 | q++; |
343 | } | |
b2a00c89 | 344 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 345 | |
b2a00c89 | 346 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
347 | // set invalid flag |
348 | *pfpsf |= INVALID_EXCEPTION; | |
349 | // return Integer Indefinite | |
350 | res = 0x8000000000000000ull; | |
351 | BID_RETURN (res); | |
b2a00c89 | 352 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
353 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
354 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
355 | // the cases that do not fit are identified here; the ones that fit | |
356 | // fall through and will be handled with other cases further, | |
357 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 358 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
359 | // if n < -2^63 - 1/2 then n is too large |
360 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) > 2^63+1/2 | |
361 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 > 5*(2^64+1), 1<=q<=16 | |
362 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 > 0x50000000000000005, 1<=q<=16 | |
363 | // <=> C * 10^(20-q) > 0x50000000000000005, 1<=q<=16 | |
364 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 365 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
366 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x50000000000000005, has 20 |
367 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] > 0x05ull)) { | |
368 | // set invalid flag | |
369 | *pfpsf |= INVALID_EXCEPTION; | |
370 | // return Integer Indefinite | |
371 | res = 0x8000000000000000ull; | |
372 | BID_RETURN (res); | |
373 | } | |
374 | // else cases that can be rounded to a 64-bit int fall through | |
375 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 376 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
377 | // if n >= 2^63 - 1/2 then n is too large |
378 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63-1/2 | |
379 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64-1), 1<=q<=16 | |
380 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 0x4fffffffffffffffb, 1<=q<=16 | |
381 | // <=> if C * 10^(20-q) >= 0x4fffffffffffffffb, 1<=q<=16 | |
382 | C.w[1] = 0x0000000000000004ull; | |
383 | C.w[0] = 0xfffffffffffffffbull; | |
384 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 385 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
386 | if (C.w[1] > 0x04ull || |
387 | (C.w[1] == 0x04ull && C.w[0] >= 0xfffffffffffffffbull)) { | |
388 | // set invalid flag | |
389 | *pfpsf |= INVALID_EXCEPTION; | |
390 | // return Integer Indefinite | |
391 | res = 0x8000000000000000ull; | |
392 | BID_RETURN (res); | |
393 | } | |
394 | // else cases that can be rounded to a 64-bit int fall through | |
395 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
396 | } // end else if n > 0 and q + exp = 19 |
397 | } // end else if ((q + exp) == 19) | |
200359e8 L |
398 | |
399 | // n is not too large to be converted to int64: -2^63-1/2 <= n < 2^63-1/2 | |
400 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 401 | if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
402 | // set inexact flag |
403 | *pfpsf |= INEXACT_EXCEPTION; | |
404 | // return 0 | |
405 | res = 0x0000000000000000ull; | |
406 | BID_RETURN (res); | |
b2a00c89 | 407 | } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1) |
200359e8 L |
408 | // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1) |
409 | // res = 0 | |
410 | // else | |
411 | // res = +/-1 | |
b2a00c89 L |
412 | ind = q - 1; // 0 <= ind <= 15 |
413 | if (C1 <= midpoint64[ind]) { | |
414 | res = 0x0000000000000000ull; // return 0 | |
415 | } else if (x_sign) { // n < 0 | |
416 | res = 0xffffffffffffffffull; // return -1 | |
417 | } else { // n > 0 | |
418 | res = 0x0000000000000001ull; // return +1 | |
200359e8 L |
419 | } |
420 | // set inexact flag | |
421 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 | 422 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
423 | // -2^63-1/2 <= x <= -1 or 1 <= x < 2^63-1/2 so x can be rounded |
424 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
425 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
426 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
427 | // chop off ind digits from the lower part of C1 |
428 | // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 64 bits | |
b2a00c89 | 429 | C1 = C1 + midpoint64[ind - 1]; |
200359e8 L |
430 | // calculate C* and f* |
431 | // C* is actually floor(C*) in this case | |
432 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 433 | // shiftright128[] and maskhigh128[] |
200359e8 | 434 | // 1 <= x <= 15 |
b2a00c89 | 435 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
436 | // C* = (C1 + 1/2 * 10^x) * 10^(-x) |
437 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 438 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 439 | Cstar = P128.w[1]; |
b2a00c89 | 440 | fstar.w[1] = P128.w[1] & maskhigh128[ind - 1]; |
200359e8 | 441 | fstar.w[0] = P128.w[0]; |
b2a00c89 L |
442 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
443 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
444 | // if (0 < f* < 10^(-x)) then the result is a midpoint |
445 | // if floor(C*) is even then C* = floor(C*) - logical right | |
446 | // shift; C* has p decimal digits, correct by Prop. 1) | |
447 | // else if floor(C*) is odd C* = floor(C*)-1 (logical right | |
448 | // shift; C* has p decimal digits, correct by Pr. 1) | |
449 | // else | |
450 | // C* = floor(C*) (logical right shift; C has p decimal digits, | |
451 | // correct by Property 1) | |
452 | // n = C* * 10^(e+x) | |
453 | ||
b2a00c89 L |
454 | // shift right C* by Ex-64 = shiftright128[ind] |
455 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
456 | Cstar = Cstar >> shift; |
457 | // determine inexactness of the rounding of C* | |
458 | // if (0 < f* - 1/2 < 10^(-x)) then | |
459 | // the result is exact | |
460 | // else // if (f* - 1/2 > T*) then | |
461 | // the result is inexact | |
462 | if (ind - 1 <= 2) { | |
463 | if (fstar.w[0] > 0x8000000000000000ull) { | |
464 | // f* > 1/2 and the result may be exact | |
b2a00c89 L |
465 | tmp64 = fstar.w[0] - 0x8000000000000000ull; // f* - 1/2 |
466 | if ((tmp64 > ten2mk128trunc[ind - 1].w[1])) { | |
467 | // ten2mk128trunc[ind -1].w[1] is identical to | |
468 | // ten2mk128[ind -1].w[1] | |
200359e8 L |
469 | // set the inexact flag |
470 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 L |
471 | } // else the result is exact |
472 | } else { // the result is inexact; f2* <= 1/2 | |
200359e8 L |
473 | // set the inexact flag |
474 | *pfpsf |= INEXACT_EXCEPTION; | |
475 | } | |
b2a00c89 L |
476 | } else { // if 3 <= ind - 1 <= 14 |
477 | if (fstar.w[1] > onehalf128[ind - 1] || | |
478 | (fstar.w[1] == onehalf128[ind - 1] && fstar.w[0])) { | |
200359e8 L |
479 | // f2* > 1/2 and the result may be exact |
480 | // Calculate f2* - 1/2 | |
b2a00c89 L |
481 | tmp64 = fstar.w[1] - onehalf128[ind - 1]; |
482 | if (tmp64 || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
483 | // ten2mk128trunc[ind -1].w[1] is identical to | |
484 | // ten2mk128[ind -1].w[1] | |
200359e8 L |
485 | // set the inexact flag |
486 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 L |
487 | } // else the result is exact |
488 | } else { // the result is inexact; f2* <= 1/2 | |
200359e8 L |
489 | // set the inexact flag |
490 | *pfpsf |= INEXACT_EXCEPTION; | |
491 | } | |
492 | } | |
493 | ||
494 | // if the result was a midpoint it was rounded away from zero, so | |
495 | // it will need a correction | |
496 | // check for midpoints | |
497 | if ((fstar.w[1] == 0) && fstar.w[0] && | |
b2a00c89 L |
498 | (fstar.w[0] <= ten2mk128trunc[ind - 1].w[1])) { |
499 | // ten2mk128trunc[ind -1].w[1] is identical to | |
500 | // ten2mk128[ind -1].w[1] | |
200359e8 | 501 | // the result is a midpoint; round to nearest |
b2a00c89 | 502 | if (Cstar & 0x01) { // Cstar is odd; MP in [EVEN, ODD] |
200359e8 | 503 | // if floor(C*) is odd C = floor(C*) - 1; the result >= 1 |
b2a00c89 L |
504 | Cstar--; // Cstar is now even |
505 | } // else MP in [ODD, EVEN] | |
200359e8 L |
506 | } |
507 | if (x_sign) | |
508 | res = -Cstar; | |
509 | else | |
510 | res = Cstar; | |
511 | } else if (exp == 0) { | |
512 | // 1 <= q <= 16 | |
513 | // res = +/-C (exact) | |
514 | if (x_sign) | |
515 | res = -C1; | |
516 | else | |
517 | res = C1; | |
b2a00c89 | 518 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
519 | // (the upper limit of 20 on q + exp is due to the fact that |
520 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
521 | // res = +/-C * 10^exp (exact) | |
522 | if (x_sign) | |
b2a00c89 | 523 | res = -C1 * ten2k64[exp]; |
200359e8 | 524 | else |
b2a00c89 | 525 | res = C1 * ten2k64[exp]; |
200359e8 L |
526 | } |
527 | } | |
528 | BID_RETURN (res); | |
529 | } | |
530 | ||
531 | /***************************************************************************** | |
532 | * BID64_to_int64_floor | |
533 | ****************************************************************************/ | |
534 | ||
535 | #if DECIMAL_CALL_BY_REFERENCE | |
536 | void | |
b2a00c89 | 537 | bid64_to_int64_floor (SINT64 * pres, UINT64 * px |
200359e8 L |
538 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
539 | _EXC_INFO_PARAM) { | |
540 | UINT64 x = *px; | |
541 | #else | |
542 | SINT64 | |
b2a00c89 | 543 | bid64_to_int64_floor (UINT64 x |
200359e8 L |
544 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
545 | _EXC_INFO_PARAM) { | |
546 | #endif | |
547 | SINT64 res; | |
548 | UINT64 x_sign; | |
549 | UINT64 x_exp; | |
b2a00c89 | 550 | int exp; // unbiased exponent |
200359e8 L |
551 | // Note: C1 represents x_significand (UINT64) |
552 | BID_UI64DOUBLE tmp1; | |
553 | unsigned int x_nr_bits; | |
554 | int q, ind, shift; | |
555 | UINT64 C1; | |
556 | UINT128 C; | |
b2a00c89 | 557 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
558 | UINT128 fstar; |
559 | UINT128 P128; | |
560 | ||
561 | // check for NaN or Infinity | |
562 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
563 | // set invalid flag | |
564 | *pfpsf |= INVALID_EXCEPTION; | |
565 | // return Integer Indefinite | |
566 | res = 0x8000000000000000ull; | |
567 | BID_RETURN (res); | |
568 | } | |
569 | // unpack x | |
b2a00c89 | 570 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
571 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
572 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 573 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 574 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 575 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
576 | x_exp = 0; |
577 | C1 = 0; | |
578 | } | |
579 | } else { | |
b2a00c89 | 580 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
581 | C1 = x & MASK_BINARY_SIG1; |
582 | } | |
583 | ||
584 | // check for zeros (possibly from non-canonical values) | |
585 | if (C1 == 0x0ull) { | |
586 | // x is 0 | |
587 | res = 0x0000000000000000ull; | |
588 | BID_RETURN (res); | |
589 | } | |
590 | // x is not special and is not zero | |
591 | ||
592 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
593 | // determine first the nr. of bits in x | |
b2a00c89 | 594 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 595 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
596 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
597 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
598 | x_nr_bits = |
599 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
600 | } else { // x < 2^32 |
601 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
602 | x_nr_bits = |
603 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
604 | } | |
b2a00c89 L |
605 | } else { // if x < 2^53 |
606 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
607 | x_nr_bits = |
608 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
609 | } | |
b2a00c89 | 610 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 611 | if (q == 0) { |
b2a00c89 L |
612 | q = nr_digits[x_nr_bits - 1].digits1; |
613 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
614 | q++; |
615 | } | |
b2a00c89 | 616 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 617 | |
b2a00c89 | 618 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
619 | // set invalid flag |
620 | *pfpsf |= INVALID_EXCEPTION; | |
621 | // return Integer Indefinite | |
622 | res = 0x8000000000000000ull; | |
623 | BID_RETURN (res); | |
b2a00c89 | 624 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
625 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
626 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
627 | // the cases that do not fit are identified here; the ones that fit | |
628 | // fall through and will be handled with other cases further, | |
629 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 630 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
631 | // if n < -2^63 then n is too large |
632 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) > 2^63 | |
633 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 5*2^64, 1<=q<=16 | |
634 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 > 0x50000000000000000, 1<=q<=16 | |
635 | // <=> C * 10^(20-q) > 0x50000000000000000, 1<=q<=16 | |
636 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 637 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
638 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x5000000000000000a, has 20 |
639 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] != 0)) { | |
640 | // set invalid flag | |
641 | *pfpsf |= INVALID_EXCEPTION; | |
642 | // return Integer Indefinite | |
643 | res = 0x8000000000000000ull; | |
644 | BID_RETURN (res); | |
645 | } | |
646 | // else cases that can be rounded to a 64-bit int fall through | |
647 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 648 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
649 | // if n >= 2^63 then n is too large |
650 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63 | |
651 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 5*2^64, 1<=q<=16 | |
652 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 0x50000000000000000, 1<=q<=16 | |
653 | // <=> if C * 10^(20-q) >= 0x50000000000000000, 1<=q<=16 | |
654 | C.w[1] = 0x0000000000000005ull; | |
655 | C.w[0] = 0x0000000000000000ull; | |
656 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 657 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
658 | if (C.w[1] >= 0x05ull) { |
659 | // actually C.w[1] == 0x05ull && C.w[0] >= 0x0000000000000000ull) { | |
660 | // set invalid flag | |
661 | *pfpsf |= INVALID_EXCEPTION; | |
662 | // return Integer Indefinite | |
663 | res = 0x8000000000000000ull; | |
664 | BID_RETURN (res); | |
665 | } | |
666 | // else cases that can be rounded to a 64-bit int fall through | |
667 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
668 | } // end else if n > 0 and q + exp = 19 |
669 | } // end else if ((q + exp) == 19) | |
200359e8 L |
670 | |
671 | // n is not too large to be converted to int64: -2^63 <= n < 2^63 | |
672 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 673 | if ((q + exp) <= 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
674 | // return -1 or 0 |
675 | if (x_sign) | |
676 | res = 0xffffffffffffffffull; | |
677 | else | |
678 | res = 0x0000000000000000ull; | |
679 | BID_RETURN (res); | |
b2a00c89 | 680 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
681 | // -2^63 <= x <= -1 or 1 <= x < 2^63 so x can be rounded |
682 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
683 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
684 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
685 | // chop off ind digits from the lower part of C1 |
686 | // C1 fits in 64 bits | |
687 | // calculate C* and f* | |
688 | // C* is actually floor(C*) in this case | |
689 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 690 | // shiftright128[] and maskhigh128[] |
200359e8 | 691 | // 1 <= x <= 15 |
b2a00c89 | 692 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
693 | // C* = C1 * 10^(-x) |
694 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 695 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 696 | Cstar = P128.w[1]; |
b2a00c89 | 697 | fstar.w[1] = P128.w[1] & maskhigh128[ind - 1]; |
200359e8 | 698 | fstar.w[0] = P128.w[0]; |
b2a00c89 L |
699 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
700 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
701 | // C* = floor(C*) (logical right shift; C has p decimal digits, |
702 | // correct by Property 1) | |
703 | // n = C* * 10^(e+x) | |
704 | ||
b2a00c89 L |
705 | // shift right C* by Ex-64 = shiftright128[ind] |
706 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
707 | Cstar = Cstar >> shift; |
708 | // determine inexactness of the rounding of C* | |
709 | // if (0 < f* < 10^(-x)) then | |
710 | // the result is exact | |
711 | // else // if (f* > T*) then | |
712 | // the result is inexact | |
b2a00c89 L |
713 | if (ind - 1 <= 2) { // fstar.w[1] is 0 |
714 | if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
715 | // ten2mk128trunc[ind -1].w[1] is identical to | |
716 | // ten2mk128[ind -1].w[1] | |
717 | if (x_sign) { // negative and inexact | |
200359e8 L |
718 | Cstar++; |
719 | } | |
b2a00c89 L |
720 | } // else the result is exact |
721 | } else { // if 3 <= ind - 1 <= 14 | |
722 | if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
723 | // ten2mk128trunc[ind -1].w[1] is identical to | |
724 | // ten2mk128[ind -1].w[1] | |
725 | if (x_sign) { // negative and inexact | |
200359e8 L |
726 | Cstar++; |
727 | } | |
b2a00c89 | 728 | } // else the result is exact |
200359e8 L |
729 | } |
730 | ||
731 | if (x_sign) | |
732 | res = -Cstar; | |
733 | else | |
734 | res = Cstar; | |
735 | } else if (exp == 0) { | |
736 | // 1 <= q <= 16 | |
737 | // res = +/-C (exact) | |
738 | if (x_sign) | |
739 | res = -C1; | |
740 | else | |
741 | res = C1; | |
b2a00c89 | 742 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
743 | // (the upper limit of 20 on q + exp is due to the fact that |
744 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
745 | // res = +/-C * 10^exp (exact) | |
746 | if (x_sign) | |
b2a00c89 | 747 | res = -C1 * ten2k64[exp]; |
200359e8 | 748 | else |
b2a00c89 | 749 | res = C1 * ten2k64[exp]; |
200359e8 L |
750 | } |
751 | } | |
752 | BID_RETURN (res); | |
753 | } | |
754 | ||
755 | /***************************************************************************** | |
756 | * BID64_to_int64_xfloor | |
757 | ****************************************************************************/ | |
758 | ||
759 | #if DECIMAL_CALL_BY_REFERENCE | |
760 | void | |
b2a00c89 | 761 | bid64_to_int64_xfloor (SINT64 * pres, UINT64 * px |
200359e8 L |
762 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
763 | _EXC_INFO_PARAM) { | |
764 | UINT64 x = *px; | |
765 | #else | |
766 | SINT64 | |
b2a00c89 | 767 | bid64_to_int64_xfloor (UINT64 x |
200359e8 L |
768 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
769 | _EXC_INFO_PARAM) { | |
770 | #endif | |
771 | SINT64 res; | |
772 | UINT64 x_sign; | |
773 | UINT64 x_exp; | |
b2a00c89 | 774 | int exp; // unbiased exponent |
200359e8 L |
775 | // Note: C1 represents x_significand (UINT64) |
776 | BID_UI64DOUBLE tmp1; | |
777 | unsigned int x_nr_bits; | |
778 | int q, ind, shift; | |
779 | UINT64 C1; | |
780 | UINT128 C; | |
b2a00c89 | 781 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
782 | UINT128 fstar; |
783 | UINT128 P128; | |
784 | ||
785 | // check for NaN or Infinity | |
786 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
787 | // set invalid flag | |
788 | *pfpsf |= INVALID_EXCEPTION; | |
789 | // return Integer Indefinite | |
790 | res = 0x8000000000000000ull; | |
791 | BID_RETURN (res); | |
792 | } | |
793 | // unpack x | |
b2a00c89 | 794 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
795 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
796 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 797 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 798 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 799 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
800 | x_exp = 0; |
801 | C1 = 0; | |
802 | } | |
803 | } else { | |
b2a00c89 | 804 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
805 | C1 = x & MASK_BINARY_SIG1; |
806 | } | |
807 | ||
808 | // check for zeros (possibly from non-canonical values) | |
809 | if (C1 == 0x0ull) { | |
810 | // x is 0 | |
811 | res = 0x0000000000000000ull; | |
812 | BID_RETURN (res); | |
813 | } | |
814 | // x is not special and is not zero | |
815 | ||
816 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
817 | // determine first the nr. of bits in x | |
b2a00c89 | 818 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 819 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
820 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
821 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
822 | x_nr_bits = |
823 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
824 | } else { // x < 2^32 |
825 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
826 | x_nr_bits = |
827 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
828 | } | |
b2a00c89 L |
829 | } else { // if x < 2^53 |
830 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
831 | x_nr_bits = |
832 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
833 | } | |
b2a00c89 | 834 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 835 | if (q == 0) { |
b2a00c89 L |
836 | q = nr_digits[x_nr_bits - 1].digits1; |
837 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
838 | q++; |
839 | } | |
b2a00c89 | 840 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 841 | |
b2a00c89 | 842 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
843 | // set invalid flag |
844 | *pfpsf |= INVALID_EXCEPTION; | |
845 | // return Integer Indefinite | |
846 | res = 0x8000000000000000ull; | |
847 | BID_RETURN (res); | |
b2a00c89 | 848 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
849 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
850 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
851 | // the cases that do not fit are identified here; the ones that fit | |
852 | // fall through and will be handled with other cases further, | |
853 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 854 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
855 | // if n < -2^63 then n is too large |
856 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) > 2^63 | |
857 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 5*2^64, 1<=q<=16 | |
858 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 > 0x50000000000000000, 1<=q<=16 | |
859 | // <=> C * 10^(20-q) > 0x50000000000000000, 1<=q<=16 | |
860 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 861 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
862 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x5000000000000000a, has 20 |
863 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] != 0)) { | |
864 | // set invalid flag | |
865 | *pfpsf |= INVALID_EXCEPTION; | |
866 | // return Integer Indefinite | |
867 | res = 0x8000000000000000ull; | |
868 | BID_RETURN (res); | |
869 | } | |
870 | // else cases that can be rounded to a 64-bit int fall through | |
871 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 872 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
873 | // if n >= 2^63 then n is too large |
874 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63 | |
875 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 5*2^64, 1<=q<=16 | |
876 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 0x50000000000000000, 1<=q<=16 | |
877 | // <=> if C * 10^(20-q) >= 0x50000000000000000, 1<=q<=16 | |
878 | C.w[1] = 0x0000000000000005ull; | |
879 | C.w[0] = 0x0000000000000000ull; | |
880 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 881 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
882 | if (C.w[1] >= 0x05ull) { |
883 | // actually C.w[1] == 0x05ull && C.w[0] >= 0x0000000000000000ull) { | |
884 | // set invalid flag | |
885 | *pfpsf |= INVALID_EXCEPTION; | |
886 | // return Integer Indefinite | |
887 | res = 0x8000000000000000ull; | |
888 | BID_RETURN (res); | |
889 | } | |
890 | // else cases that can be rounded to a 64-bit int fall through | |
891 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
892 | } // end else if n > 0 and q + exp = 19 |
893 | } // end else if ((q + exp) == 19) | |
200359e8 L |
894 | |
895 | // n is not too large to be converted to int64: -2^63 <= n < 2^63 | |
896 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 897 | if ((q + exp) <= 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
898 | // set inexact flag |
899 | *pfpsf |= INEXACT_EXCEPTION; | |
900 | // return -1 or 0 | |
901 | if (x_sign) | |
902 | res = 0xffffffffffffffffull; | |
903 | else | |
904 | res = 0x0000000000000000ull; | |
905 | BID_RETURN (res); | |
b2a00c89 | 906 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
907 | // -2^63 <= x <= -1 or 1 <= x < 2^63 so x can be rounded |
908 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
909 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
910 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
911 | // chop off ind digits from the lower part of C1 |
912 | // C1 fits in 64 bits | |
913 | // calculate C* and f* | |
914 | // C* is actually floor(C*) in this case | |
915 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 916 | // shiftright128[] and maskhigh128[] |
200359e8 | 917 | // 1 <= x <= 15 |
b2a00c89 | 918 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
919 | // C* = C1 * 10^(-x) |
920 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 921 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 922 | Cstar = P128.w[1]; |
b2a00c89 | 923 | fstar.w[1] = P128.w[1] & maskhigh128[ind - 1]; |
200359e8 | 924 | fstar.w[0] = P128.w[0]; |
b2a00c89 L |
925 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
926 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
927 | // C* = floor(C*) (logical right shift; C has p decimal digits, |
928 | // correct by Property 1) | |
929 | // n = C* * 10^(e+x) | |
930 | ||
b2a00c89 L |
931 | // shift right C* by Ex-64 = shiftright128[ind] |
932 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
933 | Cstar = Cstar >> shift; |
934 | // determine inexactness of the rounding of C* | |
935 | // if (0 < f* < 10^(-x)) then | |
936 | // the result is exact | |
937 | // else // if (f* > T*) then | |
938 | // the result is inexact | |
b2a00c89 L |
939 | if (ind - 1 <= 2) { // fstar.w[1] is 0 |
940 | if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
941 | // ten2mk128trunc[ind -1].w[1] is identical to | |
942 | // ten2mk128[ind -1].w[1] | |
943 | if (x_sign) { // negative and inexact | |
200359e8 L |
944 | Cstar++; |
945 | } | |
946 | // set the inexact flag | |
947 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 L |
948 | } // else the result is exact |
949 | } else { // if 3 <= ind - 1 <= 14 | |
950 | if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
951 | // ten2mk128trunc[ind -1].w[1] is identical to | |
952 | // ten2mk128[ind -1].w[1] | |
953 | if (x_sign) { // negative and inexact | |
200359e8 L |
954 | Cstar++; |
955 | } | |
956 | // set the inexact flag | |
957 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 | 958 | } // else the result is exact |
200359e8 L |
959 | } |
960 | ||
961 | if (x_sign) | |
962 | res = -Cstar; | |
963 | else | |
964 | res = Cstar; | |
965 | } else if (exp == 0) { | |
966 | // 1 <= q <= 16 | |
967 | // res = +/-C (exact) | |
968 | if (x_sign) | |
969 | res = -C1; | |
970 | else | |
971 | res = C1; | |
b2a00c89 | 972 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
973 | // (the upper limit of 20 on q + exp is due to the fact that |
974 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
975 | // res = +/-C * 10^exp (exact) | |
976 | if (x_sign) | |
b2a00c89 | 977 | res = -C1 * ten2k64[exp]; |
200359e8 | 978 | else |
b2a00c89 | 979 | res = C1 * ten2k64[exp]; |
200359e8 L |
980 | } |
981 | } | |
982 | BID_RETURN (res); | |
983 | } | |
984 | ||
985 | /***************************************************************************** | |
986 | * BID64_to_int64_ceil | |
987 | ****************************************************************************/ | |
988 | ||
989 | #if DECIMAL_CALL_BY_REFERENCE | |
990 | void | |
b2a00c89 | 991 | bid64_to_int64_ceil (SINT64 * pres, UINT64 * px |
200359e8 L |
992 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) |
993 | { | |
994 | UINT64 x = *px; | |
995 | #else | |
996 | SINT64 | |
b2a00c89 | 997 | bid64_to_int64_ceil (UINT64 x |
200359e8 L |
998 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) |
999 | { | |
1000 | #endif | |
1001 | SINT64 res; | |
1002 | UINT64 x_sign; | |
1003 | UINT64 x_exp; | |
b2a00c89 | 1004 | int exp; // unbiased exponent |
200359e8 L |
1005 | // Note: C1 represents x_significand (UINT64) |
1006 | BID_UI64DOUBLE tmp1; | |
1007 | unsigned int x_nr_bits; | |
1008 | int q, ind, shift; | |
1009 | UINT64 C1; | |
1010 | UINT128 C; | |
b2a00c89 | 1011 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
1012 | UINT128 fstar; |
1013 | UINT128 P128; | |
1014 | ||
1015 | // check for NaN or Infinity | |
1016 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
1017 | // set invalid flag | |
1018 | *pfpsf |= INVALID_EXCEPTION; | |
1019 | // return Integer Indefinite | |
1020 | res = 0x8000000000000000ull; | |
1021 | BID_RETURN (res); | |
1022 | } | |
1023 | // unpack x | |
b2a00c89 | 1024 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
1025 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
1026 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 1027 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 1028 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 1029 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
1030 | x_exp = 0; |
1031 | C1 = 0; | |
1032 | } | |
1033 | } else { | |
b2a00c89 | 1034 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
1035 | C1 = x & MASK_BINARY_SIG1; |
1036 | } | |
1037 | ||
1038 | // check for zeros (possibly from non-canonical values) | |
1039 | if (C1 == 0x0ull) { | |
1040 | // x is 0 | |
1041 | res = 0x00000000; | |
1042 | BID_RETURN (res); | |
1043 | } | |
1044 | // x is not special and is not zero | |
1045 | ||
1046 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
1047 | // determine first the nr. of bits in x | |
b2a00c89 | 1048 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 1049 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
1050 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
1051 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
1052 | x_nr_bits = |
1053 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
1054 | } else { // x < 2^32 |
1055 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1056 | x_nr_bits = |
1057 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1058 | } | |
b2a00c89 L |
1059 | } else { // if x < 2^53 |
1060 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1061 | x_nr_bits = |
1062 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1063 | } | |
b2a00c89 | 1064 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 1065 | if (q == 0) { |
b2a00c89 L |
1066 | q = nr_digits[x_nr_bits - 1].digits1; |
1067 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
1068 | q++; |
1069 | } | |
b2a00c89 | 1070 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 1071 | |
b2a00c89 | 1072 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
1073 | // set invalid flag |
1074 | *pfpsf |= INVALID_EXCEPTION; | |
1075 | // return Integer Indefinite | |
1076 | res = 0x8000000000000000ull; | |
1077 | BID_RETURN (res); | |
b2a00c89 | 1078 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
1079 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
1080 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
1081 | // the cases that do not fit are identified here; the ones that fit | |
1082 | // fall through and will be handled with other cases further, | |
1083 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 1084 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
1085 | // if n <= -2^63 - 1 then n is too large |
1086 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63+1 | |
1087 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64+2), 1<=q<=16 | |
1088 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 0x5000000000000000a, 1<=q<=16 | |
1089 | // <=> C * 10^(20-q) >= 0x5000000000000000a, 1<=q<=16 | |
1090 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1091 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1092 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x5000000000000000a, has 20 |
1093 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] >= 0x0aull)) { | |
1094 | // set invalid flag | |
1095 | *pfpsf |= INVALID_EXCEPTION; | |
1096 | // return Integer Indefinite | |
1097 | res = 0x8000000000000000ull; | |
1098 | BID_RETURN (res); | |
1099 | } | |
1100 | // else cases that can be rounded to a 64-bit int fall through | |
1101 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 1102 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
1103 | // if n > 2^63 - 1 then n is too large |
1104 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) > 2^63 - 1 | |
1105 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 > 5*(2^64-2), 1<=q<=16 | |
1106 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 > 0x4fffffffffffffff6, 1<=q<=16 | |
1107 | // <=> if C * 10^(20-q) > 0x4fffffffffffffff6, 1<=q<=16 | |
1108 | C.w[1] = 0x0000000000000004ull; | |
1109 | C.w[0] = 0xfffffffffffffff6ull; | |
1110 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1111 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1112 | if (C.w[1] > 0x04ull || |
1113 | (C.w[1] == 0x04ull && C.w[0] > 0xfffffffffffffff6ull)) { | |
1114 | // set invalid flag | |
1115 | *pfpsf |= INVALID_EXCEPTION; | |
1116 | // return Integer Indefinite | |
1117 | res = 0x8000000000000000ull; | |
1118 | BID_RETURN (res); | |
1119 | } | |
1120 | // else cases that can be rounded to a 64-bit int fall through | |
1121 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
1122 | } // end else if n > 0 and q + exp = 19 |
1123 | } // end else if ((q + exp) == 19) | |
200359e8 L |
1124 | |
1125 | // n is not too large to be converted to int64: -2^63-1 < n < 2^63 | |
1126 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 1127 | if ((q + exp) <= 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
1128 | // return 0 or 1 |
1129 | if (x_sign) | |
1130 | res = 0x00000000; | |
1131 | else | |
1132 | res = 0x00000001; | |
1133 | BID_RETURN (res); | |
b2a00c89 | 1134 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
1135 | // -2^63-1 < x <= -1 or 1 <= x <= 2^63 - 1 so x can be rounded |
1136 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
1137 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
1138 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
1139 | // chop off ind digits from the lower part of C1 |
1140 | // C1 fits in 64 bits | |
1141 | // calculate C* and f* | |
1142 | // C* is actually floor(C*) in this case | |
1143 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 1144 | // shiftright128[] and maskhigh128[] |
200359e8 | 1145 | // 1 <= x <= 15 |
b2a00c89 | 1146 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
1147 | // C* = C1 * 10^(-x) |
1148 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 1149 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 1150 | Cstar = P128.w[1]; |
b2a00c89 | 1151 | fstar.w[1] = P128.w[1] & maskhigh128[ind - 1]; |
200359e8 | 1152 | fstar.w[0] = P128.w[0]; |
b2a00c89 L |
1153 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
1154 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
1155 | // C* = floor(C*) (logical right shift; C has p decimal digits, |
1156 | // correct by Property 1) | |
1157 | // n = C* * 10^(e+x) | |
1158 | ||
b2a00c89 L |
1159 | // shift right C* by Ex-64 = shiftright128[ind] |
1160 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
1161 | Cstar = Cstar >> shift; |
1162 | // determine inexactness of the rounding of C* | |
1163 | // if (0 < f* < 10^(-x)) then | |
1164 | // the result is exact | |
1165 | // else // if (f* > T*) then | |
1166 | // the result is inexact | |
b2a00c89 L |
1167 | if (ind - 1 <= 2) { // fstar.w[1] is 0 |
1168 | if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
1169 | // ten2mk128trunc[ind -1].w[1] is identical to | |
1170 | // ten2mk128[ind -1].w[1] | |
1171 | if (!x_sign) { // positive and inexact | |
200359e8 L |
1172 | Cstar++; |
1173 | } | |
b2a00c89 L |
1174 | } // else the result is exact |
1175 | } else { // if 3 <= ind - 1 <= 14 | |
1176 | if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
1177 | // ten2mk128trunc[ind -1].w[1] is identical to | |
1178 | // ten2mk128[ind -1].w[1] | |
1179 | if (!x_sign) { // positive and inexact | |
200359e8 L |
1180 | Cstar++; |
1181 | } | |
b2a00c89 | 1182 | } // else the result is exact |
200359e8 L |
1183 | } |
1184 | ||
1185 | if (x_sign) | |
1186 | res = -Cstar; | |
1187 | else | |
1188 | res = Cstar; | |
1189 | } else if (exp == 0) { | |
1190 | // 1 <= q <= 16 | |
1191 | // res = +/-C (exact) | |
1192 | if (x_sign) | |
1193 | res = -C1; | |
1194 | else | |
1195 | res = C1; | |
b2a00c89 | 1196 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
1197 | // (the upper limit of 20 on q + exp is due to the fact that |
1198 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
1199 | // res = +/-C * 10^exp (exact) | |
1200 | if (x_sign) | |
b2a00c89 | 1201 | res = -C1 * ten2k64[exp]; |
200359e8 | 1202 | else |
b2a00c89 | 1203 | res = C1 * ten2k64[exp]; |
200359e8 L |
1204 | } |
1205 | } | |
1206 | BID_RETURN (res); | |
1207 | } | |
1208 | ||
1209 | /***************************************************************************** | |
1210 | * BID64_to_int64_xceil | |
1211 | ****************************************************************************/ | |
1212 | ||
1213 | #if DECIMAL_CALL_BY_REFERENCE | |
1214 | void | |
b2a00c89 | 1215 | bid64_to_int64_xceil (SINT64 * pres, UINT64 * px |
200359e8 L |
1216 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
1217 | _EXC_INFO_PARAM) { | |
1218 | UINT64 x = *px; | |
1219 | #else | |
1220 | SINT64 | |
b2a00c89 | 1221 | bid64_to_int64_xceil (UINT64 x |
200359e8 L |
1222 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
1223 | _EXC_INFO_PARAM) { | |
1224 | #endif | |
1225 | SINT64 res; | |
1226 | UINT64 x_sign; | |
1227 | UINT64 x_exp; | |
b2a00c89 | 1228 | int exp; // unbiased exponent |
200359e8 L |
1229 | // Note: C1 represents x_significand (UINT64) |
1230 | BID_UI64DOUBLE tmp1; | |
1231 | unsigned int x_nr_bits; | |
1232 | int q, ind, shift; | |
1233 | UINT64 C1; | |
1234 | UINT128 C; | |
b2a00c89 | 1235 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
1236 | UINT128 fstar; |
1237 | UINT128 P128; | |
1238 | ||
1239 | // check for NaN or Infinity | |
1240 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
1241 | // set invalid flag | |
1242 | *pfpsf |= INVALID_EXCEPTION; | |
1243 | // return Integer Indefinite | |
1244 | res = 0x8000000000000000ull; | |
1245 | BID_RETURN (res); | |
1246 | } | |
1247 | // unpack x | |
b2a00c89 | 1248 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
1249 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
1250 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 1251 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 1252 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 1253 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
1254 | x_exp = 0; |
1255 | C1 = 0; | |
1256 | } | |
1257 | } else { | |
b2a00c89 | 1258 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
1259 | C1 = x & MASK_BINARY_SIG1; |
1260 | } | |
1261 | ||
1262 | // check for zeros (possibly from non-canonical values) | |
1263 | if (C1 == 0x0ull) { | |
1264 | // x is 0 | |
1265 | res = 0x00000000; | |
1266 | BID_RETURN (res); | |
1267 | } | |
1268 | // x is not special and is not zero | |
1269 | ||
1270 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
1271 | // determine first the nr. of bits in x | |
b2a00c89 | 1272 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 1273 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
1274 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
1275 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
1276 | x_nr_bits = |
1277 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
1278 | } else { // x < 2^32 |
1279 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1280 | x_nr_bits = |
1281 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1282 | } | |
b2a00c89 L |
1283 | } else { // if x < 2^53 |
1284 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1285 | x_nr_bits = |
1286 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1287 | } | |
b2a00c89 | 1288 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 1289 | if (q == 0) { |
b2a00c89 L |
1290 | q = nr_digits[x_nr_bits - 1].digits1; |
1291 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
1292 | q++; |
1293 | } | |
b2a00c89 | 1294 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 1295 | |
b2a00c89 | 1296 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
1297 | // set invalid flag |
1298 | *pfpsf |= INVALID_EXCEPTION; | |
1299 | // return Integer Indefinite | |
1300 | res = 0x8000000000000000ull; | |
1301 | BID_RETURN (res); | |
b2a00c89 | 1302 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
1303 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
1304 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
1305 | // the cases that do not fit are identified here; the ones that fit | |
1306 | // fall through and will be handled with other cases further, | |
1307 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 1308 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
1309 | // if n <= -2^63 - 1 then n is too large |
1310 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63+1 | |
1311 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64+2), 1<=q<=16 | |
1312 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 0x5000000000000000a, 1<=q<=16 | |
1313 | // <=> C * 10^(20-q) >= 0x5000000000000000a, 1<=q<=16 | |
1314 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1315 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1316 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x5000000000000000a, has 20 |
1317 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] >= 0x0aull)) { | |
1318 | // set invalid flag | |
1319 | *pfpsf |= INVALID_EXCEPTION; | |
1320 | // return Integer Indefinite | |
1321 | res = 0x8000000000000000ull; | |
1322 | BID_RETURN (res); | |
1323 | } | |
1324 | // else cases that can be rounded to a 64-bit int fall through | |
1325 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 1326 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
1327 | // if n > 2^63 - 1 then n is too large |
1328 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) > 2^63 - 1 | |
1329 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 > 5*(2^64-2), 1<=q<=16 | |
1330 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 > 0x4fffffffffffffff6, 1<=q<=16 | |
1331 | // <=> if C * 10^(20-q) > 0x4fffffffffffffff6, 1<=q<=16 | |
1332 | C.w[1] = 0x0000000000000004ull; | |
1333 | C.w[0] = 0xfffffffffffffff6ull; | |
1334 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1335 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1336 | if (C.w[1] > 0x04ull || |
1337 | (C.w[1] == 0x04ull && C.w[0] > 0xfffffffffffffff6ull)) { | |
1338 | // set invalid flag | |
1339 | *pfpsf |= INVALID_EXCEPTION; | |
1340 | // return Integer Indefinite | |
1341 | res = 0x8000000000000000ull; | |
1342 | BID_RETURN (res); | |
1343 | } | |
1344 | // else cases that can be rounded to a 64-bit int fall through | |
1345 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
1346 | } // end else if n > 0 and q + exp = 19 |
1347 | } // end else if ((q + exp) == 19) | |
200359e8 L |
1348 | |
1349 | // n is not too large to be converted to int64: -2^63-1 < n < 2^63 | |
1350 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 1351 | if ((q + exp) <= 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
1352 | // set inexact flag |
1353 | *pfpsf |= INEXACT_EXCEPTION; | |
1354 | // return 0 or 1 | |
1355 | if (x_sign) | |
1356 | res = 0x00000000; | |
1357 | else | |
1358 | res = 0x00000001; | |
1359 | BID_RETURN (res); | |
b2a00c89 | 1360 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
1361 | // -2^63-1 < x <= -1 or 1 <= x <= 2^63 - 1 so x can be rounded |
1362 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
1363 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
1364 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
1365 | // chop off ind digits from the lower part of C1 |
1366 | // C1 fits in 64 bits | |
1367 | // calculate C* and f* | |
1368 | // C* is actually floor(C*) in this case | |
1369 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 1370 | // shiftright128[] and maskhigh128[] |
200359e8 | 1371 | // 1 <= x <= 15 |
b2a00c89 | 1372 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
1373 | // C* = C1 * 10^(-x) |
1374 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 1375 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 1376 | Cstar = P128.w[1]; |
b2a00c89 | 1377 | fstar.w[1] = P128.w[1] & maskhigh128[ind - 1]; |
200359e8 | 1378 | fstar.w[0] = P128.w[0]; |
b2a00c89 L |
1379 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
1380 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
1381 | // C* = floor(C*) (logical right shift; C has p decimal digits, |
1382 | // correct by Property 1) | |
1383 | // n = C* * 10^(e+x) | |
1384 | ||
b2a00c89 L |
1385 | // shift right C* by Ex-64 = shiftright128[ind] |
1386 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
1387 | Cstar = Cstar >> shift; |
1388 | // determine inexactness of the rounding of C* | |
1389 | // if (0 < f* < 10^(-x)) then | |
1390 | // the result is exact | |
1391 | // else // if (f* > T*) then | |
1392 | // the result is inexact | |
b2a00c89 L |
1393 | if (ind - 1 <= 2) { // fstar.w[1] is 0 |
1394 | if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
1395 | // ten2mk128trunc[ind -1].w[1] is identical to | |
1396 | // ten2mk128[ind -1].w[1] | |
1397 | if (!x_sign) { // positive and inexact | |
200359e8 L |
1398 | Cstar++; |
1399 | } | |
1400 | // set the inexact flag | |
1401 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 L |
1402 | } // else the result is exact |
1403 | } else { // if 3 <= ind - 1 <= 14 | |
1404 | if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
1405 | // ten2mk128trunc[ind -1].w[1] is identical to | |
1406 | // ten2mk128[ind -1].w[1] | |
1407 | if (!x_sign) { // positive and inexact | |
200359e8 L |
1408 | Cstar++; |
1409 | } | |
1410 | // set the inexact flag | |
1411 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 | 1412 | } // else the result is exact |
200359e8 L |
1413 | } |
1414 | ||
1415 | if (x_sign) | |
1416 | res = -Cstar; | |
1417 | else | |
1418 | res = Cstar; | |
1419 | } else if (exp == 0) { | |
1420 | // 1 <= q <= 16 | |
1421 | // res = +/-C (exact) | |
1422 | if (x_sign) | |
1423 | res = -C1; | |
1424 | else | |
1425 | res = C1; | |
b2a00c89 | 1426 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
1427 | // (the upper limit of 20 on q + exp is due to the fact that |
1428 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
1429 | // res = +/-C * 10^exp (exact) | |
1430 | if (x_sign) | |
b2a00c89 | 1431 | res = -C1 * ten2k64[exp]; |
200359e8 | 1432 | else |
b2a00c89 | 1433 | res = C1 * ten2k64[exp]; |
200359e8 L |
1434 | } |
1435 | } | |
1436 | BID_RETURN (res); | |
1437 | } | |
1438 | ||
1439 | /***************************************************************************** | |
1440 | * BID64_to_int64_int | |
1441 | ****************************************************************************/ | |
1442 | ||
1443 | #if DECIMAL_CALL_BY_REFERENCE | |
1444 | void | |
b2a00c89 | 1445 | bid64_to_int64_int (SINT64 * pres, UINT64 * px |
200359e8 L |
1446 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) { |
1447 | UINT64 x = *px; | |
1448 | #else | |
1449 | SINT64 | |
b2a00c89 | 1450 | bid64_to_int64_int (UINT64 x |
200359e8 L |
1451 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) { |
1452 | #endif | |
1453 | SINT64 res; | |
1454 | UINT64 x_sign; | |
1455 | UINT64 x_exp; | |
b2a00c89 | 1456 | int exp; // unbiased exponent |
200359e8 L |
1457 | // Note: C1 represents x_significand (UINT64) |
1458 | BID_UI64DOUBLE tmp1; | |
1459 | unsigned int x_nr_bits; | |
1460 | int q, ind, shift; | |
1461 | UINT64 C1; | |
1462 | UINT128 C; | |
b2a00c89 | 1463 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
1464 | UINT128 P128; |
1465 | ||
1466 | // check for NaN or Infinity | |
1467 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
1468 | // set invalid flag | |
1469 | *pfpsf |= INVALID_EXCEPTION; | |
1470 | // return Integer Indefinite | |
1471 | res = 0x8000000000000000ull; | |
1472 | BID_RETURN (res); | |
1473 | } | |
1474 | // unpack x | |
b2a00c89 | 1475 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
1476 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
1477 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 1478 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 1479 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 1480 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
1481 | x_exp = 0; |
1482 | C1 = 0; | |
1483 | } | |
1484 | } else { | |
b2a00c89 | 1485 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
1486 | C1 = x & MASK_BINARY_SIG1; |
1487 | } | |
1488 | ||
1489 | // check for zeros (possibly from non-canonical values) | |
1490 | if (C1 == 0x0ull) { | |
1491 | // x is 0 | |
1492 | res = 0x00000000; | |
1493 | BID_RETURN (res); | |
1494 | } | |
1495 | // x is not special and is not zero | |
1496 | ||
1497 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
1498 | // determine first the nr. of bits in x | |
b2a00c89 | 1499 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 1500 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
1501 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
1502 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
1503 | x_nr_bits = |
1504 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
1505 | } else { // x < 2^32 |
1506 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1507 | x_nr_bits = |
1508 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1509 | } | |
b2a00c89 L |
1510 | } else { // if x < 2^53 |
1511 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1512 | x_nr_bits = |
1513 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1514 | } | |
b2a00c89 | 1515 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 1516 | if (q == 0) { |
b2a00c89 L |
1517 | q = nr_digits[x_nr_bits - 1].digits1; |
1518 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
1519 | q++; |
1520 | } | |
b2a00c89 | 1521 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 1522 | |
b2a00c89 | 1523 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
1524 | // set invalid flag |
1525 | *pfpsf |= INVALID_EXCEPTION; | |
1526 | // return Integer Indefinite | |
1527 | res = 0x8000000000000000ull; | |
1528 | BID_RETURN (res); | |
b2a00c89 | 1529 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
1530 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
1531 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
1532 | // the cases that do not fit are identified here; the ones that fit | |
1533 | // fall through and will be handled with other cases further, | |
1534 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 1535 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
1536 | // if n <= -2^63 - 1 then n is too large |
1537 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63+1 | |
1538 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64+2), 1<=q<=16 | |
1539 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 0x5000000000000000a, 1<=q<=16 | |
1540 | // <=> C * 10^(20-q) >= 0x5000000000000000a, 1<=q<=16 | |
1541 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1542 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1543 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x5000000000000000a, has 20 |
1544 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] >= 0x0aull)) { | |
1545 | // set invalid flag | |
1546 | *pfpsf |= INVALID_EXCEPTION; | |
1547 | // return Integer Indefinite | |
1548 | res = 0x8000000000000000ull; | |
1549 | BID_RETURN (res); | |
1550 | } | |
1551 | // else cases that can be rounded to a 64-bit int fall through | |
1552 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 1553 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
1554 | // if n >= 2^63 then n is too large |
1555 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63 | |
1556 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 5*2^64, 1<=q<=16 | |
1557 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 0x50000000000000000, 1<=q<=16 | |
1558 | // <=> if C * 10^(20-q) >= 0x50000000000000000, 1<=q<=16 | |
1559 | C.w[1] = 0x0000000000000005ull; | |
1560 | C.w[0] = 0x0000000000000000ull; | |
1561 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1562 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1563 | if (C.w[1] >= 0x05ull) { |
1564 | // actually C.w[1] == 0x05ull && C.w[0] >= 0x0000000000000000ull) { | |
1565 | // set invalid flag | |
1566 | *pfpsf |= INVALID_EXCEPTION; | |
1567 | // return Integer Indefinite | |
1568 | res = 0x8000000000000000ull; | |
1569 | BID_RETURN (res); | |
1570 | } | |
1571 | // else cases that can be rounded to a 64-bit int fall through | |
1572 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
1573 | } // end else if n > 0 and q + exp = 19 |
1574 | } // end else if ((q + exp) == 19) | |
200359e8 L |
1575 | |
1576 | // n is not too large to be converted to int64: -2^63-1 < n < 2^63 | |
1577 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 1578 | if ((q + exp) <= 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
1579 | // return 0 |
1580 | res = 0x0000000000000000ull; | |
1581 | BID_RETURN (res); | |
b2a00c89 | 1582 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
1583 | // -2^63-1 < x <= -1 or 1 <= x < 2^63 so x can be rounded |
1584 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
1585 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
1586 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
1587 | // chop off ind digits from the lower part of C1 |
1588 | // C1 fits in 64 bits | |
1589 | // calculate C* and f* | |
1590 | // C* is actually floor(C*) in this case | |
1591 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 1592 | // shiftright128[] and maskhigh128[] |
200359e8 | 1593 | // 1 <= x <= 15 |
b2a00c89 | 1594 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
1595 | // C* = C1 * 10^(-x) |
1596 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 1597 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 1598 | Cstar = P128.w[1]; |
b2a00c89 L |
1599 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
1600 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
1601 | // C* = floor(C*) (logical right shift; C has p decimal digits, |
1602 | // correct by Property 1) | |
1603 | // n = C* * 10^(e+x) | |
1604 | ||
b2a00c89 L |
1605 | // shift right C* by Ex-64 = shiftright128[ind] |
1606 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
1607 | Cstar = Cstar >> shift; |
1608 | ||
1609 | if (x_sign) | |
1610 | res = -Cstar; | |
1611 | else | |
1612 | res = Cstar; | |
1613 | } else if (exp == 0) { | |
1614 | // 1 <= q <= 16 | |
1615 | // res = +/-C (exact) | |
1616 | if (x_sign) | |
1617 | res = -C1; | |
1618 | else | |
1619 | res = C1; | |
b2a00c89 | 1620 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
1621 | // (the upper limit of 20 on q + exp is due to the fact that |
1622 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
1623 | // res = +/-C * 10^exp (exact) | |
1624 | if (x_sign) | |
b2a00c89 | 1625 | res = -C1 * ten2k64[exp]; |
200359e8 | 1626 | else |
b2a00c89 | 1627 | res = C1 * ten2k64[exp]; |
200359e8 L |
1628 | } |
1629 | } | |
1630 | BID_RETURN (res); | |
1631 | } | |
1632 | ||
1633 | /***************************************************************************** | |
1634 | * BID64_to_int64_xint | |
1635 | ****************************************************************************/ | |
1636 | ||
1637 | #if DECIMAL_CALL_BY_REFERENCE | |
1638 | void | |
b2a00c89 | 1639 | bid64_to_int64_xint (SINT64 * pres, UINT64 * px |
200359e8 L |
1640 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) |
1641 | { | |
1642 | UINT64 x = *px; | |
1643 | #else | |
1644 | SINT64 | |
b2a00c89 | 1645 | bid64_to_int64_xint (UINT64 x |
200359e8 L |
1646 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) |
1647 | { | |
1648 | #endif | |
1649 | SINT64 res; | |
1650 | UINT64 x_sign; | |
1651 | UINT64 x_exp; | |
b2a00c89 | 1652 | int exp; // unbiased exponent |
200359e8 L |
1653 | // Note: C1 represents x_significand (UINT64) |
1654 | BID_UI64DOUBLE tmp1; | |
1655 | unsigned int x_nr_bits; | |
1656 | int q, ind, shift; | |
1657 | UINT64 C1; | |
1658 | UINT128 C; | |
b2a00c89 | 1659 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
1660 | UINT128 fstar; |
1661 | UINT128 P128; | |
1662 | ||
1663 | // check for NaN or Infinity | |
1664 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
1665 | // set invalid flag | |
1666 | *pfpsf |= INVALID_EXCEPTION; | |
1667 | // return Integer Indefinite | |
1668 | res = 0x8000000000000000ull; | |
1669 | BID_RETURN (res); | |
1670 | } | |
1671 | // unpack x | |
b2a00c89 | 1672 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
1673 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
1674 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 1675 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 1676 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 1677 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
1678 | x_exp = 0; |
1679 | C1 = 0; | |
1680 | } | |
1681 | } else { | |
b2a00c89 | 1682 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
1683 | C1 = x & MASK_BINARY_SIG1; |
1684 | } | |
1685 | ||
1686 | // check for zeros (possibly from non-canonical values) | |
1687 | if (C1 == 0x0ull) { | |
1688 | // x is 0 | |
1689 | res = 0x00000000; | |
1690 | BID_RETURN (res); | |
1691 | } | |
1692 | // x is not special and is not zero | |
1693 | ||
1694 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
1695 | // determine first the nr. of bits in x | |
b2a00c89 | 1696 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 1697 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
1698 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
1699 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
1700 | x_nr_bits = |
1701 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
1702 | } else { // x < 2^32 |
1703 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1704 | x_nr_bits = |
1705 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1706 | } | |
b2a00c89 L |
1707 | } else { // if x < 2^53 |
1708 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1709 | x_nr_bits = |
1710 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1711 | } | |
b2a00c89 | 1712 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 1713 | if (q == 0) { |
b2a00c89 L |
1714 | q = nr_digits[x_nr_bits - 1].digits1; |
1715 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
1716 | q++; |
1717 | } | |
b2a00c89 | 1718 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 1719 | |
b2a00c89 | 1720 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
1721 | // set invalid flag |
1722 | *pfpsf |= INVALID_EXCEPTION; | |
1723 | // return Integer Indefinite | |
1724 | res = 0x8000000000000000ull; | |
1725 | BID_RETURN (res); | |
b2a00c89 | 1726 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
1727 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
1728 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
1729 | // the cases that do not fit are identified here; the ones that fit | |
1730 | // fall through and will be handled with other cases further, | |
1731 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 1732 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
1733 | // if n <= -2^63 - 1 then n is too large |
1734 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63+1 | |
1735 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64+2), 1<=q<=16 | |
1736 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 0x5000000000000000a, 1<=q<=16 | |
1737 | // <=> C * 10^(20-q) >= 0x5000000000000000a, 1<=q<=16 | |
1738 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1739 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1740 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x5000000000000000a, has 20 |
1741 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] >= 0x0aull)) { | |
1742 | // set invalid flag | |
1743 | *pfpsf |= INVALID_EXCEPTION; | |
1744 | // return Integer Indefinite | |
1745 | res = 0x8000000000000000ull; | |
1746 | BID_RETURN (res); | |
1747 | } | |
1748 | // else cases that can be rounded to a 64-bit int fall through | |
1749 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 1750 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
1751 | // if n >= 2^63 then n is too large |
1752 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63 | |
1753 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 5*2^64, 1<=q<=16 | |
1754 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 0x50000000000000000, 1<=q<=16 | |
1755 | // <=> if C * 10^(20-q) >= 0x50000000000000000, 1<=q<=16 | |
1756 | C.w[1] = 0x0000000000000005ull; | |
1757 | C.w[0] = 0x0000000000000000ull; | |
1758 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1759 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1760 | if (C.w[1] >= 0x05ull) { |
1761 | // actually C.w[1] == 0x05ull && C.w[0] >= 0x0000000000000000ull) { | |
1762 | // set invalid flag | |
1763 | *pfpsf |= INVALID_EXCEPTION; | |
1764 | // return Integer Indefinite | |
1765 | res = 0x8000000000000000ull; | |
1766 | BID_RETURN (res); | |
1767 | } | |
1768 | // else cases that can be rounded to a 64-bit int fall through | |
1769 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
1770 | } // end else if n > 0 and q + exp = 19 |
1771 | } // end else if ((q + exp) == 19) | |
200359e8 L |
1772 | |
1773 | // n is not too large to be converted to int64: -2^63-1 < n < 2^63 | |
1774 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 1775 | if ((q + exp) <= 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
1776 | // set inexact flag |
1777 | *pfpsf |= INEXACT_EXCEPTION; | |
1778 | // return 0 | |
1779 | res = 0x0000000000000000ull; | |
1780 | BID_RETURN (res); | |
b2a00c89 | 1781 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
1782 | // -2^63-1 < x <= -1 or 1 <= x < 2^63 so x can be rounded |
1783 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
1784 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
1785 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
1786 | // chop off ind digits from the lower part of C1 |
1787 | // C1 fits in 64 bits | |
1788 | // calculate C* and f* | |
1789 | // C* is actually floor(C*) in this case | |
1790 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 1791 | // shiftright128[] and maskhigh128[] |
200359e8 | 1792 | // 1 <= x <= 15 |
b2a00c89 | 1793 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
1794 | // C* = C1 * 10^(-x) |
1795 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 1796 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 1797 | Cstar = P128.w[1]; |
b2a00c89 | 1798 | fstar.w[1] = P128.w[1] & maskhigh128[ind - 1]; |
200359e8 | 1799 | fstar.w[0] = P128.w[0]; |
b2a00c89 L |
1800 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
1801 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
1802 | // C* = floor(C*) (logical right shift; C has p decimal digits, |
1803 | // correct by Property 1) | |
1804 | // n = C* * 10^(e+x) | |
1805 | ||
b2a00c89 L |
1806 | // shift right C* by Ex-64 = shiftright128[ind] |
1807 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
1808 | Cstar = Cstar >> shift; |
1809 | // determine inexactness of the rounding of C* | |
1810 | // if (0 < f* < 10^(-x)) then | |
1811 | // the result is exact | |
1812 | // else // if (f* > T*) then | |
1813 | // the result is inexact | |
b2a00c89 L |
1814 | if (ind - 1 <= 2) { // fstar.w[1] is 0 |
1815 | if (fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
1816 | // ten2mk128trunc[ind -1].w[1] is identical to | |
1817 | // ten2mk128[ind -1].w[1] | |
200359e8 L |
1818 | // set the inexact flag |
1819 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 L |
1820 | } // else the result is exact |
1821 | } else { // if 3 <= ind - 1 <= 14 | |
1822 | if (fstar.w[1] || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
1823 | // ten2mk128trunc[ind -1].w[1] is identical to | |
1824 | // ten2mk128[ind -1].w[1] | |
200359e8 L |
1825 | // set the inexact flag |
1826 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 | 1827 | } // else the result is exact |
200359e8 L |
1828 | } |
1829 | if (x_sign) | |
1830 | res = -Cstar; | |
1831 | else | |
1832 | res = Cstar; | |
1833 | } else if (exp == 0) { | |
1834 | // 1 <= q <= 16 | |
1835 | // res = +/-C (exact) | |
1836 | if (x_sign) | |
1837 | res = -C1; | |
1838 | else | |
1839 | res = C1; | |
b2a00c89 | 1840 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
1841 | // (the upper limit of 20 on q + exp is due to the fact that |
1842 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
1843 | // res = +/-C * 10^exp (exact) | |
1844 | if (x_sign) | |
b2a00c89 | 1845 | res = -C1 * ten2k64[exp]; |
200359e8 | 1846 | else |
b2a00c89 | 1847 | res = C1 * ten2k64[exp]; |
200359e8 L |
1848 | } |
1849 | } | |
1850 | BID_RETURN (res); | |
1851 | } | |
1852 | ||
1853 | /***************************************************************************** | |
1854 | * BID64_to_int64_rninta | |
1855 | ****************************************************************************/ | |
1856 | ||
1857 | #if DECIMAL_CALL_BY_REFERENCE | |
1858 | void | |
b2a00c89 | 1859 | bid64_to_int64_rninta (SINT64 * pres, UINT64 * px |
200359e8 L |
1860 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
1861 | _EXC_INFO_PARAM) { | |
1862 | UINT64 x = *px; | |
1863 | #else | |
1864 | SINT64 | |
b2a00c89 | 1865 | bid64_to_int64_rninta (UINT64 x |
200359e8 L |
1866 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
1867 | _EXC_INFO_PARAM) { | |
1868 | #endif | |
1869 | SINT64 res; | |
1870 | UINT64 x_sign; | |
1871 | UINT64 x_exp; | |
b2a00c89 | 1872 | int exp; // unbiased exponent |
200359e8 L |
1873 | // Note: C1 represents x_significand (UINT64) |
1874 | BID_UI64DOUBLE tmp1; | |
1875 | unsigned int x_nr_bits; | |
1876 | int q, ind, shift; | |
1877 | UINT64 C1; | |
1878 | UINT128 C; | |
b2a00c89 | 1879 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
1880 | UINT128 P128; |
1881 | ||
1882 | // check for NaN or Infinity | |
1883 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
1884 | // set invalid flag | |
1885 | *pfpsf |= INVALID_EXCEPTION; | |
1886 | // return Integer Indefinite | |
1887 | res = 0x8000000000000000ull; | |
1888 | BID_RETURN (res); | |
1889 | } | |
1890 | // unpack x | |
b2a00c89 | 1891 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
1892 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
1893 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 1894 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 1895 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 1896 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
1897 | x_exp = 0; |
1898 | C1 = 0; | |
1899 | } | |
1900 | } else { | |
b2a00c89 | 1901 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
1902 | C1 = x & MASK_BINARY_SIG1; |
1903 | } | |
1904 | ||
1905 | // check for zeros (possibly from non-canonical values) | |
1906 | if (C1 == 0x0ull) { | |
1907 | // x is 0 | |
1908 | res = 0x00000000; | |
1909 | BID_RETURN (res); | |
1910 | } | |
1911 | // x is not special and is not zero | |
1912 | ||
1913 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
1914 | // determine first the nr. of bits in x | |
b2a00c89 | 1915 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 1916 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
1917 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
1918 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
1919 | x_nr_bits = |
1920 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
1921 | } else { // x < 2^32 |
1922 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1923 | x_nr_bits = |
1924 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1925 | } | |
b2a00c89 L |
1926 | } else { // if x < 2^53 |
1927 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
1928 | x_nr_bits = |
1929 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
1930 | } | |
b2a00c89 | 1931 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 1932 | if (q == 0) { |
b2a00c89 L |
1933 | q = nr_digits[x_nr_bits - 1].digits1; |
1934 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
1935 | q++; |
1936 | } | |
b2a00c89 | 1937 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 1938 | |
b2a00c89 | 1939 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
1940 | // set invalid flag |
1941 | *pfpsf |= INVALID_EXCEPTION; | |
1942 | // return Integer Indefinite | |
1943 | res = 0x8000000000000000ull; | |
1944 | BID_RETURN (res); | |
b2a00c89 | 1945 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
1946 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
1947 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
1948 | // the cases that do not fit are identified here; the ones that fit | |
1949 | // fall through and will be handled with other cases further, | |
1950 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 1951 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
1952 | // if n <= -2^63 - 1/2 then n is too large |
1953 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63+1/2 | |
1954 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64+1), 1<=q<=16 | |
1955 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 0x50000000000000005, 1<=q<=16 | |
1956 | // <=> C * 10^(20-q) >= 0x50000000000000005, 1<=q<=16 | |
1957 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1958 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1959 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x50000000000000005, has 20 |
1960 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] >= 0x05ull)) { | |
1961 | // set invalid flag | |
1962 | *pfpsf |= INVALID_EXCEPTION; | |
1963 | // return Integer Indefinite | |
1964 | res = 0x8000000000000000ull; | |
1965 | BID_RETURN (res); | |
1966 | } | |
1967 | // else cases that can be rounded to a 64-bit int fall through | |
1968 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 1969 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
1970 | // if n >= 2^63 - 1/2 then n is too large |
1971 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63-1/2 | |
1972 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64-1), 1<=q<=16 | |
1973 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 0x4fffffffffffffffb, 1<=q<=16 | |
1974 | // <=> if C * 10^(20-q) >= 0x4fffffffffffffffb, 1<=q<=16 | |
1975 | C.w[1] = 0x0000000000000004ull; | |
1976 | C.w[0] = 0xfffffffffffffffbull; | |
1977 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 1978 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
1979 | if (C.w[1] > 0x04ull || |
1980 | (C.w[1] == 0x04ull && C.w[0] >= 0xfffffffffffffffbull)) { | |
1981 | // set invalid flag | |
1982 | *pfpsf |= INVALID_EXCEPTION; | |
1983 | // return Integer Indefinite | |
1984 | res = 0x8000000000000000ull; | |
1985 | BID_RETURN (res); | |
1986 | } | |
1987 | // else cases that can be rounded to a 64-bit int fall through | |
1988 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
1989 | } // end else if n > 0 and q + exp = 19 |
1990 | } // end else if ((q + exp) == 19) | |
200359e8 L |
1991 | |
1992 | // n is not too large to be converted to int64: -2^63-1/2 < n < 2^63-1/2 | |
1993 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 1994 | if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
1995 | // return 0 |
1996 | res = 0x0000000000000000ull; | |
1997 | BID_RETURN (res); | |
b2a00c89 | 1998 | } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1) |
200359e8 L |
1999 | // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1) |
2000 | // res = 0 | |
2001 | // else | |
2002 | // res = +/-1 | |
b2a00c89 L |
2003 | ind = q - 1; // 0 <= ind <= 15 |
2004 | if (C1 < midpoint64[ind]) { | |
2005 | res = 0x0000000000000000ull; // return 0 | |
2006 | } else if (x_sign) { // n < 0 | |
2007 | res = 0xffffffffffffffffull; // return -1 | |
2008 | } else { // n > 0 | |
2009 | res = 0x0000000000000001ull; // return +1 | |
200359e8 | 2010 | } |
b2a00c89 | 2011 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
2012 | // -2^63-1/2 < x <= -1 or 1 <= x < 2^63-1/2 so x can be rounded |
2013 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
2014 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
2015 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
2016 | // chop off ind digits from the lower part of C1 |
2017 | // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 64 bits | |
b2a00c89 | 2018 | C1 = C1 + midpoint64[ind - 1]; |
200359e8 L |
2019 | // calculate C* and f* |
2020 | // C* is actually floor(C*) in this case | |
2021 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 2022 | // shiftright128[] and maskhigh128[] |
200359e8 | 2023 | // 1 <= x <= 15 |
b2a00c89 | 2024 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
2025 | // C* = (C1 + 1/2 * 10^x) * 10^(-x) |
2026 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 2027 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 2028 | Cstar = P128.w[1]; |
b2a00c89 L |
2029 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
2030 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
2031 | // if (0 < f* < 10^(-x)) then the result is a midpoint |
2032 | // if floor(C*) is even then C* = floor(C*) - logical right | |
2033 | // shift; C* has p decimal digits, correct by Prop. 1) | |
2034 | // else if floor(C*) is odd C* = floor(C*)-1 (logical right | |
2035 | // shift; C* has p decimal digits, correct by Pr. 1) | |
2036 | // else | |
2037 | // C* = floor(C*) (logical right shift; C has p decimal digits, | |
2038 | // correct by Property 1) | |
2039 | // n = C* * 10^(e+x) | |
2040 | ||
b2a00c89 L |
2041 | // shift right C* by Ex-64 = shiftright128[ind] |
2042 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
2043 | Cstar = Cstar >> shift; |
2044 | ||
2045 | // if the result was a midpoint it was rounded away from zero | |
2046 | if (x_sign) | |
2047 | res = -Cstar; | |
2048 | else | |
2049 | res = Cstar; | |
2050 | } else if (exp == 0) { | |
2051 | // 1 <= q <= 16 | |
2052 | // res = +/-C (exact) | |
2053 | if (x_sign) | |
2054 | res = -C1; | |
2055 | else | |
2056 | res = C1; | |
b2a00c89 | 2057 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
2058 | // (the upper limit of 20 on q + exp is due to the fact that |
2059 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
2060 | // res = +/-C * 10^exp (exact) | |
2061 | if (x_sign) | |
b2a00c89 | 2062 | res = -C1 * ten2k64[exp]; |
200359e8 | 2063 | else |
b2a00c89 | 2064 | res = C1 * ten2k64[exp]; |
200359e8 L |
2065 | } |
2066 | } | |
2067 | BID_RETURN (res); | |
2068 | } | |
2069 | ||
2070 | /***************************************************************************** | |
2071 | * BID64_to_int64_xrninta | |
2072 | ****************************************************************************/ | |
2073 | ||
2074 | #if DECIMAL_CALL_BY_REFERENCE | |
2075 | void | |
b2a00c89 | 2076 | bid64_to_int64_xrninta (SINT64 * pres, UINT64 * px |
200359e8 L |
2077 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
2078 | _EXC_INFO_PARAM) { | |
2079 | UINT64 x = *px; | |
2080 | #else | |
2081 | SINT64 | |
b2a00c89 | 2082 | bid64_to_int64_xrninta (UINT64 x |
200359e8 L |
2083 | _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
2084 | _EXC_INFO_PARAM) { | |
2085 | #endif | |
2086 | SINT64 res; | |
2087 | UINT64 x_sign; | |
2088 | UINT64 x_exp; | |
b2a00c89 | 2089 | int exp; // unbiased exponent |
200359e8 L |
2090 | // Note: C1 represents x_significand (UINT64) |
2091 | UINT64 tmp64; | |
2092 | BID_UI64DOUBLE tmp1; | |
2093 | unsigned int x_nr_bits; | |
2094 | int q, ind, shift; | |
2095 | UINT64 C1; | |
2096 | UINT128 C; | |
b2a00c89 | 2097 | UINT64 Cstar; // C* represents up to 16 decimal digits ~ 54 bits |
200359e8 L |
2098 | UINT128 fstar; |
2099 | UINT128 P128; | |
2100 | ||
2101 | // check for NaN or Infinity | |
2102 | if ((x & MASK_NAN) == MASK_NAN || (x & MASK_INF) == MASK_INF) { | |
2103 | // set invalid flag | |
2104 | *pfpsf |= INVALID_EXCEPTION; | |
2105 | // return Integer Indefinite | |
2106 | res = 0x8000000000000000ull; | |
2107 | BID_RETURN (res); | |
2108 | } | |
2109 | // unpack x | |
b2a00c89 | 2110 | x_sign = x & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
200359e8 L |
2111 | // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] => |
2112 | if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) { | |
b2a00c89 | 2113 | x_exp = (x & MASK_BINARY_EXPONENT2) >> 51; // biased |
200359e8 | 2114 | C1 = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2; |
b2a00c89 | 2115 | if (C1 > 9999999999999999ull) { // non-canonical |
200359e8 L |
2116 | x_exp = 0; |
2117 | C1 = 0; | |
2118 | } | |
2119 | } else { | |
b2a00c89 | 2120 | x_exp = (x & MASK_BINARY_EXPONENT1) >> 53; // biased |
200359e8 L |
2121 | C1 = x & MASK_BINARY_SIG1; |
2122 | } | |
2123 | ||
2124 | // check for zeros (possibly from non-canonical values) | |
2125 | if (C1 == 0x0ull) { | |
2126 | // x is 0 | |
2127 | res = 0x00000000; | |
2128 | BID_RETURN (res); | |
2129 | } | |
2130 | // x is not special and is not zero | |
2131 | ||
2132 | // q = nr. of decimal digits in x (1 <= q <= 54) | |
2133 | // determine first the nr. of bits in x | |
b2a00c89 | 2134 | if (C1 >= 0x0020000000000000ull) { // x >= 2^53 |
200359e8 | 2135 | // split the 64-bit value in two 32-bit halves to avoid rounding errors |
b2a00c89 L |
2136 | if (C1 >= 0x0000000100000000ull) { // x >= 2^32 |
2137 | tmp1.d = (double) (C1 >> 32); // exact conversion | |
200359e8 L |
2138 | x_nr_bits = |
2139 | 33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
b2a00c89 L |
2140 | } else { // x < 2^32 |
2141 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
2142 | x_nr_bits = |
2143 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
2144 | } | |
b2a00c89 L |
2145 | } else { // if x < 2^53 |
2146 | tmp1.d = (double) C1; // exact conversion | |
200359e8 L |
2147 | x_nr_bits = |
2148 | 1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); | |
2149 | } | |
b2a00c89 | 2150 | q = nr_digits[x_nr_bits - 1].digits; |
200359e8 | 2151 | if (q == 0) { |
b2a00c89 L |
2152 | q = nr_digits[x_nr_bits - 1].digits1; |
2153 | if (C1 >= nr_digits[x_nr_bits - 1].threshold_lo) | |
200359e8 L |
2154 | q++; |
2155 | } | |
b2a00c89 | 2156 | exp = x_exp - 398; // unbiased exponent |
200359e8 | 2157 | |
b2a00c89 | 2158 | if ((q + exp) > 19) { // x >= 10^19 ~= 2^63.11... (cannot fit in SINT64) |
200359e8 L |
2159 | // set invalid flag |
2160 | *pfpsf |= INVALID_EXCEPTION; | |
2161 | // return Integer Indefinite | |
2162 | res = 0x8000000000000000ull; | |
2163 | BID_RETURN (res); | |
b2a00c89 | 2164 | } else if ((q + exp) == 19) { // x = c(0)c(1)...c(18).c(19)...c(q-1) |
200359e8 L |
2165 | // in this case 2^63.11... ~= 10^19 <= x < 10^20 ~= 2^66.43... |
2166 | // so x rounded to an integer may or may not fit in a signed 64-bit int | |
2167 | // the cases that do not fit are identified here; the ones that fit | |
2168 | // fall through and will be handled with other cases further, | |
2169 | // under '1 <= q + exp <= 19' | |
b2a00c89 | 2170 | if (x_sign) { // if n < 0 and q + exp = 19 |
200359e8 L |
2171 | // if n <= -2^63 - 1/2 then n is too large |
2172 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63+1/2 | |
2173 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64+1), 1<=q<=16 | |
2174 | // <=> 0.c(0)c(1)...c(q-1) * 10^20 >= 0x50000000000000005, 1<=q<=16 | |
2175 | // <=> C * 10^(20-q) >= 0x50000000000000005, 1<=q<=16 | |
2176 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 2177 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
2178 | // Note: C1 * 10^(11-q) has 19 or 20 digits; 0x50000000000000005, has 20 |
2179 | if (C.w[1] > 0x05ull || (C.w[1] == 0x05ull && C.w[0] >= 0x05ull)) { | |
2180 | // set invalid flag | |
2181 | *pfpsf |= INVALID_EXCEPTION; | |
2182 | // return Integer Indefinite | |
2183 | res = 0x8000000000000000ull; | |
2184 | BID_RETURN (res); | |
2185 | } | |
2186 | // else cases that can be rounded to a 64-bit int fall through | |
2187 | // to '1 <= q + exp <= 19' | |
b2a00c89 | 2188 | } else { // if n > 0 and q + exp = 19 |
200359e8 L |
2189 | // if n >= 2^63 - 1/2 then n is too large |
2190 | // too large if c(0)c(1)...c(18).c(19)...c(q-1) >= 2^63-1/2 | |
2191 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 5*(2^64-1), 1<=q<=16 | |
2192 | // <=> if 0.c(0)c(1)...c(q-1) * 10^20 >= 0x4fffffffffffffffb, 1<=q<=16 | |
2193 | // <=> if C * 10^(20-q) >= 0x4fffffffffffffffb, 1<=q<=16 | |
2194 | C.w[1] = 0x0000000000000004ull; | |
2195 | C.w[0] = 0xfffffffffffffffbull; | |
2196 | // 1 <= q <= 16 => 4 <= 20-q <= 19 => 10^(20-q) is 64-bit, and so is C1 | |
b2a00c89 | 2197 | __mul_64x64_to_128MACH (C, C1, ten2k64[20 - q]); |
200359e8 L |
2198 | if (C.w[1] > 0x04ull || |
2199 | (C.w[1] == 0x04ull && C.w[0] >= 0xfffffffffffffffbull)) { | |
2200 | // set invalid flag | |
2201 | *pfpsf |= INVALID_EXCEPTION; | |
2202 | // return Integer Indefinite | |
2203 | res = 0x8000000000000000ull; | |
2204 | BID_RETURN (res); | |
2205 | } | |
2206 | // else cases that can be rounded to a 64-bit int fall through | |
2207 | // to '1 <= q + exp <= 19' | |
b2a00c89 L |
2208 | } // end else if n > 0 and q + exp = 19 |
2209 | } // end else if ((q + exp) == 19) | |
200359e8 L |
2210 | |
2211 | // n is not too large to be converted to int64: -2^63-1/2 < n < 2^63-1/2 | |
2212 | // Note: some of the cases tested for above fall through to this point | |
b2a00c89 | 2213 | if ((q + exp) < 0) { // n = +/-0.0...c(0)c(1)...c(q-1) |
200359e8 L |
2214 | // set inexact flag |
2215 | *pfpsf |= INEXACT_EXCEPTION; | |
2216 | // return 0 | |
2217 | res = 0x0000000000000000ull; | |
2218 | BID_RETURN (res); | |
b2a00c89 | 2219 | } else if ((q + exp) == 0) { // n = +/-0.c(0)c(1)...c(q-1) |
200359e8 L |
2220 | // if 0.c(0)c(1)...c(q-1) <= 0.5 <=> c(0)c(1)...c(q-1) <= 5 * 10^(q-1) |
2221 | // res = 0 | |
2222 | // else | |
2223 | // res = +/-1 | |
b2a00c89 L |
2224 | ind = q - 1; // 0 <= ind <= 15 |
2225 | if (C1 < midpoint64[ind]) { | |
2226 | res = 0x0000000000000000ull; // return 0 | |
2227 | } else if (x_sign) { // n < 0 | |
2228 | res = 0xffffffffffffffffull; // return -1 | |
2229 | } else { // n > 0 | |
2230 | res = 0x0000000000000001ull; // return +1 | |
200359e8 L |
2231 | } |
2232 | // set inexact flag | |
2233 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 | 2234 | } else { // if (1 <= q + exp <= 19, 1 <= q <= 16, -15 <= exp <= 18) |
200359e8 L |
2235 | // -2^63-1/2 < x <= -1 or 1 <= x < 2^63-1/2 so x can be rounded |
2236 | // to nearest to a 64-bit signed integer | |
b2a00c89 L |
2237 | if (exp < 0) { // 2 <= q <= 16, -15 <= exp <= -1, 1 <= q + exp <= 19 |
2238 | ind = -exp; // 1 <= ind <= 15; ind is a synonym for 'x' | |
200359e8 L |
2239 | // chop off ind digits from the lower part of C1 |
2240 | // C1 = C1 + 1/2 * 10^ind where the result C1 fits in 64 bits | |
b2a00c89 | 2241 | C1 = C1 + midpoint64[ind - 1]; |
200359e8 L |
2242 | // calculate C* and f* |
2243 | // C* is actually floor(C*) in this case | |
2244 | // C* and f* need shifting and masking, as shown by | |
b2a00c89 | 2245 | // shiftright128[] and maskhigh128[] |
200359e8 | 2246 | // 1 <= x <= 15 |
b2a00c89 | 2247 | // kx = 10^(-x) = ten2mk64[ind - 1] |
200359e8 L |
2248 | // C* = (C1 + 1/2 * 10^x) * 10^(-x) |
2249 | // the approximation of 10^(-x) was rounded up to 54 bits | |
b2a00c89 | 2250 | __mul_64x64_to_128MACH (P128, C1, ten2mk64[ind - 1]); |
200359e8 | 2251 | Cstar = P128.w[1]; |
b2a00c89 | 2252 | fstar.w[1] = P128.w[1] & maskhigh128[ind - 1]; |
200359e8 | 2253 | fstar.w[0] = P128.w[0]; |
b2a00c89 L |
2254 | // the top Ex bits of 10^(-x) are T* = ten2mk128trunc[ind].w[0], e.g. |
2255 | // if x=1, T*=ten2mk128trunc[0].w[0]=0x1999999999999999 | |
200359e8 L |
2256 | // if (0 < f* < 10^(-x)) then the result is a midpoint |
2257 | // if floor(C*) is even then C* = floor(C*) - logical right | |
2258 | // shift; C* has p decimal digits, correct by Prop. 1) | |
2259 | // else if floor(C*) is odd C* = floor(C*)-1 (logical right | |
2260 | // shift; C* has p decimal digits, correct by Pr. 1) | |
2261 | // else | |
2262 | // C* = floor(C*) (logical right shift; C has p decimal digits, | |
2263 | // correct by Property 1) | |
2264 | // n = C* * 10^(e+x) | |
2265 | ||
b2a00c89 L |
2266 | // shift right C* by Ex-64 = shiftright128[ind] |
2267 | shift = shiftright128[ind - 1]; // 0 <= shift <= 39 | |
200359e8 L |
2268 | Cstar = Cstar >> shift; |
2269 | // determine inexactness of the rounding of C* | |
2270 | // if (0 < f* - 1/2 < 10^(-x)) then | |
2271 | // the result is exact | |
2272 | // else // if (f* - 1/2 > T*) then | |
2273 | // the result is inexact | |
2274 | if (ind - 1 <= 2) { | |
2275 | if (fstar.w[0] > 0x8000000000000000ull) { | |
2276 | // f* > 1/2 and the result may be exact | |
b2a00c89 L |
2277 | tmp64 = fstar.w[0] - 0x8000000000000000ull; // f* - 1/2 |
2278 | if ((tmp64 > ten2mk128trunc[ind - 1].w[1])) { | |
2279 | // ten2mk128trunc[ind -1].w[1] is identical to | |
2280 | // ten2mk128[ind -1].w[1] | |
200359e8 L |
2281 | // set the inexact flag |
2282 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 L |
2283 | } // else the result is exact |
2284 | } else { // the result is inexact; f2* <= 1/2 | |
200359e8 L |
2285 | // set the inexact flag |
2286 | *pfpsf |= INEXACT_EXCEPTION; | |
2287 | } | |
b2a00c89 L |
2288 | } else { // if 3 <= ind - 1 <= 14 |
2289 | if (fstar.w[1] > onehalf128[ind - 1] || | |
2290 | (fstar.w[1] == onehalf128[ind - 1] && fstar.w[0])) { | |
200359e8 L |
2291 | // f2* > 1/2 and the result may be exact |
2292 | // Calculate f2* - 1/2 | |
b2a00c89 L |
2293 | tmp64 = fstar.w[1] - onehalf128[ind - 1]; |
2294 | if (tmp64 || fstar.w[0] > ten2mk128trunc[ind - 1].w[1]) { | |
2295 | // ten2mk128trunc[ind -1].w[1] is identical to | |
2296 | // ten2mk128[ind -1].w[1] | |
200359e8 L |
2297 | // set the inexact flag |
2298 | *pfpsf |= INEXACT_EXCEPTION; | |
b2a00c89 L |
2299 | } // else the result is exact |
2300 | } else { // the result is inexact; f2* <= 1/2 | |
200359e8 L |
2301 | // set the inexact flag |
2302 | *pfpsf |= INEXACT_EXCEPTION; | |
2303 | } | |
2304 | } | |
2305 | ||
2306 | // if the result was a midpoint it was rounded away from zero | |
2307 | if (x_sign) | |
2308 | res = -Cstar; | |
2309 | else | |
2310 | res = Cstar; | |
2311 | } else if (exp == 0) { | |
2312 | // 1 <= q <= 16 | |
2313 | // res = +/-C (exact) | |
2314 | if (x_sign) | |
2315 | res = -C1; | |
2316 | else | |
2317 | res = C1; | |
b2a00c89 | 2318 | } else { // if (exp > 0) => 1 <= exp <= 18, 1 <= q <= 16, 2 <= q + exp <= 20 |
200359e8 L |
2319 | // (the upper limit of 20 on q + exp is due to the fact that |
2320 | // +/-C * 10^exp is guaranteed to fit in 64 bits) | |
2321 | // res = +/-C * 10^exp (exact) | |
2322 | if (x_sign) | |
b2a00c89 | 2323 | res = -C1 * ten2k64[exp]; |
200359e8 | 2324 | else |
b2a00c89 | 2325 | res = C1 * ten2k64[exp]; |
200359e8 L |
2326 | } |
2327 | } | |
2328 | BID_RETURN (res); | |
2329 | } |