This is the mail archive of the libc-alpha@sourceware.org mailing list for the glibc project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

Fix expm1 missing underflows (bug 16353) [committed]


Similar to various other bugs in this area, some expm1 implementations
do not raise the underflow exception for subnormal arguments, when the
result is tiny and inexact.  This patch forces the exception in a
similar way to previous fixes.

(The issue does not apply to the ldbl-* implementations or to those
for x86 / x86_64 long double.  The change to
sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c is one I missed when
previously fixing bug 16354; the bug in that implementation was
previously latent, but the expm1 fixes stopped it being latent and so
required it to be fixed to avoid spurious underflows from cosh.)

Tested for x86_64 and x86.  Committed.

(auto-libm-test-out diffs omitted below.)

2015-06-22  Joseph Myers  <joseph@codesourcery.com>

	[BZ #16353]
	* sysdeps/i386/fpu/s_expm1.S (dbl_min): New object.
	(__expm1): Force underflow exception for arguments with small
	absolute value.
	* sysdeps/i386/fpu/s_expm1f.S (flt_min): New object.
	(__expm1f): Force underflow exception for arguments with small
	absolute value.
	* sysdeps/ieee754/dbl-64/s_expm1.c: Include <float.h>.
	(__expm1): Force underflow exception for arguments with small
	absolute value.
	* sysdeps/ieee754/flt-32/s_expm1f.c: Include <float.h>.
	(__expm1f): Force underflow exception for arguments with small
	absolute value.
	* sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c (__ieee754_cosh):
	Check for small arguments before calling __expm1.
	* math/auto-libm-test-in: Do not mark underflow exceptions as
	possibly missing for bug 16353.
	* math/auto-libm-test-out: Regenerated.

diff --git a/math/auto-libm-test-in b/math/auto-libm-test-in
index 6b05459..9bae9d5 100644
--- a/math/auto-libm-test-in
+++ b/math/auto-libm-test-in
@@ -1362,12 +1362,11 @@ expm1 0x6.289a78p-4
 expm1 0x6.1b4d318238d4a2a8p-4
 expm1 0x5.fb8dc64e91a74p-4
 expm1 0x3.735f497c4e67535cp-4
-# Bug 16353: underflow exception may be missing
-expm1 0x4.0000000000000028p-16384 missing-underflow
-expm1 min missing-underflow
-expm1 -min missing-underflow
-expm1 min_subnorm missing-underflow
-expm1 -min_subnorm missing-underflow
+expm1 0x4.0000000000000028p-16384
+expm1 min
+expm1 -min
+expm1 min_subnorm
+expm1 -min_subnorm
 
 fma 1.0 2.0 3.0
 fma 1.25 0.75 0.0625
diff --git a/sysdeps/i386/fpu/s_expm1.S b/sysdeps/i386/fpu/s_expm1.S
index 65426c5..05e5285 100644
--- a/sysdeps/i386/fpu/s_expm1.S
+++ b/sysdeps/i386/fpu/s_expm1.S
@@ -37,6 +37,13 @@ one:	.double 1.0
 l2e:	.tfloat 1.442695040888963407359924681002
 	ASM_SIZE_DIRECTIVE(l2e)
 
+	.section .rodata.cst8,"aM",@progbits,8
+
+	.p2align 3
+	.type dbl_min,@object
+dbl_min:	.byte 0, 0, 0, 0, 0, 0, 0x10, 0
+	ASM_SIZE_DIRECTIVE(dbl_min)
+
 #ifdef PIC
 #define MO(op) op##@GOTOFF(%edx)
 #else
@@ -74,6 +81,21 @@ ENTRY(__expm1)
 #ifdef	PIC
 	LOAD_PIC_REG (dx)
 #endif
+	fld	%st
+	fabs
+	fcoml	MO(dbl_min)
+	fstp	%st
+	fnstsw
+	sahf
+	jae	5f
+	subl	$8, %esp
+	cfi_adjust_cfa_offset (8)
+	fld	%st(0)
+	fmul	%st(0)
+	fstpl	(%esp)
+	addl	$8, %esp
+	cfi_adjust_cfa_offset (-8)
+	ret
 
 5:	fldt	MO(l2e)		// log2(e) : x
 	fmulp			// log2(e)*x
diff --git a/sysdeps/i386/fpu/s_expm1f.S b/sysdeps/i386/fpu/s_expm1f.S
index 748b7b8..a83e435 100644
--- a/sysdeps/i386/fpu/s_expm1f.S
+++ b/sysdeps/i386/fpu/s_expm1f.S
@@ -37,6 +37,13 @@ one:	.double 1.0
 l2e:	.tfloat 1.442695040888963407359924681002
 	ASM_SIZE_DIRECTIVE(l2e)
 
+	.section .rodata.cst4,"aM",@progbits,4
+
+	.p2align 2
+	.type flt_min,@object
+flt_min:	.byte 0, 0, 0x80, 0
+	ASM_SIZE_DIRECTIVE(flt_min)
+
 #ifdef PIC
 #define MO(op) op##@GOTOFF(%edx)
 #else
@@ -74,6 +81,21 @@ ENTRY(__expm1f)
 #ifdef	PIC
 	LOAD_PIC_REG (dx)
 #endif
+	fld	%st
+	fabs
+	fcoms	MO(flt_min)
+	fstp	%st
+	fnstsw
+	sahf
+	jae	5f
+	subl	$4, %esp
+	cfi_adjust_cfa_offset (4)
+	fld	%st(0)
+	fmul	%st(0)
+	fstps	(%esp)
+	addl	$4, %esp
+	cfi_adjust_cfa_offset (-4)
+	ret
 
 5:	fldt	MO(l2e)		// log2(e) : x
 	fmulp			// log2(e)*x
diff --git a/sysdeps/ieee754/dbl-64/s_expm1.c b/sysdeps/ieee754/dbl-64/s_expm1.c
index 9c9bb34..41ef63a 100644
--- a/sysdeps/ieee754/dbl-64/s_expm1.c
+++ b/sysdeps/ieee754/dbl-64/s_expm1.c
@@ -109,6 +109,7 @@
  */
 
 #include <errno.h>
+#include <float.h>
 #include <math.h>
 #include <math_private.h>
 #define one Q[0]
@@ -194,6 +195,11 @@ __expm1 (double x)
     }
   else if (hx < 0x3c900000)             /* when |x|<2**-54, return x */
     {
+      if (fabs (x) < DBL_MIN)
+	{
+	  double force_underflow = x * x;
+	  math_force_eval (force_underflow);
+	}
       t = huge + x;     /* return x with inexact flags when x!=0 */
       return x - (t - (huge + x));
     }
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c b/sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c
index 8459352..fca80b1 100644
--- a/sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/e_cosh.c
@@ -50,9 +50,10 @@ __ieee754_cosh (double x)
 	if (ix < 0x40360000) {
 	    /* |x| in [0,0.5*ln2], return 1+expm1(|x|)^2/(2*exp(|x|)) */
 		if(ix<0x3fd62e43) {
+		    if (ix<0x3c800000)			/* cosh(tiny) = 1 */
+		      return one;
 		    t = __expm1(fabs(x));
 		    w = one+t;
-		    if (ix<0x3c800000) return w;	/* cosh(tiny) = 1 */
 		    return one+(t*t)/(w+w);
 		}
 
diff --git a/sysdeps/ieee754/flt-32/s_expm1f.c b/sysdeps/ieee754/flt-32/s_expm1f.c
index ca8839f..c81b057 100644
--- a/sysdeps/ieee754/flt-32/s_expm1f.c
+++ b/sysdeps/ieee754/flt-32/s_expm1f.c
@@ -14,6 +14,7 @@
  */
 
 #include <errno.h>
+#include <float.h>
 #include <math.h>
 #include <math_private.h>
 
@@ -80,6 +81,11 @@ __expm1f(float x)
 	    c  = (hi-x)-lo;
 	}
 	else if(hx < 0x33000000) {	/* when |x|<2**-25, return x */
+	    if (fabsf (x) < FLT_MIN)
+	      {
+		float force_underflow = x * x;
+		math_force_eval (force_underflow);
+	      }
 	    t = huge+x;	/* return x with inexact flags when x!=0 */
 	    return x - (t-(huge+x));
 	}

-- 
Joseph S. Myers
joseph@codesourcery.com


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]