This is the mail archive of the
glibc-bugs@sourceware.org
mailing list for the glibc project.
[Bug math/15359] Inaccurate cos for IBM long double
- From: "carlos at redhat dot com" <sourceware-bugzilla at sourceware dot org>
- To: glibc-bugs at sourceware dot org
- Date: Fri, 12 Apr 2013 06:53:33 +0000
- Subject: [Bug math/15359] Inaccurate cos for IBM long double
- Auto-submitted: auto-generated
- References: <bug-15359-131 at http dot sourceware dot org/bugzilla/>
http://sourceware.org/bugzilla/show_bug.cgi?id=15359
Carlos O'Donell <carlos at redhat dot com> changed:
What |Removed |Added
----------------------------------------------------------------------------
Status|NEW |WAITING
--- Comment #3 from Carlos O'Donell <carlos at redhat dot com> 2013-04-12 06:53:33 UTC ---
OK I have a PPC system to test this with setup now.
pi/2 rounded to an IBM long double is:
1.570796326794896619231321691639740613531845480547870673025773690344164545962257761857472360134124755859375L
An infinite precision cos(x) of this value yields:
cos(1.570796326794896619231321691639740613531845480547870673025773690344164545962257761857472360134124755859375)=
1.0828566739219139682237461698605809743657180846737456545052536933566509276210257011310354810728258977667687179060601794311587878461976409993637130916602872593217973444675992784631595571310528753415771366618685685565006009288560082552130...
à 10^-32
The highest precision answer in an IBM long double is:
1.082856673921913968223746169860583539424442830597581817048602989306514317306816102445404953931906744377790888342946371404131442949282377978537095231104103731922805309295654296875e-32
Notice that the infinite precision answer differs from rounded IBM long double
value around ~1e-64? That's no mistake, that's around 1ulp. So it should be
possible to attain this value.
However the ppc cos(x) implemetnation yields:
1.0828566739219139682237461698548710058685937356711019144992345158825349842613103435536910382894240130262349699496720470547777605274664125545314163900911808013916015625e-32
They differ at around the 29th digit, or a difference of
5.71253355584909492648e-62.
Comparing hex digits:
Is: 0x1.c1cd129024e088a67cc7402000p-107
Should be: 0x1.c1cd129024e088a67cc74020bcp-107
It shows the answer is not as precise as the IBM long double type would allow.
A rough guess is that we have 8-12 bits lost.
Therefore there is nothing wrong with the test, it is now correctly showing
that there are 376 ulps of error in the IBM long double implementation of
cos(x). Previously we compared to zero and used a huge value for ulp, ignoring
a lot of inaccuracy around zero.
Note: At the expected value for cos(~pi/2) the value of 1 ulp is 1.5193e-64
Do you see any flaws in that analysis?
The patch I'll probably apply is to convert all the results to hex literals
with the maximum precision for the type.
e.g.
diff --git a/math/libm-test.inc b/math/libm-test.inc
index 0049fcd..f0eb81a 100644
--- a/math/libm-test.inc
+++ b/math/libm-test.inc
@@ -5337,23 +5337,23 @@ cos_test (void)
to each type. */
#ifdef TEST_FLOAT
/* 32-bit float. */
- TEST_f_f (cos, M_PI_2l, -4.371139000186241438857289400265215e-8L);
+ TEST_f_f (cos, M_PI_2l, -0x1.777a5cp-25L);
#endif
#if defined TEST_DOUBLE || (defined TEST_LDOUBLE && LDBL_MANT_DIG == 53)
/* 64-bit double or 64-bit long double. */
- TEST_f_f (cos, M_PI_2l, 6.123233995736765886130329661375001e-17L);
+ TEST_f_f (cos, M_PI_2l, 0x1.1a62633145c07p-54L);
#endif
#if defined TEST_LDOUBLE && LDBL_MANT_DIG == 64
/* 96-bit long double. */
- TEST_f_f (cos, M_PI_2l, -2.50827880633416601177866354016537e-20L);
+ TEST_f_f (cos, M_PI_2l, -0xe.ce675d1fc8f8cbbp-69L);
#endif
#if defined TEST_LDOUBLE && LDBL_MANT_DIG == 106
/* 128-bit IBM long double. */
- TEST_f_f (cos, M_PI_2l, 1.082856673921913968223746169860580e-32L);
+ TEST_f_f (cos, M_PI_2l, 0x1.c1cd129024e088a67cc74020bcp-107L);
#endif
#if defined TEST_LDOUBLE && LDBL_MANT_DIG == 113
/* 128-bit long double. */
- TEST_f_f (cos, M_PI_2l, 4.335905065061890512398522013021675e-35L);
+ TEST_f_f (cos, M_PI_2l, 0x1.cd129024e088a67cc74020bbea64p-115L);
#endif
TEST_f_f (cos, 0.75L, 0.731688868873820886311838753000084544L);
...
--
Configure bugmail: http://sourceware.org/bugzilla/userprefs.cgi?tab=email
------- You are receiving this mail because: -------
You are on the CC list for the bug.