From patchwork Mon Mar 21 22:41:38 2016 Content-Type: text/plain; charset="utf-8" MIME-Version: 1.0 Content-Transfer-Encoding: 7bit X-Patchwork-Submitter: "Paul E. Murphy" X-Patchwork-Id: 11467 X-Patchwork-Delegate: joseph@codesourcery.com Received: (qmail 79336 invoked by alias); 21 Mar 2016 22:41:46 -0000 Mailing-List: contact libc-alpha-help@sourceware.org; run by ezmlm Precedence: bulk List-Id: List-Unsubscribe: List-Subscribe: List-Archive: List-Post: List-Help: , Sender: libc-alpha-owner@sourceware.org Delivered-To: mailing list libc-alpha@sourceware.org Received: (qmail 79233 invoked by uid 89); 21 Mar 2016 22:41:45 -0000 Authentication-Results: sourceware.org; auth=none X-Virus-Found: No X-Spam-SWARE-Status: No, score=0.0 required=5.0 tests=AWL, BAYES_20, KAM_LAZY_DOMAIN_SECURITY, RP_MATCHES_RCVD autolearn=no version=3.3.2 spammy=quantities, ud, aide, ul X-HELO: e31.co.us.ibm.com X-IBM-Helo: d03dlp03.boulder.ibm.com X-IBM-MailFrom: murphyp@linux.vnet.ibm.com X-IBM-RcptTo: libc-alpha@sourceware.org From: "Paul E. Murphy" Subject: [PATCHv2] Increase internal precision of ldbl-128ibm decimal printf To: "libc-alpha@sourceware.org" Cc: Tulio Magno Quites Machado Filho Message-ID: <56F078A2.2070602@linux.vnet.ibm.com> Date: Mon, 21 Mar 2016 17:41:38 -0500 User-Agent: Mozilla/5.0 (X11; Linux x86_64; rv:38.0) Gecko/20100101 Thunderbird/38.6.0 MIME-Version: 1.0 X-TM-AS-MML: disable X-Content-Scanned: Fidelis XPS MAILER x-cbid: 16032122-8236-0000-0000-000018741171 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 * 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(-) diff --git a/stdio-common/tst-sprintf3.c b/stdio-common/tst-sprintf3.c index cffd1b6..dad216b 100644 --- a/stdio-common/tst-sprintf3.c +++ b/stdio-common/tst-sprintf3.c @@ -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); diff --git a/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c b/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c index 4f550ef..30d9bc3 100644 --- a/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c +++ b/sysdeps/ieee754/ldbl-128ibm/ldbl2mpn.c @@ -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