[PATCHv2] Increase internal precision of ldbl-128ibm decimal printf
Commit Message
Changes from V1
* Macrotize the internal precision constant
* Modify tst-sprintf3 to catch this case
* Fix my incorrect overshift of a denormal
Realistically, should we only add 1 bit of extra precision? That
should be sufficient to know which representable number we fall
closer to. 7 is consistent with how the number is presented in
hex though.
----8<----
When the signs differ, the precision of the conversion sometimes
drops below 106 bits. This strategy is identical to the
hexadecimal variant.
I've refactored tst-sprintf3 to enable testing a value with more
than 30 significant digits in order to demonstrate this failure
and its solution.
Additionally, this implicitly fixes a typo in the shift
quantities when subtracting from the high mantissa to compute
the difference.
2016-03-17 Paul E. Murphy <murphyp@linux.vnet.ibm.com>
* stdio-common/tst-sprintf3.c (TEST_N): Refactor
TEST to take significant digits as second parameter.
(TEST): Redefine in terms of TEST_N taking 30
significant digits.
(do_test): Add test case to demonstrate precision
failure in the ldbl-128ibm printf.
* sysdeps/ieee754/ldbl-128ibm/ldbl2pm.c:
(__mpn_extract_long_double): Carry 7 extra intermediate
bits of precision to aide computing difference when
signs differ.
---
stdio-common/tst-sprintf3.c | 14 ++++++++++----
sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c | 25 ++++++++++++++++++-------
2 files changed, 28 insertions(+), 11 deletions(-)
Comments
On Mon, 21 Mar 2016, Paul E. Murphy wrote:
> 2016-03-17 Paul E. Murphy <murphyp@linux.vnet.ibm.com>
>
> * stdio-common/tst-sprintf3.c (TEST_N): Refactor
Please file a bug in Bugzilla for this issue (as something that was
user-visible in releases) and include the bug number in the patch
submission ([BZ #N] in the ChangeLog entry).
I've opened and assigned BZ 19853.
The updated changelog:
2016-03-17 Paul E. Murphy <murphyp@linux.vnet.ibm.com>
[BZ #19853]
* stdio-common/tst-sprintf3.c [TEST_N]: Refactor
TEST to take significant digits as second parameter.
[TEST]: Redefine in terms of TEST_N taking 30
significant digits.
(do_test): Add test case to demonstrate precision
failure in the ldbl-128ibm printf.
* sysdeps/ieee754/ldbl-128ibm/ldbl2pm.c:
(__mpn_extract_long_double): Carry 7 extra intermediate
bits of precision to aide computing difference when
signs differ.
Ping?
On 03/22/2016 09:16 AM, Paul E. Murphy wrote:
> I've opened and assigned BZ 19853.
>
> The updated changelog:
>
> 2016-03-17 Paul E. Murphy <murphyp@linux.vnet.ibm.com>
>
> [BZ #19853]
> * stdio-common/tst-sprintf3.c [TEST_N]: Refactor
> TEST to take significant digits as second parameter.
> [TEST]: Redefine in terms of TEST_N taking 30
> significant digits.
> (do_test): Add test case to demonstrate precision
> failure in the ldbl-128ibm printf.
> * sysdeps/ieee754/ldbl-128ibm/ldbl2pm.c:
> (__mpn_extract_long_double): Carry 7 extra intermediate
> bits of precision to aide computing difference when
> signs differ.
>
This patch is OK.
It's arguable how different this bug actually is from 5268 (the original
submission of bug 5268 doesn't clearly define the scope of the bug, and
the one-bit-loss issue here gets mentioned later in the bug 5268
comments). But I suppose this bug can be considered more specific to the
one-bit-loss, with 5268 left for the more general case where there are
more than 106 bits of precision in the long double value.
Thanks Joseph. Committed as 37a4c70bd4c5c74ac562072e450dc02e8cb4c150.
On 03/30/2016 06:04 PM, Joseph Myers wrote:
> This patch is OK.
>
> It's arguable how different this bug actually is from 5268 (the original
> submission of bug 5268 doesn't clearly define the scope of the bug, and
> the one-bit-loss issue here gets mentioned later in the bug 5268
> comments). But I suppose this bug can be considered more specific to the
> one-bit-loss, with 5268 left for the more general case where there are
> more than 106 bits of precision in the long double value.
>
@@ -38,11 +38,11 @@ do_test (void)
# define COMPARE_LDBL(u, v) ((u).l == (v).l)
#endif
-#define TEST(val) \
+#define TEST_N(val, n) \
do \
{ \
u.l = (val); \
- snprintf (buf, sizeof buf, "%.30LgL", u.l); \
+ snprintf (buf, sizeof buf, "%." #n "LgL", u.l); \
if (strcmp (buf, #val) != 0) \
{ \
printf ("Error on line %d: %s != %s\n", __LINE__, buf, #val); \
@@ -50,19 +50,25 @@ do_test (void)
} \
if (sscanf (#val, "%Lg", &v.l) != 1 || !COMPARE_LDBL (u, v)) \
{ \
- printf ("Error sscanf on line %d: %.30Lg != %.30Lg\n", __LINE__, \
- u.l, v.l); \
+ printf ("Error sscanf on line %d: %." #n "Lg != %." #n "Lg\n", \
+ __LINE__, u.l, v.l); \
result = 1; \
} \
/* printf ("%s %Lg %016Lx %016Lx\n", #val, u.l, u.x[0], u.x[1]); */ \
} \
while (0)
+#define TEST(val) TEST_N (val,30)
+
#if LDBL_MANT_DIG >= 106
# if LDBL_MANT_DIG == 106
TEST (2.22507385850719347803989925739e-308L);
TEST (2.22507385850719397210554509863e-308L);
TEST (2.22507385850720088902458687609e-308L);
+
+ /* Verify precision is not lost for long doubles
+ of the form +1.pN,-1.pM. */
+ TEST_N (3.32306998946228968225951765070082e+35L, 34);
# endif
TEST (2.22507385850720138309023271733e-308L);
TEST (2.22507385850720187715587855858e-308L);
@@ -28,6 +28,12 @@
bits (106 for long double) and an integral power of two (MPN
frexpl). */
+
+/* When signs differ, the actual value is the difference between the
+ significant double and the less significant double. Sometimes a
+ bit can be lost when we borrow from the significant mantissa. */
+#define EXTRA_INTERNAL_PRECISION (7)
+
mp_size_t
__mpn_extract_long_double (mp_ptr res_ptr, mp_size_t size,
int *expt, int *is_neg,
@@ -45,10 +51,15 @@ __mpn_extract_long_double (mp_ptr res_ptr, mp_size_t size,
lo = ((long long) u.d[1].ieee.mantissa0 << 32) | u.d[1].ieee.mantissa1;
hi = ((long long) u.d[0].ieee.mantissa0 << 32) | u.d[0].ieee.mantissa1;
+ /* Hold 7 extra bits of precision in the mantissa. This allows
+ the normalizing shifts below to prevent losing precision when
+ the signs differ and the exponents are sufficiently far apart. */
+ lo <<= EXTRA_INTERNAL_PRECISION;
+
/* If the lower double is not a denormal or zero then set the hidden
53rd bit. */
if (u.d[1].ieee.exponent != 0)
- lo |= 1ULL << 52;
+ lo |= 1ULL << (52 + EXTRA_INTERNAL_PRECISION);
else
lo = lo << 1;
@@ -72,12 +83,12 @@ __mpn_extract_long_double (mp_ptr res_ptr, mp_size_t size,
if (u.d[0].ieee.negative != u.d[1].ieee.negative
&& lo != 0)
{
- lo = (1ULL << 53) - lo;
+ lo = (1ULL << (53 + EXTRA_INTERNAL_PRECISION)) - lo;
if (hi == 0)
{
/* we have a borrow from the hidden bit, so shift left 1. */
- hi = 0x0ffffffffffffeLL | (lo >> 51);
- lo = 0x1fffffffffffffLL & (lo << 1);
+ hi = 0x000ffffffffffffeLL | (lo >> (52 + EXTRA_INTERNAL_PRECISION));
+ lo = 0x0fffffffffffffffLL & (lo << 1);
(*expt)--;
}
else
@@ -85,14 +96,14 @@ __mpn_extract_long_double (mp_ptr res_ptr, mp_size_t size,
}
#if BITS_PER_MP_LIMB == 32
/* Combine the mantissas to be contiguous. */
- res_ptr[0] = lo;
- res_ptr[1] = (hi << (53 - 32)) | (lo >> 32);
+ res_ptr[0] = lo >> EXTRA_INTERNAL_PRECISION;
+ res_ptr[1] = (hi << (53 - 32)) | (lo >> (32 + EXTRA_INTERNAL_PRECISION));
res_ptr[2] = hi >> 11;
res_ptr[3] = hi >> (32 + 11);
#define N 4
#elif BITS_PER_MP_LIMB == 64
/* Combine the two mantissas to be contiguous. */
- res_ptr[0] = (hi << 53) | lo;
+ res_ptr[0] = (hi << 53) | (lo >> EXTRA_INTERNAL_PRECISION);
res_ptr[1] = hi >> 11;
#define N 2
#else