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