Fix x86/x86_64 expm1 inaccuracy near 0 in directed rounding modes (bug 16293).

Bug 16293 is inaccuracy of x86/x86_64 versions of expm1, near 0 in
directed rounding modes, that arises from frndint rounding the
exponent to 1 or -1 instead of 0, resulting in large cancellation
error.  This inaccuracy in turn affects other functions such as sinh
that use expm1.  This patch fixes the problem by setting
round-to-nearest mode temporarily around the affected calls to
frndint.  I don't think this is needed for other uses of frndint, such
as in exp itself, as only for expm1 is the cancellation error
significant.

Tested x86_64 and x86 and ulps updated accordingly.

	* sysdeps/i386/fpu/e_expl.S (IEEE754_EXPL) [USE_AS_EXPM1L]: Set
	round-to-nearest mode when using frndint.
	* sysdeps/i386/fpu/s_expm1.S (__expm1): Likewise.
	* sysdeps/i386/fpu/s_expm1f.S (__expm1f): Likewise.
	* sysdeps/x86_64/fpu/e_expl.S (IEEE754_EXPL) [USE_AS_EXPM1L]:
	Likewise.
	* math/auto-libm-test-in: Add more tests of expm1.  Do not expect
	sinh test to fail.
	* math/auto-libm-test-out: Regenerated.
	* math/libm-test.inc (TEST_COND_x86_64): Remove macro.
	(TEST_COND_x86): Likewise.
	(expm1_tonearest_test_data): New array.
	(expm1_test_tonearest): New function.
	(expm1_towardzero_test_data): New array.
	(expm1_test_towardzero): New function.
	(expm1_downward_test_data): New array.
	(expm1_test_downward): New function.
	(expm1_upward_test_data): New array.
	(expm1_test_upward): New function.
	(main): Run the new test functions.
	* sysdeps/i386/fpu/libm-test-ulps: Update.
	* sysdeps/x86_64/fpu/libm-test-ulps: Likewise.

1 files changed, 375 insertions, 0 deletions

diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps
index a471e5fe9f..bb3adcbd67 100644
--- a/sysdeps/i386/fpu/libm-test-ulps
+++ b/sysdeps/i386/fpu/libm-test-ulps
@@ -6092,9 +6092,15 @@ idouble: 1
 ifloat: 1
 
 # expm1
+Test "expm1 (-0x1p-64)":
+ildouble: 1
+ldouble: 1
 Test "expm1 (-0x2.dp+4)":
 ildouble: 1
 ldouble: 1
+Test "expm1 (-0x4p-12)":
+ildouble: 1
+ldouble: 1
 Test "expm1 (-45.0)":
 ildouble: 1
 ldouble: 1
@@ -6108,6 +6114,338 @@ Test "expm1 (11356.25)":
 ildouble: 1
 ldouble: 1
 
+# expm1_downward
+Test "expm1_downward (-0x1p-100)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (-0x2.ep+4)":
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (-0x4.9p+4)":
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (-0x4.bp+4)":
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (-0x4p-4)":
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (-0x5p+4)":
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (-0x6.4p+4)":
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (0x1p+0)":
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (0x1p-100)":
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (0x1p-32)":
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (0x3.2p+4)":
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (0x4p-12)":
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (0x4p-52)":
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (0x7.fp+4)":
+ildouble: 1
+ldouble: 1
+Test "expm1_downward (0x8p-32)":
+ildouble: 1
+ldouble: 1
+
+# expm1_tonearest
+Test "expm1_tonearest (-0x1p-64)":
+ildouble: 1
+ldouble: 1
+Test "expm1_tonearest (-0x2.dp+4)":
+ildouble: 1
+ldouble: 1
+Test "expm1_tonearest (-0x4p-12)":
+ildouble: 1
+ldouble: 1
+Test "expm1_tonearest (0x1p+0)":
+ildouble: 1
+ldouble: 1
+Test "expm1_tonearest (0x2.c5c4p+12)":
+ildouble: 1
+ldouble: 1
+
+# expm1_towardzero
+Test "expm1_towardzero (-0x1.2p+4)":
+float: 1
+Test "expm1_towardzero (-0x1.86ap+16)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x1p-100)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x1p-20)":
+ildouble: 2
+ldouble: 2
+Test "expm1_towardzero (-0x1p-32)":
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x1p-64)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x2.4p+4)":
+float: 1
+Test "expm1_towardzero (-0x2.5p+4)":
+float: 1
+Test "expm1_towardzero (-0x2.6p+4)":
+double: 1
+float: 1
+Test "expm1_towardzero (-0x2.71p+12)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x2.cp+4)":
+double: 1
+float: 1
+Test "expm1_towardzero (-0x2.dp+4)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x2.ep+4)":
+double: 1
+float: 1
+Test "expm1_towardzero (-0x3.e8p+8)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x4.9p+4)":
+double: 1
+float: 1
+Test "expm1_towardzero (-0x4.ap+4)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x4.bp+4)":
+double: 1
+float: 1
+Test "expm1_towardzero (-0x4.ep+4)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x4.fp+4)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x4p-12)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0x4p-52)":
+idouble: 1
+ildouble: 2
+ldouble: 2
+Test "expm1_towardzero (-0x5p+4)":
+double: 1
+float: 1
+Test "expm1_towardzero (-0x6.4p+4)":
+double: 1
+float: 1
+Test "expm1_towardzero (-0x8p-32)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0xf.ffffffffffff8p+1020)":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0xf.fffffffffffffffp+16380)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (-0xf.fffffp+124)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (0x1p+0)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (0x1p-100)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (0x1p-32)":
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (0x3.2p+4)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (0x4p-12)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (0x4p-52)":
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (0x7.fp+4)":
+ildouble: 1
+ldouble: 1
+Test "expm1_towardzero (0x8p-32)":
+ildouble: 1
+ldouble: 1
+
+# expm1_upward
+Test "expm1_upward (-0x1.2p+4)":
+float: 1
+Test "expm1_upward (-0x1.86ap+16)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x1p-100)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x1p-20)":
+ildouble: 2
+ldouble: 2
+Test "expm1_upward (-0x1p-32)":
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x1p-64)":
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x2.4p+4)":
+float: 1
+Test "expm1_upward (-0x2.5p+4)":
+float: 1
+Test "expm1_upward (-0x2.6p+4)":
+double: 1
+float: 1
+Test "expm1_upward (-0x2.71p+12)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x2.cp+4)":
+double: 1
+float: 1
+Test "expm1_upward (-0x2.dp+4)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x2.ep+4)":
+double: 1
+float: 1
+Test "expm1_upward (-0x3.e8p+8)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x4.9p+4)":
+double: 1
+float: 1
+Test "expm1_upward (-0x4.ap+4)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x4.bp+4)":
+double: 1
+float: 1
+Test "expm1_upward (-0x4.ep+4)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x4.fp+4)":
+double: 1
+float: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x4p-12)":
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0x4p-52)":
+idouble: 1
+ildouble: 2
+ldouble: 2
+Test "expm1_upward (-0x5p+4)":
+double: 1
+float: 1
+Test "expm1_upward (-0x6.4p+4)":
+double: 1
+float: 1
+Test "expm1_upward (-0x8p-32)":
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0xf.ffffffffffff8p+1020)":
+double: 1
+idouble: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0xf.fffffffffffffffp+16380)":
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (-0xf.fffffp+124)":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (0x1.f4p+8)":
+ildouble: 1
+ldouble: 1
+Test "expm1_upward (0x4p-4)":
+ildouble: 1
+ldouble: 1
+
 # gamma
 Test "gamma (-0.5)":
 double: 1
@@ -7773,6 +8111,9 @@ ldouble: 1
 Test "sinh_downward (0x1.8p+4)":
 ildouble: 1
 ldouble: 1
+Test "sinh_downward (0x8p-32)":
+ildouble: 1
+ldouble: 1
 Test "sinh_downward (0xcp-4)":
 float: 1
 ildouble: 1
@@ -7810,6 +8151,9 @@ ldouble: 1
 Test "sinh_towardzero (0x1.8p+4)":
 ildouble: 1
 ldouble: 1
+Test "sinh_towardzero (0x8p-32)":
+ildouble: 1
+ldouble: 1
 Test "sinh_towardzero (0xcp-4)":
 float: 1
 ildouble: 1
@@ -7840,6 +8184,9 @@ ldouble: 1
 Test "sinh_upward (0x1.7p+4)":
 ildouble: 1
 ldouble: 1
+Test "sinh_upward (0x8p-32)":
+double: 1
+float: 1
 Test "sinh_upward (0xcp-4)":
 float: 1
 ildouble: 2
@@ -10663,6 +11010,34 @@ Function: "expm1":
 ildouble: 1
 ldouble: 1
 
+Function: "expm1_downward":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 1
+ldouble: 1
+
+Function: "expm1_tonearest":
+ildouble: 1
+ldouble: 1
+
+Function: "expm1_towardzero":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 2
+ldouble: 2
+
+Function: "expm1_upward":
+double: 1
+float: 1
+idouble: 1
+ifloat: 1
+ildouble: 2
+ldouble: 2
+
 Function: "gamma":
 double: 1
 float: 2