| 6515 |
serge |
1 |
/* Copyright (C) 2007-2015 Free Software Foundation, Inc. |
|
|
2 |
|
|
|
3 |
This file is part of GCC. |
|
|
4 |
|
|
|
5 |
GCC is free software; you can redistribute it and/or modify it under |
|
|
6 |
the terms of the GNU General Public License as published by the Free |
|
|
7 |
Software Foundation; either version 3, or (at your option) any later |
|
|
8 |
version. |
|
|
9 |
|
|
|
10 |
GCC is distributed in the hope that it will be useful, but WITHOUT ANY |
|
|
11 |
WARRANTY; without even the implied warranty of MERCHANTABILITY or |
|
|
12 |
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
|
|
13 |
for more details. |
|
|
14 |
|
|
|
15 |
Under Section 7 of GPL version 3, you are granted additional |
|
|
16 |
permissions described in the GCC Runtime Library Exception, version |
|
|
17 |
3.1, as published by the Free Software Foundation. |
|
|
18 |
|
|
|
19 |
You should have received a copy of the GNU General Public License and |
|
|
20 |
a copy of the GCC Runtime Library Exception along with this program; |
|
|
21 |
see the files COPYING3 and COPYING.RUNTIME respectively. If not, see |
|
|
22 |
. */ |
|
|
23 |
|
|
|
24 |
#define BID_128RES |
|
|
25 |
#include "bid_internal.h" |
|
|
26 |
|
|
|
27 |
/***************************************************************************** |
|
|
28 |
* BID128 nextup |
|
|
29 |
****************************************************************************/ |
|
|
30 |
|
|
|
31 |
#if DECIMAL_CALL_BY_REFERENCE |
|
|
32 |
void |
|
|
33 |
bid128_nextup (UINT128 * pres, |
|
|
34 |
UINT128 * |
|
|
35 |
px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) { |
|
|
36 |
UINT128 x = *px; |
|
|
37 |
#else |
|
|
38 |
UINT128 |
|
|
39 |
bid128_nextup (UINT128 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
|
|
40 |
_EXC_INFO_PARAM) { |
|
|
41 |
#endif |
|
|
42 |
|
|
|
43 |
UINT128 res; |
|
|
44 |
UINT64 x_sign; |
|
|
45 |
UINT64 x_exp; |
|
|
46 |
int exp; |
|
|
47 |
BID_UI64DOUBLE tmp1; |
|
|
48 |
int x_nr_bits; |
|
|
49 |
int q1, ind; |
|
|
50 |
UINT128 C1; // C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (UINT64) |
|
|
51 |
|
|
|
52 |
BID_SWAP128 (x); |
|
|
53 |
// unpack the argument |
|
|
54 |
x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
|
|
55 |
C1.w[1] = x.w[1] & MASK_COEFF; |
|
|
56 |
C1.w[0] = x.w[0]; |
|
|
57 |
|
|
|
58 |
// check for NaN or Infinity |
|
|
59 |
if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) { |
|
|
60 |
// x is special |
|
|
61 |
if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN |
|
|
62 |
// if x = NaN, then res = Q (x) |
|
|
63 |
// check first for non-canonical NaN payload |
|
|
64 |
if (((x.w[1] & 0x00003fffffffffffull) > 0x0000314dc6448d93ull) || |
|
|
65 |
(((x.w[1] & 0x00003fffffffffffull) == 0x0000314dc6448d93ull) |
|
|
66 |
&& (x.w[0] > 0x38c15b09ffffffffull))) { |
|
|
67 |
x.w[1] = x.w[1] & 0xffffc00000000000ull; |
|
|
68 |
x.w[0] = 0x0ull; |
|
|
69 |
} |
|
|
70 |
if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN |
|
|
71 |
// set invalid flag |
|
|
72 |
*pfpsf |= INVALID_EXCEPTION; |
|
|
73 |
// return quiet (x) |
|
|
74 |
res.w[1] = x.w[1] & 0xfc003fffffffffffull; // clear out also G[6]-G[16] |
|
|
75 |
res.w[0] = x.w[0]; |
|
|
76 |
} else { // x is QNaN |
|
|
77 |
// return x |
|
|
78 |
res.w[1] = x.w[1] & 0xfc003fffffffffffull; // clear out G[6]-G[16] |
|
|
79 |
res.w[0] = x.w[0]; |
|
|
80 |
} |
|
|
81 |
} else { // x is not NaN, so it must be infinity |
|
|
82 |
if (!x_sign) { // x is +inf |
|
|
83 |
res.w[1] = 0x7800000000000000ull; // +inf |
|
|
84 |
res.w[0] = 0x0000000000000000ull; |
|
|
85 |
} else { // x is -inf |
|
|
86 |
res.w[1] = 0xdfffed09bead87c0ull; // -MAXFP = -999...99 * 10^emax |
|
|
87 |
res.w[0] = 0x378d8e63ffffffffull; |
|
|
88 |
} |
|
|
89 |
} |
|
|
90 |
BID_RETURN (res); |
|
|
91 |
} |
|
|
92 |
// check for non-canonical values (treated as zero) |
|
|
93 |
if ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull) { // G0_G1=11 |
|
|
94 |
// non-canonical |
|
|
95 |
x_exp = (x.w[1] << 2) & MASK_EXP; // biased and shifted left 49 bits |
|
|
96 |
C1.w[1] = 0; // significand high |
|
|
97 |
C1.w[0] = 0; // significand low |
|
|
98 |
} else { // G0_G1 != 11 |
|
|
99 |
x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bits |
|
|
100 |
if (C1.w[1] > 0x0001ed09bead87c0ull || |
|
|
101 |
(C1.w[1] == 0x0001ed09bead87c0ull |
|
|
102 |
&& C1.w[0] > 0x378d8e63ffffffffull)) { |
|
|
103 |
// x is non-canonical if coefficient is larger than 10^34 -1 |
|
|
104 |
C1.w[1] = 0; |
|
|
105 |
C1.w[0] = 0; |
|
|
106 |
} else { // canonical |
|
|
107 |
; |
|
|
108 |
} |
|
|
109 |
} |
|
|
110 |
|
|
|
111 |
if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) { |
|
|
112 |
// x is +/-0 |
|
|
113 |
res.w[1] = 0x0000000000000000ull; // +1 * 10^emin |
|
|
114 |
res.w[0] = 0x0000000000000001ull; |
|
|
115 |
} else { // x is not special and is not zero |
|
|
116 |
if (x.w[1] == 0x5fffed09bead87c0ull |
|
|
117 |
&& x.w[0] == 0x378d8e63ffffffffull) { |
|
|
118 |
// x = +MAXFP = 999...99 * 10^emax |
|
|
119 |
res.w[1] = 0x7800000000000000ull; // +inf |
|
|
120 |
res.w[0] = 0x0000000000000000ull; |
|
|
121 |
} else if (x.w[1] == 0x8000000000000000ull |
|
|
122 |
&& x.w[0] == 0x0000000000000001ull) { |
|
|
123 |
// x = -MINFP = 1...99 * 10^emin |
|
|
124 |
res.w[1] = 0x8000000000000000ull; // -0 |
|
|
125 |
res.w[0] = 0x0000000000000000ull; |
|
|
126 |
} else { // -MAXFP <= x <= -MINFP - 1 ulp OR MINFP <= x <= MAXFP - 1 ulp |
|
|
127 |
// can add/subtract 1 ulp to the significand |
|
|
128 |
|
|
|
129 |
// Note: we could check here if x >= 10^34 to speed up the case q1 = 34 |
|
|
130 |
// q1 = nr. of decimal digits in x |
|
|
131 |
// determine first the nr. of bits in x |
|
|
132 |
if (C1.w[1] == 0) { |
|
|
133 |
if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53 |
|
|
134 |
// split the 64-bit value in two 32-bit halves to avoid rnd errors |
|
|
135 |
if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32 |
|
|
136 |
tmp1.d = (double) (C1.w[0] >> 32); // exact conversion |
|
|
137 |
x_nr_bits = |
|
|
138 |
33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - |
|
|
139 |
0x3ff); |
|
|
140 |
} else { // x < 2^32 |
|
|
141 |
tmp1.d = (double) (C1.w[0]); // exact conversion |
|
|
142 |
x_nr_bits = |
|
|
143 |
1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - |
|
|
144 |
0x3ff); |
|
|
145 |
} |
|
|
146 |
} else { // if x < 2^53 |
|
|
147 |
tmp1.d = (double) C1.w[0]; // exact conversion |
|
|
148 |
x_nr_bits = |
|
|
149 |
1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); |
|
|
150 |
} |
|
|
151 |
} else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1]) |
|
|
152 |
tmp1.d = (double) C1.w[1]; // exact conversion |
|
|
153 |
x_nr_bits = |
|
|
154 |
65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); |
|
|
155 |
} |
|
|
156 |
q1 = nr_digits[x_nr_bits - 1].digits; |
|
|
157 |
if (q1 == 0) { |
|
|
158 |
q1 = nr_digits[x_nr_bits - 1].digits1; |
|
|
159 |
if (C1.w[1] > nr_digits[x_nr_bits - 1].threshold_hi |
|
|
160 |
|| (C1.w[1] == nr_digits[x_nr_bits - 1].threshold_hi |
|
|
161 |
&& C1.w[0] >= nr_digits[x_nr_bits - 1].threshold_lo)) |
|
|
162 |
q1++; |
|
|
163 |
} |
|
|
164 |
// if q1 < P34 then pad the significand with zeros |
|
|
165 |
if (q1 < P34) { |
|
|
166 |
exp = (x_exp >> 49) - 6176; |
|
|
167 |
if (exp + 6176 > P34 - q1) { |
|
|
168 |
ind = P34 - q1; // 1 <= ind <= P34 - 1 |
|
|
169 |
// pad with P34 - q1 zeros, until exponent = emin |
|
|
170 |
// C1 = C1 * 10^ind |
|
|
171 |
if (q1 <= 19) { // 64-bit C1 |
|
|
172 |
if (ind <= 19) { // 64-bit 10^ind and 64-bit C1 |
|
|
173 |
__mul_64x64_to_128MACH (C1, C1.w[0], ten2k64[ind]); |
|
|
174 |
} else { // 128-bit 10^ind and 64-bit C1 |
|
|
175 |
__mul_128x64_to_128 (C1, C1.w[0], ten2k128[ind - 20]); |
|
|
176 |
} |
|
|
177 |
} else { // C1 is (most likely) 128-bit |
|
|
178 |
if (ind <= 14) { // 64-bit 10^ind and 128-bit C1 (most likely) |
|
|
179 |
__mul_128x64_to_128 (C1, ten2k64[ind], C1); |
|
|
180 |
} else if (ind <= 19) { // 64-bit 10^ind and 64-bit C1 (q1 <= 19) |
|
|
181 |
__mul_64x64_to_128MACH (C1, C1.w[0], ten2k64[ind]); |
|
|
182 |
} else { // 128-bit 10^ind and 64-bit C1 (C1 must be 64-bit) |
|
|
183 |
__mul_128x64_to_128 (C1, C1.w[0], ten2k128[ind - 20]); |
|
|
184 |
} |
|
|
185 |
} |
|
|
186 |
x_exp = x_exp - ((UINT64) ind << 49); |
|
|
187 |
} else { // pad with zeros until the exponent reaches emin |
|
|
188 |
ind = exp + 6176; |
|
|
189 |
// C1 = C1 * 10^ind |
|
|
190 |
if (ind <= 19) { // 1 <= P34 - q1 <= 19 <=> 15 <= q1 <= 33 |
|
|
191 |
if (q1 <= 19) { // 64-bit C1, 64-bit 10^ind |
|
|
192 |
__mul_64x64_to_128MACH (C1, C1.w[0], ten2k64[ind]); |
|
|
193 |
} else { // 20 <= q1 <= 33 => 128-bit C1, 64-bit 10^ind |
|
|
194 |
__mul_128x64_to_128 (C1, ten2k64[ind], C1); |
|
|
195 |
} |
|
|
196 |
} else { // if 20 <= P34 - q1 <= 33 <=> 1 <= q1 <= 14 => |
|
|
197 |
// 64-bit C1, 128-bit 10^ind |
|
|
198 |
__mul_128x64_to_128 (C1, C1.w[0], ten2k128[ind - 20]); |
|
|
199 |
} |
|
|
200 |
x_exp = EXP_MIN; |
|
|
201 |
} |
|
|
202 |
} |
|
|
203 |
if (!x_sign) { // x > 0 |
|
|
204 |
// add 1 ulp (add 1 to the significand) |
|
|
205 |
C1.w[0]++; |
|
|
206 |
if (C1.w[0] == 0) |
|
|
207 |
C1.w[1]++; |
|
|
208 |
if (C1.w[1] == 0x0001ed09bead87c0ull && C1.w[0] == 0x378d8e6400000000ull) { // if C1 = 10^34 |
|
|
209 |
C1.w[1] = 0x0000314dc6448d93ull; // C1 = 10^33 |
|
|
210 |
C1.w[0] = 0x38c15b0a00000000ull; |
|
|
211 |
x_exp = x_exp + EXP_P1; |
|
|
212 |
} |
|
|
213 |
} else { // x < 0 |
|
|
214 |
// subtract 1 ulp (subtract 1 from the significand) |
|
|
215 |
C1.w[0]--; |
|
|
216 |
if (C1.w[0] == 0xffffffffffffffffull) |
|
|
217 |
C1.w[1]--; |
|
|
218 |
if (x_exp != 0 && C1.w[1] == 0x0000314dc6448d93ull && C1.w[0] == 0x38c15b09ffffffffull) { // if C1 = 10^33 - 1 |
|
|
219 |
C1.w[1] = 0x0001ed09bead87c0ull; // C1 = 10^34 - 1 |
|
|
220 |
C1.w[0] = 0x378d8e63ffffffffull; |
|
|
221 |
x_exp = x_exp - EXP_P1; |
|
|
222 |
} |
|
|
223 |
} |
|
|
224 |
// assemble the result |
|
|
225 |
res.w[1] = x_sign | x_exp | C1.w[1]; |
|
|
226 |
res.w[0] = C1.w[0]; |
|
|
227 |
} // end -MAXFP <= x <= -MINFP - 1 ulp OR MINFP <= x <= MAXFP - 1 ulp |
|
|
228 |
} // end x is not special and is not zero |
|
|
229 |
BID_RETURN (res); |
|
|
230 |
} |
|
|
231 |
|
|
|
232 |
/***************************************************************************** |
|
|
233 |
* BID128 nextdown |
|
|
234 |
****************************************************************************/ |
|
|
235 |
|
|
|
236 |
#if DECIMAL_CALL_BY_REFERENCE |
|
|
237 |
void |
|
|
238 |
bid128_nextdown (UINT128 * pres, |
|
|
239 |
UINT128 * |
|
|
240 |
px _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) { |
|
|
241 |
UINT128 x = *px; |
|
|
242 |
#else |
|
|
243 |
UINT128 |
|
|
244 |
bid128_nextdown (UINT128 x _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
|
|
245 |
_EXC_INFO_PARAM) { |
|
|
246 |
#endif |
|
|
247 |
|
|
|
248 |
UINT128 res; |
|
|
249 |
UINT64 x_sign; |
|
|
250 |
UINT64 x_exp; |
|
|
251 |
int exp; |
|
|
252 |
BID_UI64DOUBLE tmp1; |
|
|
253 |
int x_nr_bits; |
|
|
254 |
int q1, ind; |
|
|
255 |
UINT128 C1; // C1.w[1], C1.w[0] represent x_signif_hi, x_signif_lo (UINT64) |
|
|
256 |
|
|
|
257 |
BID_SWAP128 (x); |
|
|
258 |
// unpack the argument |
|
|
259 |
x_sign = x.w[1] & MASK_SIGN; // 0 for positive, MASK_SIGN for negative |
|
|
260 |
C1.w[1] = x.w[1] & MASK_COEFF; |
|
|
261 |
C1.w[0] = x.w[0]; |
|
|
262 |
|
|
|
263 |
// check for NaN or Infinity |
|
|
264 |
if ((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) { |
|
|
265 |
// x is special |
|
|
266 |
if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN |
|
|
267 |
// if x = NaN, then res = Q (x) |
|
|
268 |
// check first for non-canonical NaN payload |
|
|
269 |
if (((x.w[1] & 0x00003fffffffffffull) > 0x0000314dc6448d93ull) || |
|
|
270 |
(((x.w[1] & 0x00003fffffffffffull) == 0x0000314dc6448d93ull) |
|
|
271 |
&& (x.w[0] > 0x38c15b09ffffffffull))) { |
|
|
272 |
x.w[1] = x.w[1] & 0xffffc00000000000ull; |
|
|
273 |
x.w[0] = 0x0ull; |
|
|
274 |
} |
|
|
275 |
if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN |
|
|
276 |
// set invalid flag |
|
|
277 |
*pfpsf |= INVALID_EXCEPTION; |
|
|
278 |
// return quiet (x) |
|
|
279 |
res.w[1] = x.w[1] & 0xfc003fffffffffffull; // clear out also G[6]-G[16] |
|
|
280 |
res.w[0] = x.w[0]; |
|
|
281 |
} else { // x is QNaN |
|
|
282 |
// return x |
|
|
283 |
res.w[1] = x.w[1] & 0xfc003fffffffffffull; // clear out G[6]-G[16] |
|
|
284 |
res.w[0] = x.w[0]; |
|
|
285 |
} |
|
|
286 |
} else { // x is not NaN, so it must be infinity |
|
|
287 |
if (!x_sign) { // x is +inf |
|
|
288 |
res.w[1] = 0x5fffed09bead87c0ull; // +MAXFP = +999...99 * 10^emax |
|
|
289 |
res.w[0] = 0x378d8e63ffffffffull; |
|
|
290 |
} else { // x is -inf |
|
|
291 |
res.w[1] = 0xf800000000000000ull; // -inf |
|
|
292 |
res.w[0] = 0x0000000000000000ull; |
|
|
293 |
} |
|
|
294 |
} |
|
|
295 |
BID_RETURN (res); |
|
|
296 |
} |
|
|
297 |
// check for non-canonical values (treated as zero) |
|
|
298 |
if ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull) { // G0_G1=11 |
|
|
299 |
// non-canonical |
|
|
300 |
x_exp = (x.w[1] << 2) & MASK_EXP; // biased and shifted left 49 bits |
|
|
301 |
C1.w[1] = 0; // significand high |
|
|
302 |
C1.w[0] = 0; // significand low |
|
|
303 |
} else { // G0_G1 != 11 |
|
|
304 |
x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bits |
|
|
305 |
if (C1.w[1] > 0x0001ed09bead87c0ull || |
|
|
306 |
(C1.w[1] == 0x0001ed09bead87c0ull |
|
|
307 |
&& C1.w[0] > 0x378d8e63ffffffffull)) { |
|
|
308 |
// x is non-canonical if coefficient is larger than 10^34 -1 |
|
|
309 |
C1.w[1] = 0; |
|
|
310 |
C1.w[0] = 0; |
|
|
311 |
} else { // canonical |
|
|
312 |
; |
|
|
313 |
} |
|
|
314 |
} |
|
|
315 |
|
|
|
316 |
if ((C1.w[1] == 0x0ull) && (C1.w[0] == 0x0ull)) { |
|
|
317 |
// x is +/-0 |
|
|
318 |
res.w[1] = 0x8000000000000000ull; // -1 * 10^emin |
|
|
319 |
res.w[0] = 0x0000000000000001ull; |
|
|
320 |
} else { // x is not special and is not zero |
|
|
321 |
if (x.w[1] == 0xdfffed09bead87c0ull |
|
|
322 |
&& x.w[0] == 0x378d8e63ffffffffull) { |
|
|
323 |
// x = -MAXFP = -999...99 * 10^emax |
|
|
324 |
res.w[1] = 0xf800000000000000ull; // -inf |
|
|
325 |
res.w[0] = 0x0000000000000000ull; |
|
|
326 |
} else if (x.w[1] == 0x0ull && x.w[0] == 0x0000000000000001ull) { // +MINFP |
|
|
327 |
res.w[1] = 0x0000000000000000ull; // +0 |
|
|
328 |
res.w[0] = 0x0000000000000000ull; |
|
|
329 |
} else { // -MAXFP <= x <= -MINFP - 1 ulp OR MINFP <= x <= MAXFP - 1 ulp |
|
|
330 |
// can add/subtract 1 ulp to the significand |
|
|
331 |
|
|
|
332 |
// Note: we could check here if x >= 10^34 to speed up the case q1 = 34 |
|
|
333 |
// q1 = nr. of decimal digits in x |
|
|
334 |
// determine first the nr. of bits in x |
|
|
335 |
if (C1.w[1] == 0) { |
|
|
336 |
if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53 |
|
|
337 |
// split the 64-bit value in two 32-bit halves to avoid rnd errors |
|
|
338 |
if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32 |
|
|
339 |
tmp1.d = (double) (C1.w[0] >> 32); // exact conversion |
|
|
340 |
x_nr_bits = |
|
|
341 |
33 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - |
|
|
342 |
0x3ff); |
|
|
343 |
} else { // x < 2^32 |
|
|
344 |
tmp1.d = (double) (C1.w[0]); // exact conversion |
|
|
345 |
x_nr_bits = |
|
|
346 |
1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - |
|
|
347 |
0x3ff); |
|
|
348 |
} |
|
|
349 |
} else { // if x < 2^53 |
|
|
350 |
tmp1.d = (double) C1.w[0]; // exact conversion |
|
|
351 |
x_nr_bits = |
|
|
352 |
1 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); |
|
|
353 |
} |
|
|
354 |
} else { // C1.w[1] != 0 => nr. bits = 64 + nr_bits (C1.w[1]) |
|
|
355 |
tmp1.d = (double) C1.w[1]; // exact conversion |
|
|
356 |
x_nr_bits = |
|
|
357 |
65 + ((((unsigned int) (tmp1.ui64 >> 52)) & 0x7ff) - 0x3ff); |
|
|
358 |
} |
|
|
359 |
q1 = nr_digits[x_nr_bits - 1].digits; |
|
|
360 |
if (q1 == 0) { |
|
|
361 |
q1 = nr_digits[x_nr_bits - 1].digits1; |
|
|
362 |
if (C1.w[1] > nr_digits[x_nr_bits - 1].threshold_hi |
|
|
363 |
|| (C1.w[1] == nr_digits[x_nr_bits - 1].threshold_hi |
|
|
364 |
&& C1.w[0] >= nr_digits[x_nr_bits - 1].threshold_lo)) |
|
|
365 |
q1++; |
|
|
366 |
} |
|
|
367 |
// if q1 < P then pad the significand with zeros |
|
|
368 |
if (q1 < P34) { |
|
|
369 |
exp = (x_exp >> 49) - 6176; |
|
|
370 |
if (exp + 6176 > P34 - q1) { |
|
|
371 |
ind = P34 - q1; // 1 <= ind <= P34 - 1 |
|
|
372 |
// pad with P34 - q1 zeros, until exponent = emin |
|
|
373 |
// C1 = C1 * 10^ind |
|
|
374 |
if (q1 <= 19) { // 64-bit C1 |
|
|
375 |
if (ind <= 19) { // 64-bit 10^ind and 64-bit C1 |
|
|
376 |
__mul_64x64_to_128MACH (C1, C1.w[0], ten2k64[ind]); |
|
|
377 |
} else { // 128-bit 10^ind and 64-bit C1 |
|
|
378 |
__mul_128x64_to_128 (C1, C1.w[0], ten2k128[ind - 20]); |
|
|
379 |
} |
|
|
380 |
} else { // C1 is (most likely) 128-bit |
|
|
381 |
if (ind <= 14) { // 64-bit 10^ind and 128-bit C1 (most likely) |
|
|
382 |
__mul_128x64_to_128 (C1, ten2k64[ind], C1); |
|
|
383 |
} else if (ind <= 19) { // 64-bit 10^ind and 64-bit C1 (q1 <= 19) |
|
|
384 |
__mul_64x64_to_128MACH (C1, C1.w[0], ten2k64[ind]); |
|
|
385 |
} else { // 128-bit 10^ind and 64-bit C1 (C1 must be 64-bit) |
|
|
386 |
__mul_128x64_to_128 (C1, C1.w[0], ten2k128[ind - 20]); |
|
|
387 |
} |
|
|
388 |
} |
|
|
389 |
x_exp = x_exp - ((UINT64) ind << 49); |
|
|
390 |
} else { // pad with zeros until the exponent reaches emin |
|
|
391 |
ind = exp + 6176; |
|
|
392 |
// C1 = C1 * 10^ind |
|
|
393 |
if (ind <= 19) { // 1 <= P34 - q1 <= 19 <=> 15 <= q1 <= 33 |
|
|
394 |
if (q1 <= 19) { // 64-bit C1, 64-bit 10^ind |
|
|
395 |
__mul_64x64_to_128MACH (C1, C1.w[0], ten2k64[ind]); |
|
|
396 |
} else { // 20 <= q1 <= 33 => 128-bit C1, 64-bit 10^ind |
|
|
397 |
__mul_128x64_to_128 (C1, ten2k64[ind], C1); |
|
|
398 |
} |
|
|
399 |
} else { // if 20 <= P34 - q1 <= 33 <=> 1 <= q1 <= 14 => |
|
|
400 |
// 64-bit C1, 128-bit 10^ind |
|
|
401 |
__mul_128x64_to_128 (C1, C1.w[0], ten2k128[ind - 20]); |
|
|
402 |
} |
|
|
403 |
x_exp = EXP_MIN; |
|
|
404 |
} |
|
|
405 |
} |
|
|
406 |
if (x_sign) { // x < 0 |
|
|
407 |
// add 1 ulp (add 1 to the significand) |
|
|
408 |
C1.w[0]++; |
|
|
409 |
if (C1.w[0] == 0) |
|
|
410 |
C1.w[1]++; |
|
|
411 |
if (C1.w[1] == 0x0001ed09bead87c0ull && C1.w[0] == 0x378d8e6400000000ull) { // if C1 = 10^34 |
|
|
412 |
C1.w[1] = 0x0000314dc6448d93ull; // C1 = 10^33 |
|
|
413 |
C1.w[0] = 0x38c15b0a00000000ull; |
|
|
414 |
x_exp = x_exp + EXP_P1; |
|
|
415 |
} |
|
|
416 |
} else { // x > 0 |
|
|
417 |
// subtract 1 ulp (subtract 1 from the significand) |
|
|
418 |
C1.w[0]--; |
|
|
419 |
if (C1.w[0] == 0xffffffffffffffffull) |
|
|
420 |
C1.w[1]--; |
|
|
421 |
if (x_exp != 0 && C1.w[1] == 0x0000314dc6448d93ull && C1.w[0] == 0x38c15b09ffffffffull) { // if C1 = 10^33 - 1 |
|
|
422 |
C1.w[1] = 0x0001ed09bead87c0ull; // C1 = 10^34 - 1 |
|
|
423 |
C1.w[0] = 0x378d8e63ffffffffull; |
|
|
424 |
x_exp = x_exp - EXP_P1; |
|
|
425 |
} |
|
|
426 |
} |
|
|
427 |
// assemble the result |
|
|
428 |
res.w[1] = x_sign | x_exp | C1.w[1]; |
|
|
429 |
res.w[0] = C1.w[0]; |
|
|
430 |
} // end -MAXFP <= x <= -MINFP - 1 ulp OR MINFP <= x <= MAXFP - 1 ulp |
|
|
431 |
} // end x is not special and is not zero |
|
|
432 |
BID_RETURN (res); |
|
|
433 |
} |
|
|
434 |
|
|
|
435 |
/***************************************************************************** |
|
|
436 |
* BID128 nextafter |
|
|
437 |
****************************************************************************/ |
|
|
438 |
|
|
|
439 |
#if DECIMAL_CALL_BY_REFERENCE |
|
|
440 |
void |
|
|
441 |
bid128_nextafter (UINT128 * pres, UINT128 * px, |
|
|
442 |
UINT128 * |
|
|
443 |
py _EXC_FLAGS_PARAM _EXC_MASKS_PARAM _EXC_INFO_PARAM) |
|
|
444 |
{ |
|
|
445 |
UINT128 x = *px; |
|
|
446 |
UINT128 y = *py; |
|
|
447 |
UINT128 xnswp = *px; |
|
|
448 |
UINT128 ynswp = *py; |
|
|
449 |
#else |
|
|
450 |
UINT128 |
|
|
451 |
bid128_nextafter (UINT128 x, |
|
|
452 |
UINT128 y _EXC_FLAGS_PARAM _EXC_MASKS_PARAM |
|
|
453 |
_EXC_INFO_PARAM) { |
|
|
454 |
UINT128 xnswp = x; |
|
|
455 |
UINT128 ynswp = y; |
|
|
456 |
#endif |
|
|
457 |
|
|
|
458 |
UINT128 res; |
|
|
459 |
UINT128 tmp1, tmp2, tmp3; |
|
|
460 |
FPSC tmp_fpsf = 0; // dummy fpsf for calls to comparison functions |
|
|
461 |
int res1, res2; |
|
|
462 |
UINT64 x_exp; |
|
|
463 |
|
|
|
464 |
|
|
|
465 |
BID_SWAP128 (x); |
|
|
466 |
BID_SWAP128 (y); |
|
|
467 |
// check for NaNs |
|
|
468 |
if (((x.w[1] & MASK_SPECIAL) == MASK_SPECIAL) |
|
|
469 |
|| ((y.w[1] & MASK_SPECIAL) == MASK_SPECIAL)) { |
|
|
470 |
// x is special or y is special |
|
|
471 |
if ((x.w[1] & MASK_NAN) == MASK_NAN) { // x is NAN |
|
|
472 |
// if x = NaN, then res = Q (x) |
|
|
473 |
// check first for non-canonical NaN payload |
|
|
474 |
if (((x.w[1] & 0x00003fffffffffffull) > 0x0000314dc6448d93ull) || |
|
|
475 |
(((x.w[1] & 0x00003fffffffffffull) == 0x0000314dc6448d93ull) |
|
|
476 |
&& (x.w[0] > 0x38c15b09ffffffffull))) { |
|
|
477 |
x.w[1] = x.w[1] & 0xffffc00000000000ull; |
|
|
478 |
x.w[0] = 0x0ull; |
|
|
479 |
} |
|
|
480 |
if ((x.w[1] & MASK_SNAN) == MASK_SNAN) { // x is SNAN |
|
|
481 |
// set invalid flag |
|
|
482 |
*pfpsf |= INVALID_EXCEPTION; |
|
|
483 |
// return quiet (x) |
|
|
484 |
res.w[1] = x.w[1] & 0xfc003fffffffffffull; // clear out also G[6]-G[16] |
|
|
485 |
res.w[0] = x.w[0]; |
|
|
486 |
} else { // x is QNaN |
|
|
487 |
// return x |
|
|
488 |
res.w[1] = x.w[1] & 0xfc003fffffffffffull; // clear out G[6]-G[16] |
|
|
489 |
res.w[0] = x.w[0]; |
|
|
490 |
if ((y.w[1] & MASK_SNAN) == MASK_SNAN) { // y is SNAN |
|
|
491 |
// set invalid flag |
|
|
492 |
*pfpsf |= INVALID_EXCEPTION; |
|
|
493 |
} |
|
|
494 |
} |
|
|
495 |
BID_RETURN (res) |
|
|
496 |
} else if ((y.w[1] & MASK_NAN) == MASK_NAN) { // y is NAN |
|
|
497 |
// if x = NaN, then res = Q (x) |
|
|
498 |
// check first for non-canonical NaN payload |
|
|
499 |
if (((y.w[1] & 0x00003fffffffffffull) > 0x0000314dc6448d93ull) || |
|
|
500 |
(((y.w[1] & 0x00003fffffffffffull) == 0x0000314dc6448d93ull) |
|
|
501 |
&& (y.w[0] > 0x38c15b09ffffffffull))) { |
|
|
502 |
y.w[1] = y.w[1] & 0xffffc00000000000ull; |
|
|
503 |
y.w[0] = 0x0ull; |
|
|
504 |
} |
|
|
505 |
if ((y.w[1] & MASK_SNAN) == MASK_SNAN) { // y is SNAN |
|
|
506 |
// set invalid flag |
|
|
507 |
*pfpsf |= INVALID_EXCEPTION; |
|
|
508 |
// return quiet (x) |
|
|
509 |
res.w[1] = y.w[1] & 0xfc003fffffffffffull; // clear out also G[6]-G[16] |
|
|
510 |
res.w[0] = y.w[0]; |
|
|
511 |
} else { // x is QNaN |
|
|
512 |
// return x |
|
|
513 |
res.w[1] = y.w[1] & 0xfc003fffffffffffull; // clear out G[6]-G[16] |
|
|
514 |
res.w[0] = y.w[0]; |
|
|
515 |
} |
|
|
516 |
BID_RETURN (res) |
|
|
517 |
} else { // at least one is infinity |
|
|
518 |
if ((x.w[1] & MASK_ANY_INF) == MASK_INF) { // x = inf |
|
|
519 |
x.w[1] = x.w[1] & (MASK_SIGN | MASK_INF); |
|
|
520 |
x.w[0] = 0x0ull; |
|
|
521 |
} |
|
|
522 |
if ((y.w[1] & MASK_ANY_INF) == MASK_INF) { // y = inf |
|
|
523 |
y.w[1] = y.w[1] & (MASK_SIGN | MASK_INF); |
|
|
524 |
y.w[0] = 0x0ull; |
|
|
525 |
} |
|
|
526 |
} |
|
|
527 |
} |
|
|
528 |
// neither x nor y is NaN |
|
|
529 |
|
|
|
530 |
// if not infinity, check for non-canonical values x (treated as zero) |
|
|
531 |
if ((x.w[1] & MASK_ANY_INF) != MASK_INF) { // x != inf |
|
|
532 |
if ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull) { // G0_G1=11 |
|
|
533 |
// non-canonical |
|
|
534 |
x_exp = (x.w[1] << 2) & MASK_EXP; // biased and shifted left 49 bits |
|
|
535 |
x.w[1] = (x.w[1] & MASK_SIGN) | x_exp; |
|
|
536 |
x.w[0] = 0x0ull; |
|
|
537 |
} else { // G0_G1 != 11 |
|
|
538 |
x_exp = x.w[1] & MASK_EXP; // biased and shifted left 49 bits |
|
|
539 |
if ((x.w[1] & MASK_COEFF) > 0x0001ed09bead87c0ull || |
|
|
540 |
((x.w[1] & MASK_COEFF) == 0x0001ed09bead87c0ull |
|
|
541 |
&& x.w[0] > 0x378d8e63ffffffffull)) { |
|
|
542 |
// x is non-canonical if coefficient is larger than 10^34 -1 |
|
|
543 |
x.w[1] = (x.w[1] & MASK_SIGN) | x_exp; |
|
|
544 |
x.w[0] = 0x0ull; |
|
|
545 |
} else { // canonical |
|
|
546 |
; |
|
|
547 |
} |
|
|
548 |
} |
|
|
549 |
} |
|
|
550 |
// no need to check for non-canonical y |
|
|
551 |
|
|
|
552 |
// neither x nor y is NaN |
|
|
553 |
tmp_fpsf = *pfpsf; // save fpsf |
|
|
554 |
#if DECIMAL_CALL_BY_REFERENCE |
|
|
555 |
bid128_quiet_equal (&res1, &xnswp, |
|
|
556 |
&ynswp _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
557 |
_EXC_INFO_ARG); |
|
|
558 |
bid128_quiet_greater (&res2, &xnswp, |
|
|
559 |
&ynswp _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
560 |
_EXC_INFO_ARG); |
|
|
561 |
#else |
|
|
562 |
res1 = |
|
|
563 |
bid128_quiet_equal (xnswp, |
|
|
564 |
ynswp _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
565 |
_EXC_INFO_ARG); |
|
|
566 |
res2 = |
|
|
567 |
bid128_quiet_greater (xnswp, |
|
|
568 |
ynswp _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
569 |
_EXC_INFO_ARG); |
|
|
570 |
#endif |
|
|
571 |
*pfpsf = tmp_fpsf; // restore fpsf |
|
|
572 |
|
|
|
573 |
if (res1) { // x = y |
|
|
574 |
// return x with the sign of y |
|
|
575 |
res.w[1] = |
|
|
576 |
(x.w[1] & 0x7fffffffffffffffull) | (y. |
|
|
577 |
w[1] & 0x8000000000000000ull); |
|
|
578 |
res.w[0] = x.w[0]; |
|
|
579 |
} else if (res2) { // x > y |
|
|
580 |
#if DECIMAL_CALL_BY_REFERENCE |
|
|
581 |
bid128_nextdown (&res, |
|
|
582 |
&xnswp _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
583 |
_EXC_INFO_ARG); |
|
|
584 |
#else |
|
|
585 |
res = |
|
|
586 |
bid128_nextdown (xnswp _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
587 |
_EXC_INFO_ARG); |
|
|
588 |
#endif |
|
|
589 |
BID_SWAP128 (res); |
|
|
590 |
} else { // x < y |
|
|
591 |
#if DECIMAL_CALL_BY_REFERENCE |
|
|
592 |
bid128_nextup (&res, |
|
|
593 |
&xnswp _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG); |
|
|
594 |
#else |
|
|
595 |
res = |
|
|
596 |
bid128_nextup (xnswp _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG); |
|
|
597 |
#endif |
|
|
598 |
BID_SWAP128 (res); |
|
|
599 |
} |
|
|
600 |
// if the operand x is finite but the result is infinite, signal |
|
|
601 |
// overflow and inexact |
|
|
602 |
if (((x.w[1] & MASK_SPECIAL) != MASK_SPECIAL) |
|
|
603 |
&& ((res.w[1] & MASK_SPECIAL) == MASK_SPECIAL)) { |
|
|
604 |
// set the inexact flag |
|
|
605 |
*pfpsf |= INEXACT_EXCEPTION; |
|
|
606 |
// set the overflow flag |
|
|
607 |
*pfpsf |= OVERFLOW_EXCEPTION; |
|
|
608 |
} |
|
|
609 |
// if the result is in (-10^emin, 10^emin), and is different from the |
|
|
610 |
// operand x, signal underflow and inexact |
|
|
611 |
tmp1.w[HIGH_128W] = 0x0000314dc6448d93ull; |
|
|
612 |
tmp1.w[LOW_128W] = 0x38c15b0a00000000ull; // +100...0[34] * 10^emin |
|
|
613 |
tmp2.w[HIGH_128W] = res.w[1] & 0x7fffffffffffffffull; |
|
|
614 |
tmp2.w[LOW_128W] = res.w[0]; |
|
|
615 |
tmp3.w[HIGH_128W] = res.w[1]; |
|
|
616 |
tmp3.w[LOW_128W] = res.w[0]; |
|
|
617 |
tmp_fpsf = *pfpsf; // save fpsf |
|
|
618 |
#if DECIMAL_CALL_BY_REFERENCE |
|
|
619 |
bid128_quiet_greater (&res1, &tmp1, |
|
|
620 |
&tmp2 _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
621 |
_EXC_INFO_ARG); |
|
|
622 |
bid128_quiet_not_equal (&res2, &xnswp, |
|
|
623 |
&tmp3 _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
624 |
_EXC_INFO_ARG); |
|
|
625 |
#else |
|
|
626 |
res1 = |
|
|
627 |
bid128_quiet_greater (tmp1, |
|
|
628 |
tmp2 _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
629 |
_EXC_INFO_ARG); |
|
|
630 |
res2 = |
|
|
631 |
bid128_quiet_not_equal (xnswp, |
|
|
632 |
tmp3 _EXC_FLAGS_ARG _EXC_MASKS_ARG |
|
|
633 |
_EXC_INFO_ARG); |
|
|
634 |
#endif |
|
|
635 |
*pfpsf = tmp_fpsf; // restore fpsf |
|
|
636 |
if (res1 && res2) { |
|
|
637 |
// set the inexact flag |
|
|
638 |
*pfpsf |= INEXACT_EXCEPTION; |
|
|
639 |
// set the underflow flag |
|
|
640 |
*pfpsf |= UNDERFLOW_EXCEPTION; |
|
|
641 |
} |
|
|
642 |
BID_RETURN (res); |
|
|
643 |
} |