* [PATCH] Update libbid according to the latest Intel Decimal Floating-Point Math Library.
@ 2024-04-28 5:50 liuhongt
2024-05-02 9:20 ` Richard Biener
0 siblings, 1 reply; 2+ messages in thread
From: liuhongt @ 2024-04-28 5:50 UTC (permalink / raw)
To: gcc-patches; +Cc: crazylht, hjl.tools
The Intel Decimal Floating-Point Math Library is available as open-source on Netlib[1].
[1] https://www.netlib.org/misc/intel/.
Bootstrapped and regtested on x86_64-pc-linux-gnu{-m32,}.
Ok for trunk?
libgcc/config/libbid/ChangeLog:
* bid128_fma.c (add_and_round): Fix bug: the result
of (+5E+368)*(+10E-34)+(-10E+369) was returning
-9999999999999999999999999999999999E+336 instead of expected
result -1000000000000000000000000000000000E+337.
(bid128_ext_fma): Ditto.
(bid64qqq_fma): Ditto.
* bid128_noncomp.c: Change return type of bid128_class from
int to class_t.
* bid128_round_integral.c: Add default case to avoid compiler
warning.
* bid128_string.c (bid128_to_string): Replace 0x30 with '0'
for zero digit.
(bid128_from_string): Ditto.
* bid32_to_bid128.c (bid128_to_bid32): Fix Bug. In addition
to the INEXACT flag, the UNDERFLOW flag needs to be set (and
was not) when converting an input such as
+6931674235302037148946035460357709E+1857 to +1000000E-101
* bid32_to_bid64.c (bid64_to_bid32): fix Bug, In addition to
the INEXACT flag, the UNDERFLOW flag needs to be set (and was
not) when converting an input such as +9999999000000001E-111
to +1000000E-101. Furthermore, significant bits of NaNs are
set correctly now. For example, 0x7c00003b9aca0000 was
returning 0x7c000002 instead of 0x 7c000100.
* bid64_noncomp.c: Change return type of bid64_class from int
to class_t.
* bid64_round_integral.c (bid64_round_integral_exact): Add
default case to avoid compiler warning.
* bid64_string.c (bid64_from_string): Fix bug for rounding
up. The input string "10000000000000000" was returning
+1000000000000001E+1 instead of +1000000000000000E+1.
* bid64_to_bid128.c (bid128_to_bid64): Fix bug, in addition to
the INEXACT flag, the UNDERFLOW flag needs to be set (and was
not) when converting an input such as
+9999999999999999999999999999999999E-417 to
+1000000000000000E-398.
* bid_binarydecimal.c (bid32_to_binary64): Fix bug for
conversion between binary and bid types. For example,
0x7c0F4240 was returning 0x7FFFA12000000000 instead of
expected double precision 0x7FF8000000000000.
(binary64_to_bid32): Ditto.
(binary80_to_bid32): Ditto.
(binary128_to_bid32): Ditto.
(binary80_to_bid64): Ditto.
(binary128_to_bid64): Ditto.
* bid_conf.h (BID_HIGH_128W): New macro.
(BID_LOW_128W): Ditto.
* bid_functions.h (__ENABLE_BINARY80__): Ditto.
(ALIGN): Ditto.
* bid_inline_add.h (get_add128): Add default case to avoid compiler
warning.
* bid_internal.h (get_BID64): Ditto.
(fast_get_BID64_check_OF): Ditto.
(ALIGN): New macro.
Co-authored-by: Anderson, Cristina S <cristina.s.anderson@intel.com>
Co-authored-by: Akkas, Ahmet <ahmet.akkas@intel.com>
Co-authored-by: Cornea, Marius <marius.cornea@intel.com>
---
libgcc/config/libbid/bid128_fma.c | 188 ++++++++++---------
libgcc/config/libbid/bid128_noncomp.c | 2 +-
libgcc/config/libbid/bid128_round_integral.c | 2 +
libgcc/config/libbid/bid128_string.c | 7 +-
libgcc/config/libbid/bid32_to_bid128.c | 3 -
libgcc/config/libbid/bid32_to_bid64.c | 11 +-
libgcc/config/libbid/bid64_noncomp.c | 2 +-
libgcc/config/libbid/bid64_round_integral.c | 2 +
libgcc/config/libbid/bid64_string.c | 21 ++-
libgcc/config/libbid/bid64_to_bid128.c | 3 -
libgcc/config/libbid/bid_binarydecimal.c | 167 ++++++----------
libgcc/config/libbid/bid_conf.h | 8 +
libgcc/config/libbid/bid_functions.h | 23 ++-
libgcc/config/libbid/bid_inline_add.h | 2 +
libgcc/config/libbid/bid_internal.h | 17 +-
15 files changed, 220 insertions(+), 238 deletions(-)
diff --git a/libgcc/config/libbid/bid128_fma.c b/libgcc/config/libbid/bid128_fma.c
index 67233193a42..cbcf225546f 100644
--- a/libgcc/config/libbid/bid128_fma.c
+++ b/libgcc/config/libbid/bid128_fma.c
@@ -417,13 +417,12 @@ add_and_round (int q3,
R128.w[1] = R256.w[1];
R128.w[0] = R256.w[0];
}
+ if (e4 + x0 < expmin) { // for all rounding modes
+ is_tiny = 1;
+ }
// the rounded result has p34 = 34 digits
e4 = e4 + x0 + incr_exp;
- if (rnd_mode == ROUNDING_TO_NEAREST) {
- if (e4 < expmin) {
- is_tiny = 1; // for other rounding modes apply correction
- }
- } else {
+ if (rnd_mode != ROUNDING_TO_NEAREST) {
// for RM, RP, RZ, RA apply correction in order to determine tininess
// but do not save the result; apply the correction to
// (-1)^p_sign * significand * 10^0
@@ -434,10 +433,6 @@ add_and_round (int q3,
is_inexact_gt_midpoint, is_midpoint_lt_even,
is_midpoint_gt_even, 0, &P128, ptrfpsf);
scale = ((P128.w[1] & MASK_EXP) >> 49) - 6176; // -1, 0, or +1
- // the number of digits in the significand is p34 = 34
- if (e4 + scale < expmin) {
- is_tiny = 1;
- }
}
ind = p34; // the number of decimal digits in the signifcand of res
res.w[1] = p_sign | ((UINT64) (e4 + 6176) << 49) | R128.w[1]; // RN
@@ -851,7 +846,6 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
}
}
}
-
p_sign = x_sign ^ y_sign; // sign of the product
// identify cases where at least one operand is infinity
@@ -988,15 +982,10 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
if (C1.w[1] == 0) {
if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
// split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
+
tmp.d = (double) (C1.w[0] >> 32); // exact conversion
x_nr_bits =
33 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // x < 2^32
- tmp.d = (double) (C1.w[0]); // exact conversion
- x_nr_bits =
- 1 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
} else { // if x < 2^53
tmp.d = (double) C1.w[0]; // exact conversion
x_nr_bits =
@@ -1011,42 +1000,36 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
if (q1 == 0) {
q1 = nr_digits[x_nr_bits - 1].digits1;
if (C1.w[1] > nr_digits[x_nr_bits - 1].threshold_hi ||
- (C1.w[1] == nr_digits[x_nr_bits - 1].threshold_hi &&
- C1.w[0] >= nr_digits[x_nr_bits - 1].threshold_lo))
+ (C1.w[1] == nr_digits[x_nr_bits - 1].threshold_hi &&
+ C1.w[0] >= nr_digits[x_nr_bits - 1].threshold_lo))
q1++;
}
}
-
+ // q2 = nr. of decimal digits in y
+ // determine first the nr. of bits in y
if (C2.w[1] != 0 || C2.w[0] != 0) { // y = f (non-zero finite)
if (C2.w[1] == 0) {
if (C2.w[0] >= 0x0020000000000000ull) { // y >= 2^53
// split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C2.w[0] >= 0x0000000100000000ull) { // y >= 2^32
tmp.d = (double) (C2.w[0] >> 32); // exact conversion
y_nr_bits =
- 32 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // y < 2^32
- tmp.d = (double) C2.w[0]; // exact conversion
- y_nr_bits =
- ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
+ 33 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
} else { // if y < 2^53
tmp.d = (double) C2.w[0]; // exact conversion
y_nr_bits =
- ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
+ 1 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
}
} else { // C2.w[1] != 0 => nr. bits = 64 + nr_bits (C2.w[1])
tmp.d = (double) C2.w[1]; // exact conversion
y_nr_bits =
- 64 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
+ 65 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
}
-
- q2 = nr_digits[y_nr_bits].digits;
+ q2 = nr_digits[y_nr_bits - 1].digits;
if (q2 == 0) {
- q2 = nr_digits[y_nr_bits].digits1;
- if (C2.w[1] > nr_digits[y_nr_bits].threshold_hi ||
- (C2.w[1] == nr_digits[y_nr_bits].threshold_hi &&
- C2.w[0] >= nr_digits[y_nr_bits].threshold_lo))
+ q2 = nr_digits[y_nr_bits - 1].digits1;
+ if (C2.w[1] > nr_digits[y_nr_bits - 1].threshold_hi ||
+ (C2.w[1] == nr_digits[y_nr_bits - 1].threshold_hi &&
+ C2.w[0] >= nr_digits[y_nr_bits - 1].threshold_lo))
q2++;
}
}
@@ -1055,32 +1038,25 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
if (C3.w[1] == 0) {
if (C3.w[0] >= 0x0020000000000000ull) { // z >= 2^53
// split the 64-bit value in two 32-bit halves to avoid rounding errors
- if (C3.w[0] >= 0x0000000100000000ull) { // z >= 2^32
tmp.d = (double) (C3.w[0] >> 32); // exact conversion
z_nr_bits =
- 32 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
- } else { // z < 2^32
- tmp.d = (double) C3.w[0]; // exact conversion
- z_nr_bits =
- ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
- }
+ 33 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
} else { // if z < 2^53
tmp.d = (double) C3.w[0]; // exact conversion
z_nr_bits =
- ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
+ 1 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
}
} else { // C3.w[1] != 0 => nr. bits = 64 + nr_bits (C3.w[1])
tmp.d = (double) C3.w[1]; // exact conversion
z_nr_bits =
- 64 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
+ 65 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
}
-
- q3 = nr_digits[z_nr_bits].digits;
+ q3 = nr_digits[z_nr_bits - 1].digits;
if (q3 == 0) {
- q3 = nr_digits[z_nr_bits].digits1;
- if (C3.w[1] > nr_digits[z_nr_bits].threshold_hi ||
- (C3.w[1] == nr_digits[z_nr_bits].threshold_hi &&
- C3.w[0] >= nr_digits[z_nr_bits].threshold_lo))
+ q3 = nr_digits[z_nr_bits - 1].digits1;
+ if (C3.w[1] > nr_digits[z_nr_bits - 1].threshold_hi ||
+ (C3.w[1] == nr_digits[z_nr_bits - 1].threshold_hi &&
+ C3.w[0] >= nr_digits[z_nr_bits - 1].threshold_lo))
q3++;
}
}
@@ -1128,7 +1104,6 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
} else {
; // continue with x = f, y = f, z = 0 or x = f, y = f, z = f
}
-
e1 = (x_exp >> 49) - 6176; // unbiased exponent of x
e2 = (y_exp >> 49) - 6176; // unbiased exponent of y
e3 = (z_exp >> 49) - 6176; // unbiased exponent of z
@@ -1232,22 +1207,18 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
// length of C1 * C2 rounded up to a multiple of 64 bits is len = 192;
q4 = q1 + q2; // q4 in [40, 57]
}
- } else if (q1 + q2 == 58) { // C4 = C1 * C2 fits in 192 or 256 bits
- // both C1 and C2 fit in 128 bits (actually in 113 bits); at most one
- // may fit in 64 bits
- if (C1.w[1] == 0) { // C1 * C2 will fit in 192 bits
- __mul_64x128_full (C4.w[2], C4, C1.w[0], C2); // may use 64x128_to_192
- } else if (C2.w[1] == 0) { // C1 * C2 will fit in 192 bits
- __mul_64x128_full (C4.w[2], C4, C2.w[0], C1); // may use 64x128_to_192
- } else { // C1 * C2 will fit in 192 bits or in 256 bits
- __mul_128x128_to_256 (C4, C1, C2);
- }
+ } else if (q1 + q2 == 58) { // C4 = C1 * C2 fits in 192 or 256 bits;
+ // both C1 and C2 fit in 128 bits (actually in 113 bits); none can
+ // fit in 64 bits, because each number must have at least 24 decimal
+ // digits for the sum to have 58 (as the max. nr. of digits is 34) =>
+ // C1.w[1] != 0 and C2.w[1] != 0
+ __mul_128x128_to_256 (C4, C1, C2);
// if C4 < 10^(q1+q2-1) = 10^57 then q4 = q1+q2-1 = 57 else q4 = q1+q2 = 58
if (C4.w[3] == 0 && (C4.w[2] < ten2k256[18].w[2] ||
- (C4.w[2] == ten2k256[18].w[2]
- && (C4.w[1] < ten2k256[18].w[1]
- || (C4.w[1] == ten2k256[18].w[1]
- && C4.w[0] < ten2k256[18].w[0]))))) {
+ (C4.w[2] == ten2k256[18].w[2]
+ && (C4.w[1] < ten2k256[18].w[1]
+ || (C4.w[1] == ten2k256[18].w[1]
+ && C4.w[0] < ten2k256[18].w[0]))))) {
// 18 = 57 - 39 = q1+q2-1 - 39
// length of C1 * C2 rounded up to a multiple of 64 bits is len = 192;
q4 = 57; // 57 = q1 + q2 - 1
@@ -1283,7 +1254,6 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
q4 = q1 + q2; // q4 in [59, 68]
}
}
-
if (C3.w[1] == 0x0 && C3.w[0] == 0x0) { // x = f, y = f, z = 0
save_fpsf = *pfpsf; // sticky bits - caller value must be preserved
*pfpsf = 0;
@@ -1319,10 +1289,11 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
res.w[1] = R256.w[1];
}
e4 = e4 + x0;
+ q4 = p34;
if (incr_exp) {
e4 = e4 + 1;
+ if (q4 + e4 == expmin + p34) *pfpsf |= (INEXACT_EXCEPTION | UNDERFLOW_EXCEPTION);
}
- q4 = p34;
// res is now the coefficient of the result rounded to the destination
// precision, with unbounded exponent; the exponent is e4; q4=digits(res)
} else { // if (q4 <= p34)
@@ -1648,7 +1619,6 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
delta = q3 + e3 - q4 - e4;
delta_ge_zero:
if (delta >= 0) {
-
if (p34 <= delta - 1 || // Case (1')
(p34 == delta && e3 + 6176 < p34 - q3)) { // Case (1''A)
// check for overflow, which can occur only in Case (1')
@@ -1736,7 +1706,7 @@ delta_ge_zero:
res.w[1] = z_sign | ((UINT64) (e3 + 6176) << 49) | C3.w[1];
res.w[0] = C3.w[0];
}
-
+
// use the following to avoid double rounding errors when operating on
// mixed formats in rounding to nearest, and for correcting the result
// if not rounding to nearest
@@ -1795,7 +1765,10 @@ delta_ge_zero:
R64 = 10;
}
}
- if (q4 == 1 && C4.w[0] == 5) {
+
+ if (R64 == 5 && !is_inexact_lt_midpoint && !is_inexact_gt_midpoint &&
+ !is_midpoint_lt_even && !is_midpoint_gt_even) {
+ //if (q4 == 1 && C4.w[0] == 5) {
is_inexact_lt_midpoint = 0;
is_inexact_gt_midpoint = 0;
is_midpoint_lt_even = 1;
@@ -1826,11 +1799,7 @@ delta_ge_zero:
res.w[1] = z_sign | ((UINT64) (e3 + 6176) << 49) | res.w[1];
}
if (e3 == expmin) {
- if (R64 < 5 || (R64 == 5 && !is_inexact_lt_midpoint)) {
- ; // result not tiny (in round-to-nearest mode)
- } else {
- *pfpsf |= UNDERFLOW_EXCEPTION;
- }
+ *pfpsf |= UNDERFLOW_EXCEPTION; // tiny if detected before rounding
}
} // end 10^(q3+scale-1)
// set the inexact flag
@@ -1877,10 +1846,9 @@ delta_ge_zero:
// endif
if ((e3 == expmin && (q3 + scale) < p34) ||
(e3 == expmin && (q3 + scale) == p34 &&
- (res.w[1] & MASK_COEFF) == 0x0000314dc6448d93ull && // 10^33_high
- res.w[0] == 0x38c15b0a00000000ull && // 10^33_low
- z_sign != p_sign && ((!z_sign && rnd_mode != ROUNDING_UP) ||
- (z_sign && rnd_mode != ROUNDING_DOWN)))) {
+ (res.w[1] & MASK_COEFF) == 0x0000314dc6448d93ull && // 10^33_high
+ res.w[0] == 0x38c15b0a00000000ull && // 10^33_low
+ z_sign != p_sign)) {
*pfpsf |= UNDERFLOW_EXCEPTION;
}
if (rnd_mode != ROUNDING_TO_NEAREST) {
@@ -2594,7 +2562,7 @@ delta_ge_zero:
if (e3 > expmin && ((res.w[1] < 0x0000314dc6448d93ull ||
(res.w[1] == 0x0000314dc6448d93ull &&
res.w[0] < 0x38c15b0a00000000ull)) ||
- (is_inexact_lt_midpoint
+ ((is_inexact_lt_midpoint | is_midpoint_gt_even)
&& res.w[1] == 0x0000314dc6448d93ull
&& res.w[0] == 0x38c15b0a00000000ull))
&& x0 >= 1) {
@@ -2678,6 +2646,9 @@ delta_ge_zero:
res.w[0] < 0x38c15b0a00000000ull)) {
is_tiny = 1;
}
+ if (((res.w[1] & 0x7fffffffffffffffull) == 0x0000314dc6448d93ull) &&
+ (res.w[0] == 0x38c15b0a00000000ull) && // 10^33*10^-6176
+ (z_sign != p_sign)) is_tiny = 1;
} else if (e3 < expmin) {
// the result is tiny, so we must truncate more of res
is_tiny = 1;
@@ -3328,9 +3299,6 @@ delta_ge_zero:
0, &P128, pfpsf);
scale = ((P128.w[1] & MASK_EXP) >> 49) - 6176; // -1, 0, or +1
// the number of digits in the significand is p34 = 34
- if (e4 + scale < expmin) {
- is_tiny = 1;
- }
}
// the result rounded to the destination precision with unbounded exponent
@@ -3521,6 +3489,19 @@ delta_ge_zero:
is_midpoint_lt_even, is_midpoint_gt_even,
e4, &res, pfpsf);
}
+ // correction needed for tininess detection before rounding
+ if ((((res.w[1] & 0x7fffffffffffffffull) == 0x0000314dc6448d93ull) &&
+ // 10^33*10^-6176_high
+ (res.w[0] == 0x38c15b0a00000000ull)) && // 10^33*10^-6176_low
+ (((rnd_mode == ROUNDING_TO_NEAREST ||
+ rnd_mode == ROUNDING_TIES_AWAY) &&
+ (is_midpoint_lt_even || is_inexact_gt_midpoint)) ||
+ ((((rnd_mode == ROUNDING_UP) && !(res.w[1] & MASK_SIGN)) ||
+ ((rnd_mode == ROUNDING_DOWN) && (res.w[1] & MASK_SIGN)))
+ && (is_midpoint_lt_even || is_midpoint_gt_even ||
+ is_inexact_lt_midpoint || is_inexact_gt_midpoint)))) {
+ is_tiny = 1;
+ }
if (is_midpoint_lt_even || is_midpoint_gt_even ||
is_inexact_lt_midpoint || is_inexact_gt_midpoint) {
// set the inexact flag
@@ -4162,21 +4143,34 @@ bid64qqq_fma (UINT128 x, UINT128 y, UINT128 z
// determine the unbiased exponent of the result
unbexp = ((res1 >> 53) & 0x3ff) - 398;
+ if (!((res1 & MASK_NAN) == MASK_NAN)) { // res1 not NaN
// if subnormal, res1 must have exp = -398
// if tiny and inexact set underflow and inexact status flags
- if (!((res1 & MASK_NAN) == MASK_NAN) && // res1 not NaN
- (unbexp == -398)
- && ((res1 & MASK_BINARY_SIG1) < 1000000000000000ull)
- && (is_inexact_lt_midpoint0 || is_inexact_gt_midpoint0
- || is_midpoint_lt_even0 || is_midpoint_gt_even0)) {
- // set the inexact flag and the underflow flag
- *pfpsf |= (INEXACT_EXCEPTION | UNDERFLOW_EXCEPTION);
+ if ((unbexp == -398)
+ && ((res1 & MASK_BINARY_SIG1) < 1000000000000000ull)
+ && (is_inexact_lt_midpoint0 || is_inexact_gt_midpoint0
+ || is_midpoint_lt_even0 || is_midpoint_gt_even0)) {
+ // set the inexact flag and the underflow flag
+ *pfpsf |= (INEXACT_EXCEPTION | UNDERFLOW_EXCEPTION);
} else if (is_inexact_lt_midpoint0 || is_inexact_gt_midpoint0 ||
is_midpoint_lt_even0 || is_midpoint_gt_even0) {
// set the inexact flag and the underflow flag
*pfpsf |= INEXACT_EXCEPTION;
- }
+ }
+ // correction needed for tininess detection before rounding
+ if (((res1 & 0x7fffffffffffffffull) == 1000000000000000ull) &&
+ // 10^15*10^-398
+ (((rnd_mode == ROUNDING_TO_NEAREST ||
+ rnd_mode == ROUNDING_TIES_AWAY) &&
+ (is_midpoint_lt_even || is_inexact_gt_midpoint)) ||
+ ((((rnd_mode == ROUNDING_UP) && !(res1 & MASK_SIGN)) ||
+ ((rnd_mode == ROUNDING_DOWN) && (res1 & MASK_SIGN)))
+ && (is_midpoint_lt_even || is_midpoint_gt_even ||
+ is_inexact_lt_midpoint || is_inexact_gt_midpoint)))) {
+ *pfpsf |= UNDERFLOW_EXCEPTION;
+ }
+ }
*pfpsf |= save_fpsf;
BID_RETURN (res1);
} // else continue, and use rounding to nearest to round to 16 digits
@@ -4453,6 +4447,20 @@ bid64qqq_fma (UINT128 x, UINT128 y, UINT128 z
res1 = sign | MASK_STEERING_BITS |
((UINT64) (unbexp + 398) << 51) | (res1 & MASK_BINARY_SIG2);
}
+
+ // correction needed for tininess detection before rounding
+ if (((res1 & 0x7fffffffffffffffull) == 1000000000000000ull) &&
+ // 10^15*10^-398
+ (((rnd_mode == ROUNDING_TO_NEAREST ||
+ rnd_mode == ROUNDING_TIES_AWAY) &&
+ (is_midpoint_lt_even || is_inexact_gt_midpoint)) ||
+ ((((rnd_mode == ROUNDING_UP) && !(res1 & MASK_SIGN)) ||
+ ((rnd_mode == ROUNDING_DOWN) && (res1 & MASK_SIGN)))
+ && (is_midpoint_lt_even || is_midpoint_gt_even ||
+ is_inexact_lt_midpoint || is_inexact_gt_midpoint)))) {
+ *pfpsf |= UNDERFLOW_EXCEPTION;
+ }
+
*pfpsf |= save_fpsf;
BID_RETURN (res1);
}
diff --git a/libgcc/config/libbid/bid128_noncomp.c b/libgcc/config/libbid/bid128_noncomp.c
index a79ac859ce1..4ef166c81dc 100644
--- a/libgcc/config/libbid/bid128_noncomp.c
+++ b/libgcc/config/libbid/bid128_noncomp.c
@@ -443,7 +443,7 @@ void
bid128_class (int *pres, UINT128 * px _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
UINT128 x = *px;
#else
-int
+class_t
bid128_class (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
#endif
int res;
diff --git a/libgcc/config/libbid/bid128_round_integral.c b/libgcc/config/libbid/bid128_round_integral.c
index ddaa0b818f2..051f8268aa9 100644
--- a/libgcc/config/libbid/bid128_round_integral.c
+++ b/libgcc/config/libbid/bid128_round_integral.c
@@ -177,6 +177,7 @@ case ROUNDING_TO_ZERO:
BID_RETURN (res);
}
break;
+default: break; // default added to avoid compiler warning
}
// q = nr. of decimal digits in x
@@ -804,6 +805,7 @@ case ROUNDING_TO_ZERO:
BID_RETURN (res);
}
break;
+default: break; // default added to avoid compiler warning
}
BID_RETURN (res);
diff --git a/libgcc/config/libbid/bid128_string.c b/libgcc/config/libbid/bid128_string.c
index ecd295cfbdf..8fc12ee2d76 100644
--- a/libgcc/config/libbid/bid128_string.c
+++ b/libgcc/config/libbid/bid128_string.c
@@ -56,6 +56,7 @@ bid128_to_string (char *str, UINT128 x
UINT128 C1;
unsigned int k = 0; // pointer in the string
unsigned int d0, d123;
+ unsigned int zero_digit = (unsigned int) '0';
UINT64 HI_18Dig, LO_18Dig, Tmp;
UINT32 MiDi[12], *ptr;
char *c_ptr_start, *c_ptr;
@@ -232,14 +233,14 @@ bid128_to_string (char *str, UINT128 x
d123 = exp - 1000 * d0;
if (d0) { // 1000 <= exp <= 6144 => 4 digits to return
- str[k++] = d0 + 0x30;// ASCII for decimal digit d0
+ str[k++] = d0 + zero_digit; // ASCII for decimal digit d0
ind = 3 * d123;
str[k++] = char_table3[ind];
str[k++] = char_table3[ind + 1];
str[k++] = char_table3[ind + 2];
} else { // 0 <= exp <= 999 => d0 = 0
if (d123 < 10) { // 0 <= exp <= 9 => 1 digit to return
- str[k++] = d123 + 0x30;// ASCII
+ str[k++] = d123 + zero_digit; // ASCII
} else if (d123 < 100) { // 10 <= exp <= 99 => 2 digits to return
ind = 2 * (d123 - 10);
str[k++] = char_table2[ind];
@@ -643,7 +644,7 @@ bid128_from_string (char *ps _RND_MODE_PARAM _EXC_FLAGS_PARAM
}
break;
-
+ default: break; // default added to avoid compiler warning
}
// now form the coefficient as coeff_high*10^17+coeff_low+carry
scale_high = 100000000000000000ull;
diff --git a/libgcc/config/libbid/bid32_to_bid128.c b/libgcc/config/libbid/bid32_to_bid128.c
index d1d1d3458fd..5b5ce9504e1 100644
--- a/libgcc/config/libbid/bid32_to_bid128.c
+++ b/libgcc/config/libbid/bid32_to_bid128.c
@@ -155,9 +155,6 @@ bid128_to_bid32 (UINT128 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
T128 = round_const_table_128[rmode][extra_digits];
__add_carry_out (CX1.w[0], carry, T128.w[0], CX.w[0]);
CX1.w[1] = CX.w[1] + T128.w[1] + carry;
- if (__unsigned_compare_ge_128
- (CX1, power10_table_128[extra_digits + 7]))
- uf_check = 0;
}
extra_digits =
extra_digits + DECIMAL_EXPONENT_BIAS_128 -
diff --git a/libgcc/config/libbid/bid32_to_bid64.c b/libgcc/config/libbid/bid32_to_bid64.c
index 7802346a3d1..61b24b29915 100644
--- a/libgcc/config/libbid/bid32_to_bid64.c
+++ b/libgcc/config/libbid/bid32_to_bid64.c
@@ -79,6 +79,7 @@ bid64_to_bid32 (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
UINT128 Q;
UINT64 sign_x, coefficient_x, remainder_h, carry, Stemp;
UINT32 res;
+ UINT64 t64;
int_float tempx;
int exponent_x, bin_expon_cx, extra_digits, rmode = 0, amount;
unsigned status = 0;
@@ -93,8 +94,10 @@ bid64_to_bid32 (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
// unpack arguments, check for NaN or Infinity, 0
if (!unpack_BID64 (&sign_x, &exponent_x, &coefficient_x, x)) {
if (((x) & 0x7800000000000000ull) == 0x7800000000000000ull) {
- res = (coefficient_x & 0x0003ffffffffffffull);
- res /= 1000000000ull;
+ t64 = (coefficient_x & 0x0003ffffffffffffull);
+ res = t64/1000000000ull;
+ //res = (coefficient_x & 0x0003ffffffffffffull);
+ //res /= 1000000000ull;
res |= ((coefficient_x >> 32) & 0xfc000000);
#ifdef SET_STATUS_FLAGS
if ((x & SNAN_MASK64) == SNAN_MASK64) // sNaN
@@ -139,10 +142,6 @@ bid64_to_bid32 (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
exponent_x += extra_digits;
if ((exponent_x < 0) && (exponent_x + MAX_FORMAT_DIGITS_32 >= 0)) {
status = UNDERFLOW_EXCEPTION;
- if (exponent_x == -1)
- if (coefficient_x + round_const_table[rmode][extra_digits] >=
- power10_table_128[extra_digits + 7].w[0])
- status = 0;
extra_digits -= exponent_x;
exponent_x = 0;
}
diff --git a/libgcc/config/libbid/bid64_noncomp.c b/libgcc/config/libbid/bid64_noncomp.c
index ecc9dddfd8b..ec633693f62 100644
--- a/libgcc/config/libbid/bid64_noncomp.c
+++ b/libgcc/config/libbid/bid64_noncomp.c
@@ -358,7 +358,7 @@ void
bid64_class (int *pres, UINT64 * px _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
UINT64 x = *px;
#else
-int
+class_t
bid64_class (UINT64 x _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
#endif
int res;
diff --git a/libgcc/config/libbid/bid64_round_integral.c b/libgcc/config/libbid/bid64_round_integral.c
index 99a3e6101ef..5a33c429116 100644
--- a/libgcc/config/libbid/bid64_round_integral.c
+++ b/libgcc/config/libbid/bid64_round_integral.c
@@ -142,6 +142,7 @@ bid64_round_integral_exact (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
BID_RETURN (res);
}
break;
+ default: break; // default added to avoid compiler warning
} // end switch ()
// q = nr. of decimal digits in x (1 <= q <= 54)
@@ -483,6 +484,7 @@ bid64_round_integral_exact (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
BID_RETURN (res);
}
break;
+ default: break; // default added to avoid compiler warning
} // end switch ()
BID_RETURN (res);
}
diff --git a/libgcc/config/libbid/bid64_string.c b/libgcc/config/libbid/bid64_string.c
index 55fcd9e9208..81ac5e275a4 100644
--- a/libgcc/config/libbid/bid64_string.c
+++ b/libgcc/config/libbid/bid64_string.c
@@ -251,7 +251,7 @@ bid64_from_string (char *ps
#endif
UINT64 sign_x, coefficient_x = 0, rounded = 0, res;
int expon_x = 0, sgn_expon, ndigits, add_expon = 0, midpoint =
- 0, rounded_up = 0;
+ 0, rounded_up = 0, dround = 0;
int dec_expon_scale = 0, right_radix_leading_zeros = 0, rdx_pt_enc =
0;
unsigned fpsc;
@@ -419,10 +419,10 @@ bid64_from_string (char *ps
break;
case ROUNDING_DOWN:
- if(sign_x) { coefficient_x++; rounded_up=1; }
+ if(sign_x) { if(c>'0') {coefficient_x++; rounded_up=1;} else dround=1; }
break;
case ROUNDING_UP:
- if(!sign_x) { coefficient_x++; rounded_up=1; }
+ if(!sign_x) { if(c>'0') {coefficient_x++; rounded_up=1;} else dround=1; }
break;
case ROUNDING_TIES_AWAY:
if(c>='5') { coefficient_x++; rounded_up=1; }
@@ -443,8 +443,21 @@ bid64_from_string (char *ps
midpoint = 0;
rounded_up = 1;
}
- if (c > '0')
+ if (c > '0') {
rounded = 1;
+
+ if(dround)
+ {
+ dround = 0;
+ coefficient_x ++;
+ rounded_up = 1;
+
+ if (coefficient_x == 10000000000000000ull) {
+ coefficient_x = 1000000000000000ull;
+ add_expon = 1;
+ }
+ }
+ }
}
ps++;
c = *ps;
diff --git a/libgcc/config/libbid/bid64_to_bid128.c b/libgcc/config/libbid/bid64_to_bid128.c
index 6e55ba24e02..a8daddce6bc 100644
--- a/libgcc/config/libbid/bid64_to_bid128.c
+++ b/libgcc/config/libbid/bid64_to_bid128.c
@@ -153,9 +153,6 @@ bid128_to_bid64 (UINT128 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
T128 = round_const_table_128[rmode][extra_digits];
__add_carry_out (CX1.w[0], carry, T128.w[0], CX.w[0]);
CX1.w[1] = CX.w[1] + T128.w[1] + carry;
- if (__unsigned_compare_ge_128
- (CX1, power10_table_128[extra_digits + 16]))
- uf_check = 0;
}
extra_digits =
extra_digits + DECIMAL_EXPONENT_BIAS_128 -
diff --git a/libgcc/config/libbid/bid_binarydecimal.c b/libgcc/config/libbid/bid_binarydecimal.c
index 5b5b721ffb3..6df39f69887 100644
--- a/libgcc/config/libbid/bid_binarydecimal.c
+++ b/libgcc/config/libbid/bid_binarydecimal.c
@@ -566,19 +566,19 @@ BID_BINARY80LDOUBLE;
{ if ((x & (0xFull<<27)) == (0xFull<<27)) \
{ if ((x & (0x1Full<<26)) != (0x1Full<<26)) inf; \
if ((x & (1ul<<25))!=0) *pfpsf |= INVALID_EXCEPTION; \
- nan(s,((((x) & 0xFFFFul) > 999999ul) ? 0 : \
+ nan(s,((((x) & 0xFFFFFul) > 999999ul) ? 0 : \
(((unsigned long long) x) << 44)),0ull); \
} \
e = ((x >> 21) & ((1ull<<8)-1)) - 101; \
c = (1ull<<23) + (x & ((1ull<<21)-1)); \
- if ((unsigned long)(c) > 9999999ul) c = 0; \
+ if ((unsigned long)(c) > 9999999ul) zero; \
k = 0; \
} \
else \
{ e = ((x >> 23) & ((1ull<<8)-1)) - 101; \
c = x & ((1ull<<23)-1); \
if (c == 0) zero; \
- k = clz32(c) - 8; \
+ k = clz32_nz(c) - 8; \
c = c << k; \
} \
}
@@ -594,14 +594,14 @@ BID_BINARY80LDOUBLE;
} \
e = ((x >> 51) & ((1ull<<10)-1)) - 398; \
c = (1ull<<53) + (x & ((1ull<<51)-1)); \
- if ((unsigned long long)(c) > 9999999999999999ull) c = 0; \
+ if ((unsigned long long)(c) > 9999999999999999ull) zero; \
k = 0; \
} \
else \
{ e = ((x >> 53) & ((1ull<<10)-1)) - 398; \
c = x & ((1ull<<53)-1); \
if (c == 0) zero; \
- k = clz64(c) - 10; \
+ k = clz64_nz(c) - 10; \
c = c << k; \
} \
}
@@ -144302,20 +144302,6 @@ bid32_to_binary64 (UINT32 x
// We actually check if e >= ceil((sci_emax + 1) * log_10(2))
// which in this case is e >= ceil(1024 * log_10(2)) = ceil(308.25) = 309
- if (e >= 309) {
- *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
- return_binary64_ovf (s);
- }
-// Also check for "trivial" underflow, when 10^e * 2^113 <= 2^emin * 1/4,
-// so test e <= floor((emin - 115) * log_10(2))
-// In this case just fix ourselves at that value for uniformity.
-//
-// This is important not only to keep the tables small but to maintain the
-// testing of the round/sticky words as a correct rounding method
-
- if (e <= -358)
- e = -358;
-
// Look up the breakpoint and approximate exponent
m_min = (breakpoints_binary64 + 358)[e];
@@ -144323,7 +144309,7 @@ bid32_to_binary64 (UINT32 x
// Choose provisional exponent and reciprocal multiplier based on breakpoint
- if (le128 (c.w[1], c.w[0], m_min.w[1], m_min.w[0])) {
+ if (c.w[1] < m_min.w[1]) {
r = (multipliers1_binary64 + 358)[e];
} else {
r = (multipliers2_binary64 + 358)[e];
@@ -144332,17 +144318,12 @@ bid32_to_binary64 (UINT32 x
// Do the reciprocal multiplication
- __mul_128x256_to_384 (z, c, r)
+ __mul_64x256_to_320(z, c.w[1], r);
+ z.w[5]=z.w[4]; z.w[4]=z.w[3]; z.w[3]=z.w[2]; z.w[2]=z.w[1]; z.w[1]=z.w[0]; z.w[0]=0;
+
// Check for exponent underflow and compensate by shifting the product
// Cut off the process at precision+2, since we can't really shift further
- if (e_out < 1) {
- int d;
- d = 1 - e_out;
- if (d > 55)
- d = 55;
- e_out = 1;
- srl256 (z.w[5], z.w[4], z.w[3], z.w[2], d);
- }
+
c_prov = z.w[5];
// Round using round-sticky words
@@ -144353,31 +144334,14 @@ bid32_to_binary64 (UINT32 x
w[1],
roundbound_128[(rnd_mode << 2) + ((s & 1) << 1) +
(c_prov & 1)].w[0], z.w[4], z.w[3])) {
- c_prov = c_prov + 1;
- if (c_prov == (1ull << 53)) {
- c_prov = 1ull << 52;
- e_out = e_out + 1;
- }
- }
-// Check for overflow
-
- if (e_out >= 2047) {
- *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
- return_binary64_ovf (s);
+ c_prov = c_prov + 1;
}
-// Modify exponent for a tiny result, otherwise lop the implicit bit
-
- if (c_prov < (1ull << 52))
- e_out = 0;
- else
- c_prov = c_prov & ((1ull << 52) - 1);
+ c_prov = c_prov & ((1ull << 52) - 1);
// Set the inexact and underflow flag as appropriate
if ((z.w[4] != 0) || (z.w[3] != 0)) {
*pfpsf |= INEXACT_EXCEPTION;
- if (e_out == 0)
- *pfpsf |= UNDERFLOW_EXCEPTION;
}
// Package up the result as a binary floating-point number
@@ -145756,6 +145720,14 @@ binary64_to_bid32 (double x
__mul_128x256_to_384 (z, c, r)
c_prov = z.w[5];
+// Test inexactness and underflow (when testing tininess before rounding)
+
+ if ((z.w[4] != 0) || (z.w[3] != 0)) {
+ *pfpsf |= INEXACT_EXCEPTION;
+ if (c_prov < 1000000ull)
+ *pfpsf |= UNDERFLOW_EXCEPTION;
+ }
+
// Round using round-sticky words
// If we spill over into the next decade, correct
// Flag underflow where it may be needed even for |result| = SNN
@@ -145769,27 +145741,16 @@ binary64_to_bid32 (double x
if (c_prov == 10000000ull) {
c_prov = 1000000ull;
e_out = e_out + 1;
- } else if ((c_prov == 1000000ull) && (e_out == 0)) {
- if ((((rnd_mode & 3) == 0) && (z.w[4] <= 17524406870024074035ull))
- || ((rnd_mode + (s & 1) == 2)
- && (z.w[4] <= 16602069666338596454ull)))
- *pfpsf |= UNDERFLOW_EXCEPTION;
}
}
+
// Check for overflow
if (e_out > 90 + 101) {
*pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
return_bid32_ovf (s);
}
-// Set the inexact flag as appropriate and check underflow
-// It's no doubt superfluous to check inexactness, but anyway...
- if ((z.w[4] != 0) || (z.w[3] != 0)) {
- *pfpsf |= INEXACT_EXCEPTION;
- if (c_prov < 1000000ull)
- *pfpsf |= UNDERFLOW_EXCEPTION;
- }
// Package up the result
return_bid32 (s, e_out, c_prov);
@@ -145919,6 +145880,14 @@ binary80_to_bid32 (BINARY80 x
__mul_128x256_to_384 (z, c, r)
c_prov = z.w[5];
+// Test inexactness and underflow (when testing tininess before rounding)
+
+ if ((z.w[4] != 0) || (z.w[3] != 0)) {
+ *pfpsf |= INEXACT_EXCEPTION;
+ if (c_prov < 1000000ull)
+ *pfpsf |= UNDERFLOW_EXCEPTION;
+ }
+
// Round using round-sticky words
// If we spill over into the next decade, correct
// Flag underflow where it may be needed even for |result| = SNN
@@ -145932,27 +145901,16 @@ binary80_to_bid32 (BINARY80 x
if (c_prov == 10000000ull) {
c_prov = 1000000ull;
e_out = e_out + 1;
- } else if ((c_prov == 1000000ull) && (e_out == 0)) {
- if ((((rnd_mode & 3) == 0) && (z.w[4] <= 17524406870024074035ull))
- || ((rnd_mode + (s & 1) == 2)
- && (z.w[4] <= 16602069666338596454ull)))
- *pfpsf |= UNDERFLOW_EXCEPTION;
}
}
+
// Check for overflow
if (e_out > 90 + 101) {
*pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
return_bid32_ovf (s);
}
-// Set the inexact flag as appropriate and check underflow
-// It's no doubt superfluous to check inexactness, but anyway...
- if ((z.w[4] != 0) || (z.w[3] != 0)) {
- *pfpsf |= INEXACT_EXCEPTION;
- if (c_prov < 1000000ull)
- *pfpsf |= UNDERFLOW_EXCEPTION;
- }
// Package up the result
return_bid32 (s, e_out, c_prov);
@@ -146071,6 +146029,13 @@ binary128_to_bid32 (BINARY128 x
__mul_128x256_to_384 (z, c, r)
c_prov = z.w[5];
+// Test inexactness and underflow (when testing tininess before rounding)
+ if ((z.w[4] != 0) || (z.w[3] != 0)) {
+ *pfpsf |= INEXACT_EXCEPTION;
+ if (c_prov < 1000000ull)
+ *pfpsf |= UNDERFLOW_EXCEPTION;
+ }
+
// Round using round-sticky words
// If we spill over into the next decade, correct
// Flag underflow where it may be needed even for |result| = SNN
@@ -146086,30 +146051,16 @@ binary128_to_bid32 (BINARY128 x
if (c_prov == 10000000ull) {
c_prov = 1000000ull;
e_out = e_out + 1;
- } else if ((c_prov == 1000000ull) && (e_out == 0)) {
- if ((((rnd_mode & 3) == 0) &&
- le128 (z.w[4], z.w[3],
- 17524406870024074035ull, 3689348814741910323ull)) ||
- ((rnd_mode + (s & 1) == 2) &&
- le128 (z.w[4], z.w[3],
- 16602069666338596454ull, 7378697629483820646ull)))
- *pfpsf |= UNDERFLOW_EXCEPTION;
}
}
+
// Check for overflow
if (e_out > 90 + 101) {
*pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
return_bid32_ovf (s);
}
-// Set the inexact flag as appropriate and check underflow
-// It's no doubt superfluous to check inexactness, but anyway...
- if ((z.w[4] != 0) || (z.w[3] != 0)) {
- *pfpsf |= INEXACT_EXCEPTION;
- if (c_prov < 1000000ull)
- *pfpsf |= UNDERFLOW_EXCEPTION;
- }
// Package up the result
return_bid32 (s, e_out, c_prov);
@@ -146562,6 +146513,14 @@ binary80_to_bid64 (BINARY80 x
__mul_128x256_to_384 (z, c, r)
c_prov = z.w[5];
+// Test inexactness and underflow (when testing tininess before rounding)
+
+ if ((z.w[4] != 0) || (z.w[3] != 0)) {
+ *pfpsf |= INEXACT_EXCEPTION;
+ if (c_prov < 1000000000000000ull)
+ *pfpsf |= UNDERFLOW_EXCEPTION;
+ }
+
// Round using round-sticky words
// If we spill over into the next decade, correct
// Flag underflow where it may be needed even for |result| = SNN
@@ -146575,27 +146534,16 @@ binary80_to_bid64 (BINARY80 x
if (c_prov == 10000000000000000ull) {
c_prov = 1000000000000000ull;
e_out = e_out + 1;
- } else if ((c_prov == 1000000000000000ull) && (e_out == 0)) {
- if ((((rnd_mode & 3) == 0) && (z.w[4] <= 17524406870024074035ull))
- || ((rnd_mode + (s & 1) == 2)
- && (z.w[4] <= 16602069666338596454ull)))
- *pfpsf |= UNDERFLOW_EXCEPTION;
}
}
+
// Check for overflow
if (e_out > 369 + 398) {
*pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
return_bid64_ovf (s);
}
-// Set the inexact flag as appropriate and check underflow
-// It's no doubt superfluous to check inexactness, but anyway...
- if ((z.w[4] != 0) || (z.w[3] != 0)) {
- *pfpsf |= INEXACT_EXCEPTION;
- if (c_prov < 1000000000000000ull)
- *pfpsf |= UNDERFLOW_EXCEPTION;
- }
// Package up the result
return_bid64 (s, e_out, c_prov);
@@ -146723,6 +146671,14 @@ binary128_to_bid64 (BINARY128 x
__mul_128x256_to_384 (z, c, r)
c_prov = z.w[5];
+// Test inexactness and underflow (when testing tininess before rounding)
+
+ if ((z.w[4] != 0) || (z.w[3] != 0)) {
+ *pfpsf |= INEXACT_EXCEPTION;
+ if (c_prov < 1000000000000000ull)
+ *pfpsf |= UNDERFLOW_EXCEPTION;
+ }
+
// Round using round-sticky words
// If we spill over into the next decade, correct
// Flag underflow where it may be needed even for |result| = SNN
@@ -146736,27 +146692,16 @@ binary128_to_bid64 (BINARY128 x
if (c_prov == 10000000000000000ull) {
c_prov = 1000000000000000ull;
e_out = e_out + 1;
- } else if ((c_prov == 1000000000000000ull) && (e_out == 0)) {
- if ((((rnd_mode & 3) == 0) && (z.w[4] <= 17524406870024074035ull))
- || ((rnd_mode + (s & 1) == 2)
- && (z.w[4] <= 16602069666338596454ull)))
- *pfpsf |= UNDERFLOW_EXCEPTION;
}
}
+
// Check for overflow
if (e_out > 369 + 398) {
*pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
return_bid64_ovf (s);
}
-// Set the inexact flag as appropriate and check underflow
-// It's no doubt superfluous to check inexactness, but anyway...
- if ((z.w[4] != 0) || (z.w[3] != 0)) {
- *pfpsf |= INEXACT_EXCEPTION;
- if (c_prov < 1000000000000000ull)
- *pfpsf |= UNDERFLOW_EXCEPTION;
- }
// Package up the result
return_bid64 (s, e_out, c_prov);
diff --git a/libgcc/config/libbid/bid_conf.h b/libgcc/config/libbid/bid_conf.h
index dff938edd5d..e054a3ff570 100644
--- a/libgcc/config/libbid/bid_conf.h
+++ b/libgcc/config/libbid/bid_conf.h
@@ -519,6 +519,14 @@ see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
#define BID_BIG_ENDIAN __FLOAT_WORD_ORDER__ == __ORDER_BIG_ENDIAN__
#endif
+#if BID_BIG_ENDIAN
+#define BID_HIGH_128W 0
+#define BID_LOW_128W 1
+#else
+#define BID_HIGH_128W 1
+#define BID_LOW_128W 0
+#endif
+
#ifndef BID_THREAD
#if defined (HAVE_CC_TLS) && defined (USE_TLS)
#define BID_THREAD __thread
diff --git a/libgcc/config/libbid/bid_functions.h b/libgcc/config/libbid/bid_functions.h
index 9b4fbef8d01..cce87eb4d60 100644
--- a/libgcc/config/libbid/bid_functions.h
+++ b/libgcc/config/libbid/bid_functions.h
@@ -67,9 +67,13 @@ ALIGN (16)
#endif
+#if defined __NO_BINARY80__
+#define __ENABLE_BINARY80__ 0
+#else
#if !defined _MSC_VER || defined __INTEL_COMPILER
#define __ENABLE_BINARY80__ 1
#endif
+#endif
#ifndef HPUX_OS
#define BINARY80 long double
@@ -91,6 +95,19 @@ ALIGN (16)
} UINT256;
typedef unsigned int FPSC; // floating-point status and control
+ typedef enum class_types {
+ signalingNaN,
+ quietNaN,
+ negativeInfinity,
+ negativeNormal,
+ negativeSubnormal,
+ negativeZero,
+ positiveZero,
+ positiveSubnormal,
+ positiveNormal,
+ positiveInfinity
+ } class_t;
+
// TYPE parameters
#define BID128_MAXDIGITS 34
#define BID64_MAXDIGITS 16
@@ -2948,7 +2965,7 @@ ALIGN (16)
extern UINT64 bid64_copySign (UINT64 x,
UINT64 y _EXC_MASKS_PARAM
_EXC_INFO_PARAM);
- extern int bid64_class (UINT64 x _EXC_MASKS_PARAM _EXC_INFO_PARAM);
+ extern class_t bid64_class (UINT64 x _EXC_MASKS_PARAM _EXC_INFO_PARAM);
extern int bid64_sameQuantum (UINT64 x, UINT64 y
_EXC_MASKS_PARAM _EXC_INFO_PARAM);
extern int bid64_totalOrder (UINT64 x, UINT64 y
@@ -2984,8 +3001,8 @@ ALIGN (16)
extern UINT128 bid128_copySign (UINT128 x,
UINT128 y _EXC_MASKS_PARAM
_EXC_INFO_PARAM);
- extern int bid128_class (UINT128 x _EXC_MASKS_PARAM
- _EXC_INFO_PARAM);
+ extern class_t bid128_class (UINT128 x _EXC_MASKS_PARAM
+ _EXC_INFO_PARAM);
extern int bid128_sameQuantum (UINT128 x,
UINT128 y _EXC_MASKS_PARAM
_EXC_INFO_PARAM);
diff --git a/libgcc/config/libbid/bid_inline_add.h b/libgcc/config/libbid/bid_inline_add.h
index eb16ea8e0ba..9ca41a2f89e 100644
--- a/libgcc/config/libbid/bid_inline_add.h
+++ b/libgcc/config/libbid/bid_inline_add.h
@@ -918,6 +918,7 @@ get_add128 (UINT64 sign_x, int exponent_x, UINT64 coefficient_x,
coefficient_x += D;
}
break;
+ default: break; // default added to avoid compiler warning
}
if (coefficient_x < 1000000000000000ull) {
coefficient_x -= D;
@@ -1107,6 +1108,7 @@ get_add128 (UINT64 sign_x, int exponent_x, UINT64 coefficient_x,
} else if (FS.w[1] | FS.w[0])
CYh++;
break;
+ default: break; // default added to avoid compiler warning
}
#endif
#endif
diff --git a/libgcc/config/libbid/bid_internal.h b/libgcc/config/libbid/bid_internal.h
index 764abf81057..d3c83b5c219 100644
--- a/libgcc/config/libbid/bid_internal.h
+++ b/libgcc/config/libbid/bid_internal.h
@@ -970,6 +970,8 @@ get_BID64 (UINT64 sgn, int expon, UINT64 coeff, int rmode,
// round up
if (sgn)
r = SMALLEST_BID64;
+ default:
+ break;
}
return r;
}
@@ -1086,6 +1088,8 @@ fast_get_BID64_check_OF (UINT64 sgn, int expon, UINT64 coeff, int rmode,
// round up
if (sgn)
r = SMALLEST_BID64;
+ default:
+ break;
}
return r;
}
@@ -2582,19 +2586,6 @@ ALIGN (16)
A=((tempx.i >>23) & EXPONENT_MASK32) - 0x7f;\
}
- enum class_types {
- signalingNaN,
- quietNaN,
- negativeInfinity,
- negativeNormal,
- negativeSubnormal,
- negativeZero,
- positiveZero,
- positiveSubnormal,
- positiveNormal,
- positiveInfinity
- };
-
typedef union {
UINT64 ui64;
double d;
--
2.31.1
^ permalink raw reply [flat|nested] 2+ messages in thread
* Re: [PATCH] Update libbid according to the latest Intel Decimal Floating-Point Math Library.
2024-04-28 5:50 [PATCH] Update libbid according to the latest Intel Decimal Floating-Point Math Library liuhongt
@ 2024-05-02 9:20 ` Richard Biener
0 siblings, 0 replies; 2+ messages in thread
From: Richard Biener @ 2024-05-02 9:20 UTC (permalink / raw)
To: liuhongt; +Cc: gcc-patches, crazylht, hjl.tools
On Sun, Apr 28, 2024 at 7:53 AM liuhongt <hongtao.liu@intel.com> wrote:
>
> The Intel Decimal Floating-Point Math Library is available as open-source on Netlib[1].
>
> [1] https://www.netlib.org/misc/intel/.
>
> Bootstrapped and regtested on x86_64-pc-linux-gnu{-m32,}.
> Ok for trunk?
OK for trunk.
Thanks,
Richard.
> libgcc/config/libbid/ChangeLog:
>
> * bid128_fma.c (add_and_round): Fix bug: the result
> of (+5E+368)*(+10E-34)+(-10E+369) was returning
> -9999999999999999999999999999999999E+336 instead of expected
> result -1000000000000000000000000000000000E+337.
> (bid128_ext_fma): Ditto.
> (bid64qqq_fma): Ditto.
> * bid128_noncomp.c: Change return type of bid128_class from
> int to class_t.
> * bid128_round_integral.c: Add default case to avoid compiler
> warning.
> * bid128_string.c (bid128_to_string): Replace 0x30 with '0'
> for zero digit.
> (bid128_from_string): Ditto.
> * bid32_to_bid128.c (bid128_to_bid32): Fix Bug. In addition
> to the INEXACT flag, the UNDERFLOW flag needs to be set (and
> was not) when converting an input such as
> +6931674235302037148946035460357709E+1857 to +1000000E-101
> * bid32_to_bid64.c (bid64_to_bid32): fix Bug, In addition to
> the INEXACT flag, the UNDERFLOW flag needs to be set (and was
> not) when converting an input such as +9999999000000001E-111
> to +1000000E-101. Furthermore, significant bits of NaNs are
> set correctly now. For example, 0x7c00003b9aca0000 was
> returning 0x7c000002 instead of 0x 7c000100.
> * bid64_noncomp.c: Change return type of bid64_class from int
> to class_t.
> * bid64_round_integral.c (bid64_round_integral_exact): Add
> default case to avoid compiler warning.
> * bid64_string.c (bid64_from_string): Fix bug for rounding
> up. The input string "10000000000000000" was returning
> +1000000000000001E+1 instead of +1000000000000000E+1.
> * bid64_to_bid128.c (bid128_to_bid64): Fix bug, in addition to
> the INEXACT flag, the UNDERFLOW flag needs to be set (and was
> not) when converting an input such as
> +9999999999999999999999999999999999E-417 to
> +1000000000000000E-398.
> * bid_binarydecimal.c (bid32_to_binary64): Fix bug for
> conversion between binary and bid types. For example,
> 0x7c0F4240 was returning 0x7FFFA12000000000 instead of
> expected double precision 0x7FF8000000000000.
> (binary64_to_bid32): Ditto.
> (binary80_to_bid32): Ditto.
> (binary128_to_bid32): Ditto.
> (binary80_to_bid64): Ditto.
> (binary128_to_bid64): Ditto.
> * bid_conf.h (BID_HIGH_128W): New macro.
> (BID_LOW_128W): Ditto.
> * bid_functions.h (__ENABLE_BINARY80__): Ditto.
> (ALIGN): Ditto.
> * bid_inline_add.h (get_add128): Add default case to avoid compiler
> warning.
> * bid_internal.h (get_BID64): Ditto.
> (fast_get_BID64_check_OF): Ditto.
> (ALIGN): New macro.
>
> Co-authored-by: Anderson, Cristina S <cristina.s.anderson@intel.com>
> Co-authored-by: Akkas, Ahmet <ahmet.akkas@intel.com>
> Co-authored-by: Cornea, Marius <marius.cornea@intel.com>
> ---
> libgcc/config/libbid/bid128_fma.c | 188 ++++++++++---------
> libgcc/config/libbid/bid128_noncomp.c | 2 +-
> libgcc/config/libbid/bid128_round_integral.c | 2 +
> libgcc/config/libbid/bid128_string.c | 7 +-
> libgcc/config/libbid/bid32_to_bid128.c | 3 -
> libgcc/config/libbid/bid32_to_bid64.c | 11 +-
> libgcc/config/libbid/bid64_noncomp.c | 2 +-
> libgcc/config/libbid/bid64_round_integral.c | 2 +
> libgcc/config/libbid/bid64_string.c | 21 ++-
> libgcc/config/libbid/bid64_to_bid128.c | 3 -
> libgcc/config/libbid/bid_binarydecimal.c | 167 ++++++----------
> libgcc/config/libbid/bid_conf.h | 8 +
> libgcc/config/libbid/bid_functions.h | 23 ++-
> libgcc/config/libbid/bid_inline_add.h | 2 +
> libgcc/config/libbid/bid_internal.h | 17 +-
> 15 files changed, 220 insertions(+), 238 deletions(-)
>
> diff --git a/libgcc/config/libbid/bid128_fma.c b/libgcc/config/libbid/bid128_fma.c
> index 67233193a42..cbcf225546f 100644
> --- a/libgcc/config/libbid/bid128_fma.c
> +++ b/libgcc/config/libbid/bid128_fma.c
> @@ -417,13 +417,12 @@ add_and_round (int q3,
> R128.w[1] = R256.w[1];
> R128.w[0] = R256.w[0];
> }
> + if (e4 + x0 < expmin) { // for all rounding modes
> + is_tiny = 1;
> + }
> // the rounded result has p34 = 34 digits
> e4 = e4 + x0 + incr_exp;
> - if (rnd_mode == ROUNDING_TO_NEAREST) {
> - if (e4 < expmin) {
> - is_tiny = 1; // for other rounding modes apply correction
> - }
> - } else {
> + if (rnd_mode != ROUNDING_TO_NEAREST) {
> // for RM, RP, RZ, RA apply correction in order to determine tininess
> // but do not save the result; apply the correction to
> // (-1)^p_sign * significand * 10^0
> @@ -434,10 +433,6 @@ add_and_round (int q3,
> is_inexact_gt_midpoint, is_midpoint_lt_even,
> is_midpoint_gt_even, 0, &P128, ptrfpsf);
> scale = ((P128.w[1] & MASK_EXP) >> 49) - 6176; // -1, 0, or +1
> - // the number of digits in the significand is p34 = 34
> - if (e4 + scale < expmin) {
> - is_tiny = 1;
> - }
> }
> ind = p34; // the number of decimal digits in the signifcand of res
> res.w[1] = p_sign | ((UINT64) (e4 + 6176) << 49) | R128.w[1]; // RN
> @@ -851,7 +846,6 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> }
> }
> }
> -
> p_sign = x_sign ^ y_sign; // sign of the product
>
> // identify cases where at least one operand is infinity
> @@ -988,15 +982,10 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> if (C1.w[1] == 0) {
> if (C1.w[0] >= 0x0020000000000000ull) { // x >= 2^53
> // split the 64-bit value in two 32-bit halves to avoid rounding errors
> - if (C1.w[0] >= 0x0000000100000000ull) { // x >= 2^32
> +
> tmp.d = (double) (C1.w[0] >> 32); // exact conversion
> x_nr_bits =
> 33 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> - } else { // x < 2^32
> - tmp.d = (double) (C1.w[0]); // exact conversion
> - x_nr_bits =
> - 1 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> - }
> } else { // if x < 2^53
> tmp.d = (double) C1.w[0]; // exact conversion
> x_nr_bits =
> @@ -1011,42 +1000,36 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> if (q1 == 0) {
> q1 = nr_digits[x_nr_bits - 1].digits1;
> if (C1.w[1] > nr_digits[x_nr_bits - 1].threshold_hi ||
> - (C1.w[1] == nr_digits[x_nr_bits - 1].threshold_hi &&
> - C1.w[0] >= nr_digits[x_nr_bits - 1].threshold_lo))
> + (C1.w[1] == nr_digits[x_nr_bits - 1].threshold_hi &&
> + C1.w[0] >= nr_digits[x_nr_bits - 1].threshold_lo))
> q1++;
> }
> }
> -
> + // q2 = nr. of decimal digits in y
> + // determine first the nr. of bits in y
> if (C2.w[1] != 0 || C2.w[0] != 0) { // y = f (non-zero finite)
> if (C2.w[1] == 0) {
> if (C2.w[0] >= 0x0020000000000000ull) { // y >= 2^53
> // split the 64-bit value in two 32-bit halves to avoid rounding errors
> - if (C2.w[0] >= 0x0000000100000000ull) { // y >= 2^32
> tmp.d = (double) (C2.w[0] >> 32); // exact conversion
> y_nr_bits =
> - 32 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> - } else { // y < 2^32
> - tmp.d = (double) C2.w[0]; // exact conversion
> - y_nr_bits =
> - ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> - }
> + 33 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> } else { // if y < 2^53
> tmp.d = (double) C2.w[0]; // exact conversion
> y_nr_bits =
> - ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> + 1 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> }
> } else { // C2.w[1] != 0 => nr. bits = 64 + nr_bits (C2.w[1])
> tmp.d = (double) C2.w[1]; // exact conversion
> y_nr_bits =
> - 64 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> + 65 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> }
> -
> - q2 = nr_digits[y_nr_bits].digits;
> + q2 = nr_digits[y_nr_bits - 1].digits;
> if (q2 == 0) {
> - q2 = nr_digits[y_nr_bits].digits1;
> - if (C2.w[1] > nr_digits[y_nr_bits].threshold_hi ||
> - (C2.w[1] == nr_digits[y_nr_bits].threshold_hi &&
> - C2.w[0] >= nr_digits[y_nr_bits].threshold_lo))
> + q2 = nr_digits[y_nr_bits - 1].digits1;
> + if (C2.w[1] > nr_digits[y_nr_bits - 1].threshold_hi ||
> + (C2.w[1] == nr_digits[y_nr_bits - 1].threshold_hi &&
> + C2.w[0] >= nr_digits[y_nr_bits - 1].threshold_lo))
> q2++;
> }
> }
> @@ -1055,32 +1038,25 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> if (C3.w[1] == 0) {
> if (C3.w[0] >= 0x0020000000000000ull) { // z >= 2^53
> // split the 64-bit value in two 32-bit halves to avoid rounding errors
> - if (C3.w[0] >= 0x0000000100000000ull) { // z >= 2^32
> tmp.d = (double) (C3.w[0] >> 32); // exact conversion
> z_nr_bits =
> - 32 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> - } else { // z < 2^32
> - tmp.d = (double) C3.w[0]; // exact conversion
> - z_nr_bits =
> - ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> - }
> + 33 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> } else { // if z < 2^53
> tmp.d = (double) C3.w[0]; // exact conversion
> z_nr_bits =
> - ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> + 1 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> }
> } else { // C3.w[1] != 0 => nr. bits = 64 + nr_bits (C3.w[1])
> tmp.d = (double) C3.w[1]; // exact conversion
> z_nr_bits =
> - 64 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> + 65 + ((((unsigned int) (tmp.ui64 >> 52)) & 0x7ff) - 0x3ff);
> }
> -
> - q3 = nr_digits[z_nr_bits].digits;
> + q3 = nr_digits[z_nr_bits - 1].digits;
> if (q3 == 0) {
> - q3 = nr_digits[z_nr_bits].digits1;
> - if (C3.w[1] > nr_digits[z_nr_bits].threshold_hi ||
> - (C3.w[1] == nr_digits[z_nr_bits].threshold_hi &&
> - C3.w[0] >= nr_digits[z_nr_bits].threshold_lo))
> + q3 = nr_digits[z_nr_bits - 1].digits1;
> + if (C3.w[1] > nr_digits[z_nr_bits - 1].threshold_hi ||
> + (C3.w[1] == nr_digits[z_nr_bits - 1].threshold_hi &&
> + C3.w[0] >= nr_digits[z_nr_bits - 1].threshold_lo))
> q3++;
> }
> }
> @@ -1128,7 +1104,6 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> } else {
> ; // continue with x = f, y = f, z = 0 or x = f, y = f, z = f
> }
> -
> e1 = (x_exp >> 49) - 6176; // unbiased exponent of x
> e2 = (y_exp >> 49) - 6176; // unbiased exponent of y
> e3 = (z_exp >> 49) - 6176; // unbiased exponent of z
> @@ -1232,22 +1207,18 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> // length of C1 * C2 rounded up to a multiple of 64 bits is len = 192;
> q4 = q1 + q2; // q4 in [40, 57]
> }
> - } else if (q1 + q2 == 58) { // C4 = C1 * C2 fits in 192 or 256 bits
> - // both C1 and C2 fit in 128 bits (actually in 113 bits); at most one
> - // may fit in 64 bits
> - if (C1.w[1] == 0) { // C1 * C2 will fit in 192 bits
> - __mul_64x128_full (C4.w[2], C4, C1.w[0], C2); // may use 64x128_to_192
> - } else if (C2.w[1] == 0) { // C1 * C2 will fit in 192 bits
> - __mul_64x128_full (C4.w[2], C4, C2.w[0], C1); // may use 64x128_to_192
> - } else { // C1 * C2 will fit in 192 bits or in 256 bits
> - __mul_128x128_to_256 (C4, C1, C2);
> - }
> + } else if (q1 + q2 == 58) { // C4 = C1 * C2 fits in 192 or 256 bits;
> + // both C1 and C2 fit in 128 bits (actually in 113 bits); none can
> + // fit in 64 bits, because each number must have at least 24 decimal
> + // digits for the sum to have 58 (as the max. nr. of digits is 34) =>
> + // C1.w[1] != 0 and C2.w[1] != 0
> + __mul_128x128_to_256 (C4, C1, C2);
> // if C4 < 10^(q1+q2-1) = 10^57 then q4 = q1+q2-1 = 57 else q4 = q1+q2 = 58
> if (C4.w[3] == 0 && (C4.w[2] < ten2k256[18].w[2] ||
> - (C4.w[2] == ten2k256[18].w[2]
> - && (C4.w[1] < ten2k256[18].w[1]
> - || (C4.w[1] == ten2k256[18].w[1]
> - && C4.w[0] < ten2k256[18].w[0]))))) {
> + (C4.w[2] == ten2k256[18].w[2]
> + && (C4.w[1] < ten2k256[18].w[1]
> + || (C4.w[1] == ten2k256[18].w[1]
> + && C4.w[0] < ten2k256[18].w[0]))))) {
> // 18 = 57 - 39 = q1+q2-1 - 39
> // length of C1 * C2 rounded up to a multiple of 64 bits is len = 192;
> q4 = 57; // 57 = q1 + q2 - 1
> @@ -1283,7 +1254,6 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> q4 = q1 + q2; // q4 in [59, 68]
> }
> }
> -
> if (C3.w[1] == 0x0 && C3.w[0] == 0x0) { // x = f, y = f, z = 0
> save_fpsf = *pfpsf; // sticky bits - caller value must be preserved
> *pfpsf = 0;
> @@ -1319,10 +1289,11 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> res.w[1] = R256.w[1];
> }
> e4 = e4 + x0;
> + q4 = p34;
> if (incr_exp) {
> e4 = e4 + 1;
> + if (q4 + e4 == expmin + p34) *pfpsf |= (INEXACT_EXCEPTION | UNDERFLOW_EXCEPTION);
> }
> - q4 = p34;
> // res is now the coefficient of the result rounded to the destination
> // precision, with unbounded exponent; the exponent is e4; q4=digits(res)
> } else { // if (q4 <= p34)
> @@ -1648,7 +1619,6 @@ bid128_ext_fma (int *ptr_is_midpoint_lt_even,
> delta = q3 + e3 - q4 - e4;
> delta_ge_zero:
> if (delta >= 0) {
> -
> if (p34 <= delta - 1 || // Case (1')
> (p34 == delta && e3 + 6176 < p34 - q3)) { // Case (1''A)
> // check for overflow, which can occur only in Case (1')
> @@ -1736,7 +1706,7 @@ delta_ge_zero:
> res.w[1] = z_sign | ((UINT64) (e3 + 6176) << 49) | C3.w[1];
> res.w[0] = C3.w[0];
> }
> -
> +
> // use the following to avoid double rounding errors when operating on
> // mixed formats in rounding to nearest, and for correcting the result
> // if not rounding to nearest
> @@ -1795,7 +1765,10 @@ delta_ge_zero:
> R64 = 10;
> }
> }
> - if (q4 == 1 && C4.w[0] == 5) {
> +
> + if (R64 == 5 && !is_inexact_lt_midpoint && !is_inexact_gt_midpoint &&
> + !is_midpoint_lt_even && !is_midpoint_gt_even) {
> + //if (q4 == 1 && C4.w[0] == 5) {
> is_inexact_lt_midpoint = 0;
> is_inexact_gt_midpoint = 0;
> is_midpoint_lt_even = 1;
> @@ -1826,11 +1799,7 @@ delta_ge_zero:
> res.w[1] = z_sign | ((UINT64) (e3 + 6176) << 49) | res.w[1];
> }
> if (e3 == expmin) {
> - if (R64 < 5 || (R64 == 5 && !is_inexact_lt_midpoint)) {
> - ; // result not tiny (in round-to-nearest mode)
> - } else {
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> - }
> + *pfpsf |= UNDERFLOW_EXCEPTION; // tiny if detected before rounding
> }
> } // end 10^(q3+scale-1)
> // set the inexact flag
> @@ -1877,10 +1846,9 @@ delta_ge_zero:
> // endif
> if ((e3 == expmin && (q3 + scale) < p34) ||
> (e3 == expmin && (q3 + scale) == p34 &&
> - (res.w[1] & MASK_COEFF) == 0x0000314dc6448d93ull && // 10^33_high
> - res.w[0] == 0x38c15b0a00000000ull && // 10^33_low
> - z_sign != p_sign && ((!z_sign && rnd_mode != ROUNDING_UP) ||
> - (z_sign && rnd_mode != ROUNDING_DOWN)))) {
> + (res.w[1] & MASK_COEFF) == 0x0000314dc6448d93ull && // 10^33_high
> + res.w[0] == 0x38c15b0a00000000ull && // 10^33_low
> + z_sign != p_sign)) {
> *pfpsf |= UNDERFLOW_EXCEPTION;
> }
> if (rnd_mode != ROUNDING_TO_NEAREST) {
> @@ -2594,7 +2562,7 @@ delta_ge_zero:
> if (e3 > expmin && ((res.w[1] < 0x0000314dc6448d93ull ||
> (res.w[1] == 0x0000314dc6448d93ull &&
> res.w[0] < 0x38c15b0a00000000ull)) ||
> - (is_inexact_lt_midpoint
> + ((is_inexact_lt_midpoint | is_midpoint_gt_even)
> && res.w[1] == 0x0000314dc6448d93ull
> && res.w[0] == 0x38c15b0a00000000ull))
> && x0 >= 1) {
> @@ -2678,6 +2646,9 @@ delta_ge_zero:
> res.w[0] < 0x38c15b0a00000000ull)) {
> is_tiny = 1;
> }
> + if (((res.w[1] & 0x7fffffffffffffffull) == 0x0000314dc6448d93ull) &&
> + (res.w[0] == 0x38c15b0a00000000ull) && // 10^33*10^-6176
> + (z_sign != p_sign)) is_tiny = 1;
> } else if (e3 < expmin) {
> // the result is tiny, so we must truncate more of res
> is_tiny = 1;
> @@ -3328,9 +3299,6 @@ delta_ge_zero:
> 0, &P128, pfpsf);
> scale = ((P128.w[1] & MASK_EXP) >> 49) - 6176; // -1, 0, or +1
> // the number of digits in the significand is p34 = 34
> - if (e4 + scale < expmin) {
> - is_tiny = 1;
> - }
> }
>
> // the result rounded to the destination precision with unbounded exponent
> @@ -3521,6 +3489,19 @@ delta_ge_zero:
> is_midpoint_lt_even, is_midpoint_gt_even,
> e4, &res, pfpsf);
> }
> + // correction needed for tininess detection before rounding
> + if ((((res.w[1] & 0x7fffffffffffffffull) == 0x0000314dc6448d93ull) &&
> + // 10^33*10^-6176_high
> + (res.w[0] == 0x38c15b0a00000000ull)) && // 10^33*10^-6176_low
> + (((rnd_mode == ROUNDING_TO_NEAREST ||
> + rnd_mode == ROUNDING_TIES_AWAY) &&
> + (is_midpoint_lt_even || is_inexact_gt_midpoint)) ||
> + ((((rnd_mode == ROUNDING_UP) && !(res.w[1] & MASK_SIGN)) ||
> + ((rnd_mode == ROUNDING_DOWN) && (res.w[1] & MASK_SIGN)))
> + && (is_midpoint_lt_even || is_midpoint_gt_even ||
> + is_inexact_lt_midpoint || is_inexact_gt_midpoint)))) {
> + is_tiny = 1;
> + }
> if (is_midpoint_lt_even || is_midpoint_gt_even ||
> is_inexact_lt_midpoint || is_inexact_gt_midpoint) {
> // set the inexact flag
> @@ -4162,21 +4143,34 @@ bid64qqq_fma (UINT128 x, UINT128 y, UINT128 z
> // determine the unbiased exponent of the result
> unbexp = ((res1 >> 53) & 0x3ff) - 398;
>
> + if (!((res1 & MASK_NAN) == MASK_NAN)) { // res1 not NaN
> // if subnormal, res1 must have exp = -398
> // if tiny and inexact set underflow and inexact status flags
> - if (!((res1 & MASK_NAN) == MASK_NAN) && // res1 not NaN
> - (unbexp == -398)
> - && ((res1 & MASK_BINARY_SIG1) < 1000000000000000ull)
> - && (is_inexact_lt_midpoint0 || is_inexact_gt_midpoint0
> - || is_midpoint_lt_even0 || is_midpoint_gt_even0)) {
> - // set the inexact flag and the underflow flag
> - *pfpsf |= (INEXACT_EXCEPTION | UNDERFLOW_EXCEPTION);
> + if ((unbexp == -398)
> + && ((res1 & MASK_BINARY_SIG1) < 1000000000000000ull)
> + && (is_inexact_lt_midpoint0 || is_inexact_gt_midpoint0
> + || is_midpoint_lt_even0 || is_midpoint_gt_even0)) {
> + // set the inexact flag and the underflow flag
> + *pfpsf |= (INEXACT_EXCEPTION | UNDERFLOW_EXCEPTION);
> } else if (is_inexact_lt_midpoint0 || is_inexact_gt_midpoint0 ||
> is_midpoint_lt_even0 || is_midpoint_gt_even0) {
> // set the inexact flag and the underflow flag
> *pfpsf |= INEXACT_EXCEPTION;
> - }
> + }
>
> + // correction needed for tininess detection before rounding
> + if (((res1 & 0x7fffffffffffffffull) == 1000000000000000ull) &&
> + // 10^15*10^-398
> + (((rnd_mode == ROUNDING_TO_NEAREST ||
> + rnd_mode == ROUNDING_TIES_AWAY) &&
> + (is_midpoint_lt_even || is_inexact_gt_midpoint)) ||
> + ((((rnd_mode == ROUNDING_UP) && !(res1 & MASK_SIGN)) ||
> + ((rnd_mode == ROUNDING_DOWN) && (res1 & MASK_SIGN)))
> + && (is_midpoint_lt_even || is_midpoint_gt_even ||
> + is_inexact_lt_midpoint || is_inexact_gt_midpoint)))) {
> + *pfpsf |= UNDERFLOW_EXCEPTION;
> + }
> + }
> *pfpsf |= save_fpsf;
> BID_RETURN (res1);
> } // else continue, and use rounding to nearest to round to 16 digits
> @@ -4453,6 +4447,20 @@ bid64qqq_fma (UINT128 x, UINT128 y, UINT128 z
> res1 = sign | MASK_STEERING_BITS |
> ((UINT64) (unbexp + 398) << 51) | (res1 & MASK_BINARY_SIG2);
> }
> +
> + // correction needed for tininess detection before rounding
> + if (((res1 & 0x7fffffffffffffffull) == 1000000000000000ull) &&
> + // 10^15*10^-398
> + (((rnd_mode == ROUNDING_TO_NEAREST ||
> + rnd_mode == ROUNDING_TIES_AWAY) &&
> + (is_midpoint_lt_even || is_inexact_gt_midpoint)) ||
> + ((((rnd_mode == ROUNDING_UP) && !(res1 & MASK_SIGN)) ||
> + ((rnd_mode == ROUNDING_DOWN) && (res1 & MASK_SIGN)))
> + && (is_midpoint_lt_even || is_midpoint_gt_even ||
> + is_inexact_lt_midpoint || is_inexact_gt_midpoint)))) {
> + *pfpsf |= UNDERFLOW_EXCEPTION;
> + }
> +
> *pfpsf |= save_fpsf;
> BID_RETURN (res1);
> }
> diff --git a/libgcc/config/libbid/bid128_noncomp.c b/libgcc/config/libbid/bid128_noncomp.c
> index a79ac859ce1..4ef166c81dc 100644
> --- a/libgcc/config/libbid/bid128_noncomp.c
> +++ b/libgcc/config/libbid/bid128_noncomp.c
> @@ -443,7 +443,7 @@ void
> bid128_class (int *pres, UINT128 * px _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
> UINT128 x = *px;
> #else
> -int
> +class_t
> bid128_class (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
> #endif
> int res;
> diff --git a/libgcc/config/libbid/bid128_round_integral.c b/libgcc/config/libbid/bid128_round_integral.c
> index ddaa0b818f2..051f8268aa9 100644
> --- a/libgcc/config/libbid/bid128_round_integral.c
> +++ b/libgcc/config/libbid/bid128_round_integral.c
> @@ -177,6 +177,7 @@ case ROUNDING_TO_ZERO:
> BID_RETURN (res);
> }
> break;
> +default: break; // default added to avoid compiler warning
> }
>
> // q = nr. of decimal digits in x
> @@ -804,6 +805,7 @@ case ROUNDING_TO_ZERO:
> BID_RETURN (res);
> }
> break;
> +default: break; // default added to avoid compiler warning
> }
>
> BID_RETURN (res);
> diff --git a/libgcc/config/libbid/bid128_string.c b/libgcc/config/libbid/bid128_string.c
> index ecd295cfbdf..8fc12ee2d76 100644
> --- a/libgcc/config/libbid/bid128_string.c
> +++ b/libgcc/config/libbid/bid128_string.c
> @@ -56,6 +56,7 @@ bid128_to_string (char *str, UINT128 x
> UINT128 C1;
> unsigned int k = 0; // pointer in the string
> unsigned int d0, d123;
> + unsigned int zero_digit = (unsigned int) '0';
> UINT64 HI_18Dig, LO_18Dig, Tmp;
> UINT32 MiDi[12], *ptr;
> char *c_ptr_start, *c_ptr;
> @@ -232,14 +233,14 @@ bid128_to_string (char *str, UINT128 x
> d123 = exp - 1000 * d0;
>
> if (d0) { // 1000 <= exp <= 6144 => 4 digits to return
> - str[k++] = d0 + 0x30;// ASCII for decimal digit d0
> + str[k++] = d0 + zero_digit; // ASCII for decimal digit d0
> ind = 3 * d123;
> str[k++] = char_table3[ind];
> str[k++] = char_table3[ind + 1];
> str[k++] = char_table3[ind + 2];
> } else { // 0 <= exp <= 999 => d0 = 0
> if (d123 < 10) { // 0 <= exp <= 9 => 1 digit to return
> - str[k++] = d123 + 0x30;// ASCII
> + str[k++] = d123 + zero_digit; // ASCII
> } else if (d123 < 100) { // 10 <= exp <= 99 => 2 digits to return
> ind = 2 * (d123 - 10);
> str[k++] = char_table2[ind];
> @@ -643,7 +644,7 @@ bid128_from_string (char *ps _RND_MODE_PARAM _EXC_FLAGS_PARAM
> }
> break;
>
> -
> + default: break; // default added to avoid compiler warning
> }
> // now form the coefficient as coeff_high*10^17+coeff_low+carry
> scale_high = 100000000000000000ull;
> diff --git a/libgcc/config/libbid/bid32_to_bid128.c b/libgcc/config/libbid/bid32_to_bid128.c
> index d1d1d3458fd..5b5ce9504e1 100644
> --- a/libgcc/config/libbid/bid32_to_bid128.c
> +++ b/libgcc/config/libbid/bid32_to_bid128.c
> @@ -155,9 +155,6 @@ bid128_to_bid32 (UINT128 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
> T128 = round_const_table_128[rmode][extra_digits];
> __add_carry_out (CX1.w[0], carry, T128.w[0], CX.w[0]);
> CX1.w[1] = CX.w[1] + T128.w[1] + carry;
> - if (__unsigned_compare_ge_128
> - (CX1, power10_table_128[extra_digits + 7]))
> - uf_check = 0;
> }
> extra_digits =
> extra_digits + DECIMAL_EXPONENT_BIAS_128 -
> diff --git a/libgcc/config/libbid/bid32_to_bid64.c b/libgcc/config/libbid/bid32_to_bid64.c
> index 7802346a3d1..61b24b29915 100644
> --- a/libgcc/config/libbid/bid32_to_bid64.c
> +++ b/libgcc/config/libbid/bid32_to_bid64.c
> @@ -79,6 +79,7 @@ bid64_to_bid32 (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
> UINT128 Q;
> UINT64 sign_x, coefficient_x, remainder_h, carry, Stemp;
> UINT32 res;
> + UINT64 t64;
> int_float tempx;
> int exponent_x, bin_expon_cx, extra_digits, rmode = 0, amount;
> unsigned status = 0;
> @@ -93,8 +94,10 @@ bid64_to_bid32 (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
> // unpack arguments, check for NaN or Infinity, 0
> if (!unpack_BID64 (&sign_x, &exponent_x, &coefficient_x, x)) {
> if (((x) & 0x7800000000000000ull) == 0x7800000000000000ull) {
> - res = (coefficient_x & 0x0003ffffffffffffull);
> - res /= 1000000000ull;
> + t64 = (coefficient_x & 0x0003ffffffffffffull);
> + res = t64/1000000000ull;
> + //res = (coefficient_x & 0x0003ffffffffffffull);
> + //res /= 1000000000ull;
> res |= ((coefficient_x >> 32) & 0xfc000000);
> #ifdef SET_STATUS_FLAGS
> if ((x & SNAN_MASK64) == SNAN_MASK64) // sNaN
> @@ -139,10 +142,6 @@ bid64_to_bid32 (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
> exponent_x += extra_digits;
> if ((exponent_x < 0) && (exponent_x + MAX_FORMAT_DIGITS_32 >= 0)) {
> status = UNDERFLOW_EXCEPTION;
> - if (exponent_x == -1)
> - if (coefficient_x + round_const_table[rmode][extra_digits] >=
> - power10_table_128[extra_digits + 7].w[0])
> - status = 0;
> extra_digits -= exponent_x;
> exponent_x = 0;
> }
> diff --git a/libgcc/config/libbid/bid64_noncomp.c b/libgcc/config/libbid/bid64_noncomp.c
> index ecc9dddfd8b..ec633693f62 100644
> --- a/libgcc/config/libbid/bid64_noncomp.c
> +++ b/libgcc/config/libbid/bid64_noncomp.c
> @@ -358,7 +358,7 @@ void
> bid64_class (int *pres, UINT64 * px _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
> UINT64 x = *px;
> #else
> -int
> +class_t
> bid64_class (UINT64 x _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
> #endif
> int res;
> diff --git a/libgcc/config/libbid/bid64_round_integral.c b/libgcc/config/libbid/bid64_round_integral.c
> index 99a3e6101ef..5a33c429116 100644
> --- a/libgcc/config/libbid/bid64_round_integral.c
> +++ b/libgcc/config/libbid/bid64_round_integral.c
> @@ -142,6 +142,7 @@ bid64_round_integral_exact (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
> BID_RETURN (res);
> }
> break;
> + default: break; // default added to avoid compiler warning
> } // end switch ()
>
> // q = nr. of decimal digits in x (1 <= q <= 54)
> @@ -483,6 +484,7 @@ bid64_round_integral_exact (UINT64 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
> BID_RETURN (res);
> }
> break;
> + default: break; // default added to avoid compiler warning
> } // end switch ()
> BID_RETURN (res);
> }
> diff --git a/libgcc/config/libbid/bid64_string.c b/libgcc/config/libbid/bid64_string.c
> index 55fcd9e9208..81ac5e275a4 100644
> --- a/libgcc/config/libbid/bid64_string.c
> +++ b/libgcc/config/libbid/bid64_string.c
> @@ -251,7 +251,7 @@ bid64_from_string (char *ps
> #endif
> UINT64 sign_x, coefficient_x = 0, rounded = 0, res;
> int expon_x = 0, sgn_expon, ndigits, add_expon = 0, midpoint =
> - 0, rounded_up = 0;
> + 0, rounded_up = 0, dround = 0;
> int dec_expon_scale = 0, right_radix_leading_zeros = 0, rdx_pt_enc =
> 0;
> unsigned fpsc;
> @@ -419,10 +419,10 @@ bid64_from_string (char *ps
> break;
>
> case ROUNDING_DOWN:
> - if(sign_x) { coefficient_x++; rounded_up=1; }
> + if(sign_x) { if(c>'0') {coefficient_x++; rounded_up=1;} else dround=1; }
> break;
> case ROUNDING_UP:
> - if(!sign_x) { coefficient_x++; rounded_up=1; }
> + if(!sign_x) { if(c>'0') {coefficient_x++; rounded_up=1;} else dround=1; }
> break;
> case ROUNDING_TIES_AWAY:
> if(c>='5') { coefficient_x++; rounded_up=1; }
> @@ -443,8 +443,21 @@ bid64_from_string (char *ps
> midpoint = 0;
> rounded_up = 1;
> }
> - if (c > '0')
> + if (c > '0') {
> rounded = 1;
> +
> + if(dround)
> + {
> + dround = 0;
> + coefficient_x ++;
> + rounded_up = 1;
> +
> + if (coefficient_x == 10000000000000000ull) {
> + coefficient_x = 1000000000000000ull;
> + add_expon = 1;
> + }
> + }
> + }
> }
> ps++;
> c = *ps;
> diff --git a/libgcc/config/libbid/bid64_to_bid128.c b/libgcc/config/libbid/bid64_to_bid128.c
> index 6e55ba24e02..a8daddce6bc 100644
> --- a/libgcc/config/libbid/bid64_to_bid128.c
> +++ b/libgcc/config/libbid/bid64_to_bid128.c
> @@ -153,9 +153,6 @@ bid128_to_bid64 (UINT128 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
> T128 = round_const_table_128[rmode][extra_digits];
> __add_carry_out (CX1.w[0], carry, T128.w[0], CX.w[0]);
> CX1.w[1] = CX.w[1] + T128.w[1] + carry;
> - if (__unsigned_compare_ge_128
> - (CX1, power10_table_128[extra_digits + 16]))
> - uf_check = 0;
> }
> extra_digits =
> extra_digits + DECIMAL_EXPONENT_BIAS_128 -
> diff --git a/libgcc/config/libbid/bid_binarydecimal.c b/libgcc/config/libbid/bid_binarydecimal.c
> index 5b5b721ffb3..6df39f69887 100644
> --- a/libgcc/config/libbid/bid_binarydecimal.c
> +++ b/libgcc/config/libbid/bid_binarydecimal.c
> @@ -566,19 +566,19 @@ BID_BINARY80LDOUBLE;
> { if ((x & (0xFull<<27)) == (0xFull<<27)) \
> { if ((x & (0x1Full<<26)) != (0x1Full<<26)) inf; \
> if ((x & (1ul<<25))!=0) *pfpsf |= INVALID_EXCEPTION; \
> - nan(s,((((x) & 0xFFFFul) > 999999ul) ? 0 : \
> + nan(s,((((x) & 0xFFFFFul) > 999999ul) ? 0 : \
> (((unsigned long long) x) << 44)),0ull); \
> } \
> e = ((x >> 21) & ((1ull<<8)-1)) - 101; \
> c = (1ull<<23) + (x & ((1ull<<21)-1)); \
> - if ((unsigned long)(c) > 9999999ul) c = 0; \
> + if ((unsigned long)(c) > 9999999ul) zero; \
> k = 0; \
> } \
> else \
> { e = ((x >> 23) & ((1ull<<8)-1)) - 101; \
> c = x & ((1ull<<23)-1); \
> if (c == 0) zero; \
> - k = clz32(c) - 8; \
> + k = clz32_nz(c) - 8; \
> c = c << k; \
> } \
> }
> @@ -594,14 +594,14 @@ BID_BINARY80LDOUBLE;
> } \
> e = ((x >> 51) & ((1ull<<10)-1)) - 398; \
> c = (1ull<<53) + (x & ((1ull<<51)-1)); \
> - if ((unsigned long long)(c) > 9999999999999999ull) c = 0; \
> + if ((unsigned long long)(c) > 9999999999999999ull) zero; \
> k = 0; \
> } \
> else \
> { e = ((x >> 53) & ((1ull<<10)-1)) - 398; \
> c = x & ((1ull<<53)-1); \
> if (c == 0) zero; \
> - k = clz64(c) - 10; \
> + k = clz64_nz(c) - 10; \
> c = c << k; \
> } \
> }
> @@ -144302,20 +144302,6 @@ bid32_to_binary64 (UINT32 x
> // We actually check if e >= ceil((sci_emax + 1) * log_10(2))
> // which in this case is e >= ceil(1024 * log_10(2)) = ceil(308.25) = 309
>
> - if (e >= 309) {
> - *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
> - return_binary64_ovf (s);
> - }
> -// Also check for "trivial" underflow, when 10^e * 2^113 <= 2^emin * 1/4,
> -// so test e <= floor((emin - 115) * log_10(2))
> -// In this case just fix ourselves at that value for uniformity.
> -//
> -// This is important not only to keep the tables small but to maintain the
> -// testing of the round/sticky words as a correct rounding method
> -
> - if (e <= -358)
> - e = -358;
> -
> // Look up the breakpoint and approximate exponent
>
> m_min = (breakpoints_binary64 + 358)[e];
> @@ -144323,7 +144309,7 @@ bid32_to_binary64 (UINT32 x
>
> // Choose provisional exponent and reciprocal multiplier based on breakpoint
>
> - if (le128 (c.w[1], c.w[0], m_min.w[1], m_min.w[0])) {
> + if (c.w[1] < m_min.w[1]) {
> r = (multipliers1_binary64 + 358)[e];
> } else {
> r = (multipliers2_binary64 + 358)[e];
> @@ -144332,17 +144318,12 @@ bid32_to_binary64 (UINT32 x
>
> // Do the reciprocal multiplication
>
> - __mul_128x256_to_384 (z, c, r)
> + __mul_64x256_to_320(z, c.w[1], r);
> + z.w[5]=z.w[4]; z.w[4]=z.w[3]; z.w[3]=z.w[2]; z.w[2]=z.w[1]; z.w[1]=z.w[0]; z.w[0]=0;
> +
> // Check for exponent underflow and compensate by shifting the product
> // Cut off the process at precision+2, since we can't really shift further
> - if (e_out < 1) {
> - int d;
> - d = 1 - e_out;
> - if (d > 55)
> - d = 55;
> - e_out = 1;
> - srl256 (z.w[5], z.w[4], z.w[3], z.w[2], d);
> - }
> +
> c_prov = z.w[5];
>
> // Round using round-sticky words
> @@ -144353,31 +144334,14 @@ bid32_to_binary64 (UINT32 x
> w[1],
> roundbound_128[(rnd_mode << 2) + ((s & 1) << 1) +
> (c_prov & 1)].w[0], z.w[4], z.w[3])) {
> - c_prov = c_prov + 1;
> - if (c_prov == (1ull << 53)) {
> - c_prov = 1ull << 52;
> - e_out = e_out + 1;
> - }
> - }
> -// Check for overflow
> -
> - if (e_out >= 2047) {
> - *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
> - return_binary64_ovf (s);
> + c_prov = c_prov + 1;
> }
> -// Modify exponent for a tiny result, otherwise lop the implicit bit
> -
> - if (c_prov < (1ull << 52))
> - e_out = 0;
> - else
> - c_prov = c_prov & ((1ull << 52) - 1);
> + c_prov = c_prov & ((1ull << 52) - 1);
>
> // Set the inexact and underflow flag as appropriate
>
> if ((z.w[4] != 0) || (z.w[3] != 0)) {
> *pfpsf |= INEXACT_EXCEPTION;
> - if (e_out == 0)
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> }
> // Package up the result as a binary floating-point number
>
> @@ -145756,6 +145720,14 @@ binary64_to_bid32 (double x
> __mul_128x256_to_384 (z, c, r)
> c_prov = z.w[5];
>
> +// Test inexactness and underflow (when testing tininess before rounding)
> +
> + if ((z.w[4] != 0) || (z.w[3] != 0)) {
> + *pfpsf |= INEXACT_EXCEPTION;
> + if (c_prov < 1000000ull)
> + *pfpsf |= UNDERFLOW_EXCEPTION;
> + }
> +
> // Round using round-sticky words
> // If we spill over into the next decade, correct
> // Flag underflow where it may be needed even for |result| = SNN
> @@ -145769,27 +145741,16 @@ binary64_to_bid32 (double x
> if (c_prov == 10000000ull) {
> c_prov = 1000000ull;
> e_out = e_out + 1;
> - } else if ((c_prov == 1000000ull) && (e_out == 0)) {
> - if ((((rnd_mode & 3) == 0) && (z.w[4] <= 17524406870024074035ull))
> - || ((rnd_mode + (s & 1) == 2)
> - && (z.w[4] <= 16602069666338596454ull)))
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> }
> }
> +
> // Check for overflow
>
> if (e_out > 90 + 101) {
> *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
> return_bid32_ovf (s);
> }
> -// Set the inexact flag as appropriate and check underflow
> -// It's no doubt superfluous to check inexactness, but anyway...
>
> - if ((z.w[4] != 0) || (z.w[3] != 0)) {
> - *pfpsf |= INEXACT_EXCEPTION;
> - if (c_prov < 1000000ull)
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> - }
> // Package up the result
>
> return_bid32 (s, e_out, c_prov);
> @@ -145919,6 +145880,14 @@ binary80_to_bid32 (BINARY80 x
> __mul_128x256_to_384 (z, c, r)
> c_prov = z.w[5];
>
> +// Test inexactness and underflow (when testing tininess before rounding)
> +
> + if ((z.w[4] != 0) || (z.w[3] != 0)) {
> + *pfpsf |= INEXACT_EXCEPTION;
> + if (c_prov < 1000000ull)
> + *pfpsf |= UNDERFLOW_EXCEPTION;
> + }
> +
> // Round using round-sticky words
> // If we spill over into the next decade, correct
> // Flag underflow where it may be needed even for |result| = SNN
> @@ -145932,27 +145901,16 @@ binary80_to_bid32 (BINARY80 x
> if (c_prov == 10000000ull) {
> c_prov = 1000000ull;
> e_out = e_out + 1;
> - } else if ((c_prov == 1000000ull) && (e_out == 0)) {
> - if ((((rnd_mode & 3) == 0) && (z.w[4] <= 17524406870024074035ull))
> - || ((rnd_mode + (s & 1) == 2)
> - && (z.w[4] <= 16602069666338596454ull)))
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> }
> }
> +
> // Check for overflow
>
> if (e_out > 90 + 101) {
> *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
> return_bid32_ovf (s);
> }
> -// Set the inexact flag as appropriate and check underflow
> -// It's no doubt superfluous to check inexactness, but anyway...
>
> - if ((z.w[4] != 0) || (z.w[3] != 0)) {
> - *pfpsf |= INEXACT_EXCEPTION;
> - if (c_prov < 1000000ull)
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> - }
> // Package up the result
>
> return_bid32 (s, e_out, c_prov);
> @@ -146071,6 +146029,13 @@ binary128_to_bid32 (BINARY128 x
> __mul_128x256_to_384 (z, c, r)
> c_prov = z.w[5];
>
> +// Test inexactness and underflow (when testing tininess before rounding)
> + if ((z.w[4] != 0) || (z.w[3] != 0)) {
> + *pfpsf |= INEXACT_EXCEPTION;
> + if (c_prov < 1000000ull)
> + *pfpsf |= UNDERFLOW_EXCEPTION;
> + }
> +
> // Round using round-sticky words
> // If we spill over into the next decade, correct
> // Flag underflow where it may be needed even for |result| = SNN
> @@ -146086,30 +146051,16 @@ binary128_to_bid32 (BINARY128 x
> if (c_prov == 10000000ull) {
> c_prov = 1000000ull;
> e_out = e_out + 1;
> - } else if ((c_prov == 1000000ull) && (e_out == 0)) {
> - if ((((rnd_mode & 3) == 0) &&
> - le128 (z.w[4], z.w[3],
> - 17524406870024074035ull, 3689348814741910323ull)) ||
> - ((rnd_mode + (s & 1) == 2) &&
> - le128 (z.w[4], z.w[3],
> - 16602069666338596454ull, 7378697629483820646ull)))
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> }
> }
> +
> // Check for overflow
>
> if (e_out > 90 + 101) {
> *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
> return_bid32_ovf (s);
> }
> -// Set the inexact flag as appropriate and check underflow
> -// It's no doubt superfluous to check inexactness, but anyway...
>
> - if ((z.w[4] != 0) || (z.w[3] != 0)) {
> - *pfpsf |= INEXACT_EXCEPTION;
> - if (c_prov < 1000000ull)
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> - }
> // Package up the result
>
> return_bid32 (s, e_out, c_prov);
> @@ -146562,6 +146513,14 @@ binary80_to_bid64 (BINARY80 x
> __mul_128x256_to_384 (z, c, r)
> c_prov = z.w[5];
>
> +// Test inexactness and underflow (when testing tininess before rounding)
> +
> + if ((z.w[4] != 0) || (z.w[3] != 0)) {
> + *pfpsf |= INEXACT_EXCEPTION;
> + if (c_prov < 1000000000000000ull)
> + *pfpsf |= UNDERFLOW_EXCEPTION;
> + }
> +
> // Round using round-sticky words
> // If we spill over into the next decade, correct
> // Flag underflow where it may be needed even for |result| = SNN
> @@ -146575,27 +146534,16 @@ binary80_to_bid64 (BINARY80 x
> if (c_prov == 10000000000000000ull) {
> c_prov = 1000000000000000ull;
> e_out = e_out + 1;
> - } else if ((c_prov == 1000000000000000ull) && (e_out == 0)) {
> - if ((((rnd_mode & 3) == 0) && (z.w[4] <= 17524406870024074035ull))
> - || ((rnd_mode + (s & 1) == 2)
> - && (z.w[4] <= 16602069666338596454ull)))
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> }
> }
> +
> // Check for overflow
>
> if (e_out > 369 + 398) {
> *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
> return_bid64_ovf (s);
> }
> -// Set the inexact flag as appropriate and check underflow
> -// It's no doubt superfluous to check inexactness, but anyway...
>
> - if ((z.w[4] != 0) || (z.w[3] != 0)) {
> - *pfpsf |= INEXACT_EXCEPTION;
> - if (c_prov < 1000000000000000ull)
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> - }
> // Package up the result
>
> return_bid64 (s, e_out, c_prov);
> @@ -146723,6 +146671,14 @@ binary128_to_bid64 (BINARY128 x
> __mul_128x256_to_384 (z, c, r)
> c_prov = z.w[5];
>
> +// Test inexactness and underflow (when testing tininess before rounding)
> +
> + if ((z.w[4] != 0) || (z.w[3] != 0)) {
> + *pfpsf |= INEXACT_EXCEPTION;
> + if (c_prov < 1000000000000000ull)
> + *pfpsf |= UNDERFLOW_EXCEPTION;
> + }
> +
> // Round using round-sticky words
> // If we spill over into the next decade, correct
> // Flag underflow where it may be needed even for |result| = SNN
> @@ -146736,27 +146692,16 @@ binary128_to_bid64 (BINARY128 x
> if (c_prov == 10000000000000000ull) {
> c_prov = 1000000000000000ull;
> e_out = e_out + 1;
> - } else if ((c_prov == 1000000000000000ull) && (e_out == 0)) {
> - if ((((rnd_mode & 3) == 0) && (z.w[4] <= 17524406870024074035ull))
> - || ((rnd_mode + (s & 1) == 2)
> - && (z.w[4] <= 16602069666338596454ull)))
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> }
> }
> +
> // Check for overflow
>
> if (e_out > 369 + 398) {
> *pfpsf |= (OVERFLOW_EXCEPTION | INEXACT_EXCEPTION);
> return_bid64_ovf (s);
> }
> -// Set the inexact flag as appropriate and check underflow
> -// It's no doubt superfluous to check inexactness, but anyway...
>
> - if ((z.w[4] != 0) || (z.w[3] != 0)) {
> - *pfpsf |= INEXACT_EXCEPTION;
> - if (c_prov < 1000000000000000ull)
> - *pfpsf |= UNDERFLOW_EXCEPTION;
> - }
> // Package up the result
>
> return_bid64 (s, e_out, c_prov);
> diff --git a/libgcc/config/libbid/bid_conf.h b/libgcc/config/libbid/bid_conf.h
> index dff938edd5d..e054a3ff570 100644
> --- a/libgcc/config/libbid/bid_conf.h
> +++ b/libgcc/config/libbid/bid_conf.h
> @@ -519,6 +519,14 @@ see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
> #define BID_BIG_ENDIAN __FLOAT_WORD_ORDER__ == __ORDER_BIG_ENDIAN__
> #endif
>
> +#if BID_BIG_ENDIAN
> +#define BID_HIGH_128W 0
> +#define BID_LOW_128W 1
> +#else
> +#define BID_HIGH_128W 1
> +#define BID_LOW_128W 0
> +#endif
> +
> #ifndef BID_THREAD
> #if defined (HAVE_CC_TLS) && defined (USE_TLS)
> #define BID_THREAD __thread
> diff --git a/libgcc/config/libbid/bid_functions.h b/libgcc/config/libbid/bid_functions.h
> index 9b4fbef8d01..cce87eb4d60 100644
> --- a/libgcc/config/libbid/bid_functions.h
> +++ b/libgcc/config/libbid/bid_functions.h
> @@ -67,9 +67,13 @@ ALIGN (16)
> #endif
>
>
> +#if defined __NO_BINARY80__
> +#define __ENABLE_BINARY80__ 0
> +#else
> #if !defined _MSC_VER || defined __INTEL_COMPILER
> #define __ENABLE_BINARY80__ 1
> #endif
> +#endif
>
> #ifndef HPUX_OS
> #define BINARY80 long double
> @@ -91,6 +95,19 @@ ALIGN (16)
> } UINT256;
> typedef unsigned int FPSC; // floating-point status and control
>
> + typedef enum class_types {
> + signalingNaN,
> + quietNaN,
> + negativeInfinity,
> + negativeNormal,
> + negativeSubnormal,
> + negativeZero,
> + positiveZero,
> + positiveSubnormal,
> + positiveNormal,
> + positiveInfinity
> + } class_t;
> +
> // TYPE parameters
> #define BID128_MAXDIGITS 34
> #define BID64_MAXDIGITS 16
> @@ -2948,7 +2965,7 @@ ALIGN (16)
> extern UINT64 bid64_copySign (UINT64 x,
> UINT64 y _EXC_MASKS_PARAM
> _EXC_INFO_PARAM);
> - extern int bid64_class (UINT64 x _EXC_MASKS_PARAM _EXC_INFO_PARAM);
> + extern class_t bid64_class (UINT64 x _EXC_MASKS_PARAM _EXC_INFO_PARAM);
> extern int bid64_sameQuantum (UINT64 x, UINT64 y
> _EXC_MASKS_PARAM _EXC_INFO_PARAM);
> extern int bid64_totalOrder (UINT64 x, UINT64 y
> @@ -2984,8 +3001,8 @@ ALIGN (16)
> extern UINT128 bid128_copySign (UINT128 x,
> UINT128 y _EXC_MASKS_PARAM
> _EXC_INFO_PARAM);
> - extern int bid128_class (UINT128 x _EXC_MASKS_PARAM
> - _EXC_INFO_PARAM);
> + extern class_t bid128_class (UINT128 x _EXC_MASKS_PARAM
> + _EXC_INFO_PARAM);
> extern int bid128_sameQuantum (UINT128 x,
> UINT128 y _EXC_MASKS_PARAM
> _EXC_INFO_PARAM);
> diff --git a/libgcc/config/libbid/bid_inline_add.h b/libgcc/config/libbid/bid_inline_add.h
> index eb16ea8e0ba..9ca41a2f89e 100644
> --- a/libgcc/config/libbid/bid_inline_add.h
> +++ b/libgcc/config/libbid/bid_inline_add.h
> @@ -918,6 +918,7 @@ get_add128 (UINT64 sign_x, int exponent_x, UINT64 coefficient_x,
> coefficient_x += D;
> }
> break;
> + default: break; // default added to avoid compiler warning
> }
> if (coefficient_x < 1000000000000000ull) {
> coefficient_x -= D;
> @@ -1107,6 +1108,7 @@ get_add128 (UINT64 sign_x, int exponent_x, UINT64 coefficient_x,
> } else if (FS.w[1] | FS.w[0])
> CYh++;
> break;
> + default: break; // default added to avoid compiler warning
> }
> #endif
> #endif
> diff --git a/libgcc/config/libbid/bid_internal.h b/libgcc/config/libbid/bid_internal.h
> index 764abf81057..d3c83b5c219 100644
> --- a/libgcc/config/libbid/bid_internal.h
> +++ b/libgcc/config/libbid/bid_internal.h
> @@ -970,6 +970,8 @@ get_BID64 (UINT64 sgn, int expon, UINT64 coeff, int rmode,
> // round up
> if (sgn)
> r = SMALLEST_BID64;
> + default:
> + break;
> }
> return r;
> }
> @@ -1086,6 +1088,8 @@ fast_get_BID64_check_OF (UINT64 sgn, int expon, UINT64 coeff, int rmode,
> // round up
> if (sgn)
> r = SMALLEST_BID64;
> + default:
> + break;
> }
> return r;
> }
> @@ -2582,19 +2586,6 @@ ALIGN (16)
> A=((tempx.i >>23) & EXPONENT_MASK32) - 0x7f;\
> }
>
> - enum class_types {
> - signalingNaN,
> - quietNaN,
> - negativeInfinity,
> - negativeNormal,
> - negativeSubnormal,
> - negativeZero,
> - positiveZero,
> - positiveSubnormal,
> - positiveNormal,
> - positiveInfinity
> - };
> -
> typedef union {
> UINT64 ui64;
> double d;
> --
> 2.31.1
>
^ permalink raw reply [flat|nested] 2+ messages in thread
end of thread, other threads:[~2024-05-02 9:20 UTC | newest]
Thread overview: 2+ messages (download: mbox.gz / follow: Atom feed)
-- links below jump to the message on this page --
2024-04-28 5:50 [PATCH] Update libbid according to the latest Intel Decimal Floating-Point Math Library liuhongt
2024-05-02 9:20 ` Richard Biener
This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).