This is the mail archive of the
libc-alpha@sourceware.org
mailing list for the glibc project.
Fix nexttoward bugs (bugs 2550, 2570)
- From: "Joseph S. Myers" <joseph at codesourcery dot com>
- To: libc-alpha at sourceware dot org
- Date: Tue, 1 May 2012 10:03:22 +0000 (UTC)
- Subject: Fix nexttoward bugs (bugs 2550, 2570)
Bug 2550 is a problem where nexttowardf produces incorrect results
when the float operand is subnormal. The same problem also applies to
nexttoward. Bug 2570 is a different problem with nexttoward adjusting
its double operand in the wrong direction.
Both problems arise from floating-point comparisons done by
complicated and buggy integer code working on the representations of
floating-point numbers (one in float or double, one in long double;
nexttowardl has two arguments of the same type and is just an alias
for nextafterl, which is somewhat simpler). The natural fix is just
to use floating-point "<" and ">" (NaNs have already been detected by
this point; there's no need to use unordered operations).
I propose this patch to fix these bugs that way. I adjusted all the
implementations that, from inspection, appear to have at least the bug
relating to subnormal inputs; only the x86 implementations (also used
on x86_64) have had these changes tested, but all the changes are
similar and mechanical. I didn't change the ldbl-128ibm
implementation of nexttoward (double); I don't expect that format to
suffer from subnormal-related bugs, given that the exponent range for
double and long double is the same.
Tested x86 and x86_64. If you see a failure of a nexttoward test on
some other architecture, do not add it to libm-test-ulps; nexttoward
should always produce exact results and any difference from the
expected result is a bug in either nexttoward or the testcase.
2012-05-01 Joseph Myers <joseph@codesourcery.com>
[BZ #2550]
[BZ #2570]
* math/s_nexttowardf.c (__nexttowardf): Use floating-point
comparisons to determine direction to adjust input.
* sysdeps/i386/fpu/s_nexttoward.c (__nexttoward): Likewise.
* sysdeps/i386/fpu/s_nexttowardf.c(__nexttowardf): Likewise.
* sysdeps/ieee754/ldbl-128/s_nexttoward.c (__nexttoward):
Likewise.
* sysdeps/ieee754/ldbl-128/s_nexttowardf.c (__nexttowardf):
Likewise.
* sysdeps/ieee754/ldbl-128ibm/s_nexttowardf.c (__nexttowardf):
Likewise.
* sysdeps/ieee754/ldbl-96/s_nexttoward.c (__nexttoward): Likewise.
* sysdeps/ieee754/ldbl-96/s_nexttowardf.c (__nexttowardf):
Likewise.
* sysdeps/ieee754/ldbl-opt/s_nexttowardfd.c (__nldbl_nexttowardf):
Likewise.
* math/libm-test.inc (nexttoward_test): Add more tests.
diff --git a/math/libm-test.inc b/math/libm-test.inc
index d643bad..0875e2c 100644
--- a/math/libm-test.inc
+++ b/math/libm-test.inc
@@ -5743,8 +5743,73 @@ nexttoward_test (void)
TEST_ff_f (nexttoward, 1.1L, nan_value, nan_value);
TEST_ff_f (nexttoward, nan_value, nan_value, nan_value);
- /* XXX We need the hexadecimal FP number representation here for further
- tests. */
+#ifdef TEST_FLOAT
+ TEST_ff_f (nexttoward, 1.0, 1.1L, 0x1.000002p0);
+ TEST_ff_f (nexttoward, 1.0, LDBL_MAX, 0x1.000002p0);
+ TEST_ff_f (nexttoward, 1.0, 0x1.0000000000001p0, 0x1.000002p0);
+ TEST_ff_f (nexttoward, 1.0, 0.9L, 0x0.ffffffp0);
+ TEST_ff_f (nexttoward, 1.0, -LDBL_MAX, 0x0.ffffffp0);
+ TEST_ff_f (nexttoward, 1.0, 0x0.fffffffffffff8p0, 0x0.ffffffp0);
+ TEST_ff_f (nexttoward, -1.0, -1.1L, -0x1.000002p0);
+ TEST_ff_f (nexttoward, -1.0, -LDBL_MAX, -0x1.000002p0);
+ TEST_ff_f (nexttoward, -1.0, -0x1.0000000000001p0, -0x1.000002p0);
+ TEST_ff_f (nexttoward, -1.0, -0.9L, -0x0.ffffffp0);
+ TEST_ff_f (nexttoward, -1.0, LDBL_MAX, -0x0.ffffffp0);
+ TEST_ff_f (nexttoward, -1.0, -0x0.fffffffffffff8p0, -0x0.ffffffp0);
+ TEST_ff_f (nexttoward, -0x1.3p-145, -0xap-148L, -0x1.4p-145);
+# if LDBL_MANT_DIG >= 64
+ TEST_ff_f (nexttoward, 1.0, 0x1.000000000000002p0L, 0x1.000002p0);
+ TEST_ff_f (nexttoward, 1.0, 0x0.ffffffffffffffffp0L, 0x0.ffffffp0);
+ TEST_ff_f (nexttoward, -1.0, -0x1.000000000000002p0L, -0x1.000002p0);
+ TEST_ff_f (nexttoward, -1.0, -0x0.ffffffffffffffffp0L, -0x0.ffffffp0);
+# endif
+# if LDBL_MANT_DIG >= 106
+ TEST_ff_f (nexttoward, 1.0, 0x1.000000000000000000000000008p0L, 0x1.000002p0);
+ TEST_ff_f (nexttoward, 1.0, 0x0.ffffffffffffffffffffffffffcp0L, 0x0.ffffffp0);
+ TEST_ff_f (nexttoward, -1.0, -0x1.000000000000000000000000008p0L, -0x1.000002p0);
+ TEST_ff_f (nexttoward, -1.0, -0x0.ffffffffffffffffffffffffffcp0L, -0x0.ffffffp0);
+# endif
+# if LDBL_MANT_DIG >= 113
+ TEST_ff_f (nexttoward, 1.0, 0x1.0000000000000000000000000001p0L, 0x1.000002p0);
+ TEST_ff_f (nexttoward, 1.0, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.ffffffp0);
+ TEST_ff_f (nexttoward, -1.0, -0x1.0000000000000000000000000001p0L, -0x1.000002p0);
+ TEST_ff_f (nexttoward, -1.0, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.ffffffp0);
+# endif
+#endif
+#ifdef TEST_DOUBLE
+ TEST_ff_f (nexttoward, 1.0, 1.1L, 0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, 1.0, LDBL_MAX, 0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, 1.0, 0x1.0000000000001p0, 0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, 1.0, 0.9L, 0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, 1.0, -LDBL_MAX, 0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, 1.0, 0x0.fffffffffffff8p0, 0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, -1.0, -1.1L, -0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, -1.0, -LDBL_MAX, -0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, -1.0, -0x1.0000000000001p0, -0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, -1.0, -0.9L, -0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, -1.0, LDBL_MAX, -0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, -1.0, -0x0.fffffffffffff8p0, -0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, -1.0, -0x8.00346dc5d6388p-3L, -0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, 0x1p-1074, 0x1p-1073L, 0x1p-1073);
+# if LDBL_MANT_DIG >= 64
+ TEST_ff_f (nexttoward, 1.0, 0x1.000000000000002p0L, 0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, 1.0, 0x0.ffffffffffffffffp0L, 0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, -1.0, -0x1.000000000000002p0L, -0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, -1.0, -0x0.ffffffffffffffffp0L, -0x0.fffffffffffff8p0);
+# endif
+# if LDBL_MANT_DIG >= 106
+ TEST_ff_f (nexttoward, 1.0, 0x1.000000000000000000000000008p0L, 0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, 1.0, 0x0.ffffffffffffffffffffffffffcp0L, 0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, -1.0, -0x1.000000000000000000000000008p0L, -0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, -1.0, -0x0.ffffffffffffffffffffffffffcp0L, -0x0.fffffffffffff8p0);
+# endif
+# if LDBL_MANT_DIG >= 113
+ TEST_ff_f (nexttoward, 1.0, 0x1.0000000000000000000000000001p0L, 0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, 1.0, 0x0.ffffffffffffffffffffffffffff8p0L, 0x0.fffffffffffff8p0);
+ TEST_ff_f (nexttoward, -1.0, -0x1.0000000000000000000000000001p0L, -0x1.0000000000001p0);
+ TEST_ff_f (nexttoward, -1.0, -0x0.ffffffffffffffffffffffffffff8p0L, -0x0.fffffffffffff8p0);
+# endif
+#endif
END (nexttoward);
}
diff --git a/math/s_nexttowardf.c b/math/s_nexttowardf.c
index fb008d5..e8c4dd1 100644
--- a/math/s_nexttowardf.c
+++ b/math/s_nexttowardf.c
@@ -47,16 +47,12 @@ float __nexttowardf(float x, long double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if(hy<0||(ix>>23)>(iy>>20)-0x380
- || ((ix>>23)==(iy>>20)-0x380
- && (ix&0x7fffff)>(((hy<<3)|(ly>>29))&0x7fffff))) /* x > y, x -= ulp */
+ if(x > y) /* x -= ulp */
hx -= 1;
else /* x < y, x += ulp */
hx += 1;
} else { /* x < 0 */
- if(hy>=0||(ix>>23)>(iy>>20)-0x380
- || ((ix>>23)==(iy>>20)-0x380
- && (ix&0x7fffff)>(((hy<<3)|(ly>>29))&0x7fffff))) /* x < y, x -= ulp */
+ if(x < y) /* x -= ulp */
hx -= 1;
else /* x > y, x += ulp */
hx += 1;
diff --git a/sysdeps/i386/fpu/s_nexttoward.c b/sysdeps/i386/fpu/s_nexttoward.c
index e5f0164..74147c4 100644
--- a/sysdeps/i386/fpu/s_nexttoward.c
+++ b/sysdeps/i386/fpu/s_nexttoward.c
@@ -55,11 +55,7 @@ double __nexttoward(double x, long double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if (esy>=0x8000||((ix>>20)&0x7ff)>iy-0x3c00
- || (((ix>>20)&0x7ff)==iy-0x3c00
- && (((hx<<11)|(lx>>21))>(hy&0x7fffffff)
- || (((hx<<11)|(lx>>21))==(hy&0x7fffffff)
- && (lx<<11)>ly)))) { /* x > y, x -= ulp */
+ if (x > y) { /* x -= ulp */
if(lx==0) hx -= 1;
lx -= 1;
} else { /* x < y, x += ulp */
@@ -67,11 +63,7 @@ double __nexttoward(double x, long double y)
if(lx==0) hx += 1;
}
} else { /* x < 0 */
- if (esy<0x8000||((ix>>20)&0x7ff)>iy-0x3c00
- || (((ix>>20)&0x7ff)==iy-0x3c00
- && (((hx<<11)|(lx>>21))>(hy&0x7fffffff)
- || (((hx<<11)|(lx>>21))==(hy&0x7fffffff)
- && (lx<<11)>ly)))) {/* x < y, x -= ulp */
+ if (x < y) { /* x -= ulp */
if(lx==0) hx -= 1;
lx -= 1;
} else { /* x > y, x += ulp */
diff --git a/sysdeps/i386/fpu/s_nexttowardf.c b/sysdeps/i386/fpu/s_nexttowardf.c
index 89e8771..49651be 100644
--- a/sysdeps/i386/fpu/s_nexttowardf.c
+++ b/sysdeps/i386/fpu/s_nexttowardf.c
@@ -47,17 +47,13 @@ float __nexttowardf(float x, long double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if(esy>=0x8000||((ix>>23)&0xff)>iy-0x3f80
- || (((ix>>23)&0xff)==iy-0x3f80
- && ((ix&0x7fffff)<<8)>(hy&0x7fffffff))) {/* x > y, x -= ulp */
+ if(x > y) { /* x -= ulp */
hx -= 1;
} else { /* x < y, x += ulp */
hx += 1;
}
} else { /* x < 0 */
- if(esy<0x8000||((ix>>23)&0xff)>iy-0x3f80
- || (((ix>>23)&0xff)==iy-0x3f80
- && ((ix&0x7fffff)<<8)>(hy&0x7fffffff))) {/* x < y, x -= ulp */
+ if(x < y) { /* x -= ulp */
hx -= 1;
} else { /* x > y, x += ulp */
hx += 1;
diff --git a/sysdeps/ieee754/ldbl-128/s_nexttoward.c b/sysdeps/ieee754/ldbl-128/s_nexttoward.c
index 87a9a6c..1ea0b64 100644
--- a/sysdeps/ieee754/ldbl-128/s_nexttoward.c
+++ b/sysdeps/ieee754/ldbl-128/s_nexttoward.c
@@ -55,11 +55,7 @@ double __nexttoward(double x, long double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if (hy<0||(ix>>20)>(iy>>48)-0x3c00
- || ((ix>>20)==(iy>>48)-0x3c00
- && (((((int64_t)hx)<<28)|(lx>>4))>(hy&0x0000ffffffffffffLL)
- || (((((int64_t)hx)<<28)|(lx>>4))==(hy&0x0000ffffffffffffLL)
- && (lx&0xf)>(ly>>60))))) { /* x > y, x -= ulp */
+ if (x > y) { /* x -= ulp */
if(lx==0) hx -= 1;
lx -= 1;
} else { /* x < y, x += ulp */
@@ -67,11 +63,7 @@ double __nexttoward(double x, long double y)
if(lx==0) hx += 1;
}
} else { /* x < 0 */
- if (hy>=0||(ix>>20)>(iy>>48)-0x3c00
- || ((ix>>20)==(iy>>48)-0x3c00
- && (((((int64_t)hx)<<28)|(lx>>4))>(hy&0x0000ffffffffffffLL)
- || (((((int64_t)hx)<<28)|(lx>>4))==(hy&0x0000ffffffffffffLL)
- && (lx&0xf)>(ly>>60))))) { /* x < y, x -= ulp */
+ if (x < y) { /* x -= ulp */
if(lx==0) hx -= 1;
lx -= 1;
} else { /* x > y, x += ulp */
diff --git a/sysdeps/ieee754/ldbl-128/s_nexttowardf.c b/sysdeps/ieee754/ldbl-128/s_nexttowardf.c
index 371fa80..02a1407 100644
--- a/sysdeps/ieee754/ldbl-128/s_nexttowardf.c
+++ b/sysdeps/ieee754/ldbl-128/s_nexttowardf.c
@@ -46,17 +46,13 @@ float __nexttowardf(float x, long double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if(hy<0||(ix>>23)>(iy>>48)-0x3f80
- || ((ix>>23)==(iy>>48)-0x3f80
- && (ix&0x7fffff)>((hy>>25)&0x7fffff))) {/* x > y, x -= ulp */
+ if(x > y) { /* x -= ulp */
hx -= 1;
} else { /* x < y, x += ulp */
hx += 1;
}
} else { /* x < 0 */
- if(hy>=0||(ix>>23)>(iy>>48)-0x3f80
- || ((ix>>23)==(iy>>48)-0x3f80
- && (ix&0x7fffff)>((hy>>25)&0x7fffff))) {/* x < y, x -= ulp */
+ if(x < y) { /* x < y, x -= ulp */
hx -= 1;
} else { /* x > y, x += ulp */
hx += 1;
diff --git a/sysdeps/ieee754/ldbl-128ibm/s_nexttowardf.c b/sysdeps/ieee754/ldbl-128ibm/s_nexttowardf.c
index a674583..b387a91 100644
--- a/sysdeps/ieee754/ldbl-128ibm/s_nexttowardf.c
+++ b/sysdeps/ieee754/ldbl-128ibm/s_nexttowardf.c
@@ -49,17 +49,13 @@ float __nexttowardf(float x, long double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if(hy<0||(ix>>23)>(iy>>52)-0x380
- || ((ix>>23)==(iy>>52)-0x380
- && (ix&0x7fffff)>((hy>>29)&0x7fffff))) {/* x > y, x -= ulp */
+ if(x > y) { /* x -= ulp */
hx -= 1;
} else { /* x < y, x += ulp */
hx += 1;
}
} else { /* x < 0 */
- if(hy>=0||(ix>>23)>(iy>>52)-0x380
- || ((ix>>23)==(iy>>52)-0x380
- && (ix&0x7fffff)>((hy>>29)&0x7fffff))) {/* x < y, x -= ulp */
+ if(x < y) { /* x -= ulp */
hx -= 1;
} else { /* x > y, x += ulp */
hx += 1;
diff --git a/sysdeps/ieee754/ldbl-96/s_nexttoward.c b/sysdeps/ieee754/ldbl-96/s_nexttoward.c
index 9b93eca..f7a8b21 100644
--- a/sysdeps/ieee754/ldbl-96/s_nexttoward.c
+++ b/sysdeps/ieee754/ldbl-96/s_nexttoward.c
@@ -52,11 +52,7 @@ double __nexttoward(double x, long double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if (esy>=0x8000||((ix>>20)&0x7ff)>iy-0x3c00
- || (((ix>>20)&0x7ff)==iy-0x3c00
- && (((hx<<11)|(lx>>21))>(hy&0x7fffffff)
- || (((hx<<11)|(lx>>21))==(hy&0x7fffffff)
- && (lx<<11)>ly)))) { /* x > y, x -= ulp */
+ if (x > y) { /* x -= ulp */
if(lx==0) hx -= 1;
lx -= 1;
} else { /* x < y, x += ulp */
@@ -64,11 +60,7 @@ double __nexttoward(double x, long double y)
if(lx==0) hx += 1;
}
} else { /* x < 0 */
- if (esy<0x8000||((ix>>20)&0x7ff)>iy-0x3c00
- || (((ix>>20)&0x7ff)==iy-0x3c00
- && (((hx<<11)|(lx>>21))>(hy&0x7fffffff)
- || (((hx<<11)|(lx>>21))==(hy&0x7fffffff)
- && (lx<<11)>ly)))) {/* x < y, x -= ulp */
+ if (x < y) { /* x -= ulp */
if(lx==0) hx -= 1;
lx -= 1;
} else { /* x > y, x += ulp */
diff --git a/sysdeps/ieee754/ldbl-96/s_nexttowardf.c b/sysdeps/ieee754/ldbl-96/s_nexttowardf.c
index aeb92b6..a96f9da 100644
--- a/sysdeps/ieee754/ldbl-96/s_nexttowardf.c
+++ b/sysdeps/ieee754/ldbl-96/s_nexttowardf.c
@@ -44,17 +44,13 @@ float __nexttowardf(float x, long double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if(esy>=0x8000||((ix>>23)&0xff)>iy-0x3f80
- || (((ix>>23)&0xff)==iy-0x3f80
- && ((ix&0x7fffff)<<8)>(hy&0x7fffffff))) {/* x > y, x -= ulp */
+ if(x > y) { /* x -= ulp */
hx -= 1;
} else { /* x < y, x += ulp */
hx += 1;
}
} else { /* x < 0 */
- if(esy<0x8000||((ix>>23)&0xff)>iy-0x3f80
- || (((ix>>23)&0xff)==iy-0x3f80
- && ((ix&0x7fffff)<<8)>(hy&0x7fffffff))) {/* x < y, x -= ulp */
+ if(x < y) { /* x -= ulp */
hx -= 1;
} else { /* x > y, x += ulp */
hx += 1;
diff --git a/sysdeps/ieee754/ldbl-opt/s_nexttowardfd.c b/sysdeps/ieee754/ldbl-opt/s_nexttowardfd.c
index 68027f2..7eca121 100644
--- a/sysdeps/ieee754/ldbl-opt/s_nexttowardfd.c
+++ b/sysdeps/ieee754/ldbl-opt/s_nexttowardfd.c
@@ -50,16 +50,12 @@ float __nldbl_nexttowardf(float x, double y)
return x;
}
if(hx>=0) { /* x > 0 */
- if(hy<0||(ix>>23)>(iy>>20)-0x380
- || ((ix>>23)==(iy>>20)-0x380
- && (ix&0x7fffff)>(((hy<<3)|(ly>>29))&0x7fffff))) /* x > y, x -= ulp */
+ if(x > y) /* x -= ulp */
hx -= 1;
else /* x < y, x += ulp */
hx += 1;
} else { /* x < 0 */
- if(hy>=0||(ix>>23)>(iy>>20)-0x380
- || ((ix>>23)==(iy>>20)-0x380
- && (ix&0x7fffff)>(((hy<<3)|(ly>>29))&0x7fffff))) /* x < y, x -= ulp */
+ if(x < y) /* x -= ulp */
hx -= 1;
else /* x > y, x += ulp */
hx += 1;
--
Joseph S. Myers
joseph@codesourcery.com