[PATCH] Increase internal precision of ldbl-128ibm decimal printf

[PATCH] Increase internal precision of ldbl-128ibm decimal printf

[Paul E. Murphy]

When the signs differ, the precision of the conversion
could drop below 106 bits.  This strategy is
identical to the hexadecimal variant.

0x1.ffffffffffffffffffffffffff80p111L is an example of
a number which should, but does not round trip through
printf as a decimal number when printing with 32 digits
of precision.

Note, the printf test did not pick this up, as it
only tests up to 30 digits.  Long double can have
more at nominal precision.

2016-03-17  Paul E. Murphy  <murphyp@linux.vnet.ibm.com>

	* sysdeps/ieee754/ldbl-128ibm/ldbl2pm.c:
	(__mpn_extract_long_double): Carry 7 extra intermediate
	bits of precision to aide computing difference when
	signs differ.
---
 sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c | 21 +++++++++++++--------
 1 file changed, 13 insertions(+), 8 deletions(-)

diff --git a/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c b/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c
index 4f550ef..afe6cf2 100644
--- a/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c
+++ b/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c
@@ -45,12 +45,17 @@ __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 <<= 7;
+
   /* 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 + 7);
   else
-    lo = lo << 1;
+    lo = lo << (1 + 7);
 
   /* The lower double is normalized separately from the upper.  We may
      need to adjust the lower manitissa to reflect this.  */
@@ -72,12 +77,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 + 7)) - 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 + 7));
+	  lo = 0x0fffffffffffffffLL & (lo << 1);
 	  (*expt)--;
 	}
       else
@@ -85,14 +90,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 >> 7;
+  res_ptr[1] = (hi << (53 - 32)) | (lo >> (32 + 7));
   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 >> 7);
   res_ptr[1] = hi >> 11;
   #define N 2
 #else
-- 
2.4.3

Re: [PATCH] Increase internal precision of ldbl-128ibm decimal printf

[Andreas Schwab]

"Paul E. Murphy" <murphyp@linux.vnet.ibm.com> writes:

> +  /* 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 <<= 7;

That magical number should probably get a name.

>  	  /* we have a borrow from the hidden bit, so shift left 1.  */
> -	  hi = 0x0ffffffffffffeLL | (lo >> 51);
> -	  lo = 0x1fffffffffffffLL & (lo << 1);
> +	  hi = 0x000ffffffffffffeLL | (lo >> (52 + 7));
> +	  lo = 0x0fffffffffffffffLL & (lo << 1);

Not 51 + 7?

Andreas.

-- 
Andreas Schwab, schwab@linux-m68k.org
GPG Key fingerprint = 58CA 54C7 6D53 942B 1756  01D3 44D5 214B 8276 4ED5
"And now for something completely different."

Re: [PATCH] Increase internal precision of ldbl-128ibm decimal printf

[Joseph Myers]

On Thu, 17 Mar 2016, Paul E. Murphy wrote:

> When the signs differ, the precision of the conversion
> could drop below 106 bits.  This strategy is
> identical to the hexadecimal variant.

I assume this was user-visible in a release, in which case it should have 
a bug filed in Bugzilla (if it doesn't already have one - I take it this 
is distinct from bug 5268?).

It should also have a testcase included in the patch (there are several 
tests in sysdeps/ieee754/ldbl-128ibm already for issues specific to this 
format - various of which use unions to test issues that can't be tested 
simply with 106-bit inputs GCC can represent, or where it's necessary to 
examine the two halves of a long double result separately).

-- 
Joseph S. Myers
joseph@codesourcery.com

Re: [PATCH] Increase internal precision of ldbl-128ibm decimal printf

[Paul E. Murphy]

On 03/17/2016 04:56 PM, Andreas Schwab wrote:
> "Paul E. Murphy" <murphyp@linux.vnet.ibm.com> writes:
> 
>> +  /* 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 <<= 7;
> 
> That magical number should probably get a name.

Suggestions? It was noted as an "IEEE 854 adjustment" constant in
printf_fphex.c, and otherwise integrated into the shifts later on.

INTERNAL_PRECISION_SHIFT? Setting it 1 should fix the I case noted
in the comments.

> 
>>  	  /* we have a borrow from the hidden bit, so shift left 1.  */
>> -	  hi = 0x0ffffffffffffeLL | (lo >> 51);
>> -	  lo = 0x1fffffffffffffLL & (lo << 1);
>> +	  hi = 0x000ffffffffffffeLL | (lo >> (52 + 7));
>> +	  lo = 0x0fffffffffffffffLL & (lo << 1);
> 
> Not 51 + 7?

That is intentional.  It should be moving the high bit of lo
lo into position 0.  This is identical to how the ldbl-128ibm
printf_fphex.c, which does round trip these numbers correctly.

I.e ((1<<53)-1)>>51 == 3.  The implicit bit is already factored
into lo, so the existing shift appears wrong here. 

> 
> Andreas.
>

Re: [PATCH] Increase internal precision of ldbl-128ibm decimal printf

[Paul E. Murphy]

On 03/17/2016 05:16 PM, Joseph Myers wrote:
> On Thu, 17 Mar 2016, Paul E. Murphy wrote:
> 
>> When the signs differ, the precision of the conversion
>> could drop below 106 bits.  This strategy is
>> identical to the hexadecimal variant.
> 
> I assume this was user-visible in a release, in which case it should have 
> a bug filed in Bugzilla (if it doesn't already have one - I take it this 
> is distinct from bug 5268?).
> 
> It should also have a testcase included in the patch (there are several 
> tests in sysdeps/ieee754/ldbl-128ibm already for issues specific to this 
> format - various of which use unions to test issues that can't be tested 
> simply with 106-bit inputs GCC can represent, or where it's necessary to 
> examine the two halves of a long double result separately).
> 

This came up internally. I think it is similar, the comments would indicate
it is a 1 bit internal precision loss which this patch does address.

I will figure out how to update the tests to catch this, and give the magic
shift amount a name.

Thanks,
Paul