This is the mail archive of the
libc-alpha@sourceware.org
mailing list for the glibc project.
Fix ldbl-128ibm cos range reduction near pi/2 (bug 15359)
- From: "Joseph S. Myers" <joseph at codesourcery dot com>
- To: <libc-alpha at sourceware dot org>
- Date: Thu, 9 May 2013 20:38:18 +0000
- Subject: Fix ldbl-128ibm cos range reduction near pi/2 (bug 15359)
Bug 15359 is inaccuracy of ldbl-128ibm cos near pi/2. As with the
ldbl-128 case, the cause is inaccurate range reduction, and like that
patch <http://sourceware.org/ml/libc-alpha/2013-05/msg00271.html> this
patch uses the full precision of long double for the high and low
parts of pi/2 to fix the inaccuracy. Tested powerpc to confirm the
large ulps no longer appear.
2013-05-09 Joseph Myers <joseph@codesourcery.com>
[BZ #15359]
* sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c (c): Use 106 bits for
high part of pi/2.
(__ieee754_rem_pio2l): Update comments.
diff --git a/sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c b/sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c
index 692ae24..6a72d6a 100644
--- a/sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c
+++ b/sysdeps/ieee754/ldbl-128ibm/e_rem_pio2l.c
@@ -185,13 +185,13 @@ static const int32_t two_over_pi[] = {
};
static const long double c[] = {
-/* 93 bits of pi/2 */
+/* 106 bits of pi/2 */
#define PI_2_1 c[0]
- 1.57079632679489661923132169155131424e+00L, /* 3fff921fb54442d18469898cc5100000 */
+ 0x1.921fb54442d18469898cc517018p+0L,
/* pi/2 - PI_2_1 */
#define PI_2_1t c[1]
- 8.84372056613570112025531863263659260e-29L, /* 3fa1c06e0e68948127044533e63a0106 */
+ 0x3.839a252049c1114cf98e804178p-108L,
};
int32_t __ieee754_rem_pio2l(long double x, long double *y)
@@ -216,7 +216,7 @@ int32_t __ieee754_rem_pio2l(long double x, long double *y)
{
if (hx > 0)
{
- /* 113 + 93 bit PI is ok */
+ /* 106 + 106 bit PI is ok */
z = x - PI_2_1;
y[0] = z - PI_2_1t;
y[1] = (z - y[0]) - PI_2_1t;
@@ -224,7 +224,7 @@ int32_t __ieee754_rem_pio2l(long double x, long double *y)
}
else
{
- /* 113 + 93 bit PI is ok */
+ /* 106 + 106 bit PI is ok */
z = x + PI_2_1;
y[0] = z + PI_2_1t;
y[1] = (z - y[0]) + PI_2_1t;
--
Joseph S. Myers
joseph@codesourcery.com