On Linux/mips (all 3 ABIs) and on Linux/hppa (configuring a testdir in such
a way that the configuration ignores the totalorder* functions provided by
glibc and instead uses the gnulib replacement instead:
gl_cv_func_totalorder_in_libm=no gl_cv_func_totalorder_no_libm=no \
gl_cv_func_totalorderf_in_libm=no gl_cv_func_totalorderf_no_libm=no \
gl_cv_func_totalorderl_in_libm=no gl_cv_func_totalorderl_no_libm=no \
./configure
), I see these test failures:
FAIL: test-totalorder
=====================
Failed: i=0 j=1
Failed: i=1 j=0
Failed: i=12 j=13
Failed: i=13 j=12
FAIL test-totalorder (exit status: 1)
FAIL: test-totalorderf
======================
Failed: i=0 j=1
Failed: i=1 j=0
Failed: i=12 j=13
Failed: i=13 j=12
FAIL test-totalorderf (exit status: 1)
FAIL: test-totalorderl
======================
Failed: i=0 j=1
Failed: i=1 j=0
Failed: i=12 j=13
Failed: i=13 j=12
FAIL test-totalorderl (exit status: 1)
That is, QNaN and SNaN produce a wrong comparison result. See the explanation
at the top of lib/snan.h.
This patch fixes the totalorder* replacements, to work right according to the
specification.
Possibly, on the mips side, the code should test __mips_nan2008. But so far
I have not encountered any mips machine where this would be necessary.
2023-10-14 Bruno Haible <br...@clisp.org>
totalorder*: Fix test failures on PA-RISC and MIPS CPUs.
* lib/totalorderf.c (totalorderf): On hppa and mips, invert bit 22
before comparing two NaNs.
* lib/totalorder.c (totalorder): On hppa and mips, invert bit 51 before
comparing two NaNs.
* lib/totalorderl.c: Include <float.h>.
(totalorderl): On hppa and mips, invert bit 51 or 47 of the xhi, yhi
parts before comparing two NaNs.
* modules/totalorderl (Depends-on): Add 'float'.
diff --git a/lib/totalorder.c b/lib/totalorder.c
index 4c861f971c..edbcd6998b 100644
--- a/lib/totalorder.c
+++ b/lib/totalorder.c
@@ -18,26 +18,38 @@
#include <config.h>
+/* Specification. */
#include <math.h>
int
totalorder (double const *x, double const *y)
{
- /* Although the following is not strictly portable, and won't work
- on some obsolete platforms (e.g., PA-RISC, MIPS before alignment
- to IEEE 754-2008), it should be good enough nowadays. */
+ /* If the sign bits of *X and *Y differ, the one with non-zero sign bit
+ is "smaller" than the one with sign bit == 0. */
int xs = signbit (*x);
int ys = signbit (*y);
if (!xs != !ys)
return xs;
+
+ /* If one of *X, *Y is a NaN and the other isn't, the answer is easy
+ as well: the negative NaN is "smaller", the positive NaN is "greater"
+ than the other argument. */
int xn = isnand (*x);
int yn = isnand (*y);
if (!xn != !yn)
return !xn == !xs;
+ /* If none of *X, *Y is a NaN, the '<=' operator does the job, including
+ for -Infinity and +Infinity. */
if (!xn)
return *x <= *y;
+ /* At this point, *X and *Y are NaNs with the same sign bit. */
+
unsigned long long extended_sign = -!!xs;
+#if defined __hppa || defined __mips__
+ /* Invert the most significant bit of the mantissa field. Cf. snan.h. */
+ extended_sign ^= (1ULL << 51);
+#endif
union { unsigned long long i; double f; } volatile xu = {0}, yu = {0};
xu.f = *x;
yu.f = *y;
diff --git a/lib/totalorderf.c b/lib/totalorderf.c
index d86299ff09..4cbfa06ffd 100644
--- a/lib/totalorderf.c
+++ b/lib/totalorderf.c
@@ -18,27 +18,39 @@
#include <config.h>
+/* Specification. */
#include <math.h>
int
totalorderf (float const *x, float const *y)
{
- /* Although the following is not strictly portable, and won't work
- on some obsolete platforms (e.g., PA-RISC, MIPS before alignment
- to IEEE 754-2008), it should be good enough nowadays. */
+ /* If the sign bits of *X and *Y differ, the one with non-zero sign bit
+ is "smaller" than the one with sign bit == 0. */
int xs = signbit (*x);
int ys = signbit (*y);
if (!xs != !ys)
return xs;
+
+ /* If one of *X, *Y is a NaN and the other isn't, the answer is easy
+ as well: the negative NaN is "smaller", the positive NaN is "greater"
+ than the other argument. */
int xn = isnanf (*x);
int yn = isnanf (*y);
if (!xn != !yn)
return !xn == !xs;
+ /* If none of *X, *Y is a NaN, the '<=' operator does the job, including
+ for -Infinity and +Infinity. */
if (!xn)
return *x <= *y;
- unsigned long extended_sign = -!!xs;
- union { unsigned long i; float f; } volatile xu = {0}, yu = {0};
+ /* At this point, *X and *Y are NaNs with the same sign bit. */
+
+ unsigned int extended_sign = -!!xs;
+#if defined __hppa || defined __mips__
+ /* Invert the most significant bit of the mantissa field. Cf. snan.h. */
+ extended_sign ^= (1U << 22);
+#endif
+ union { unsigned int i; float f; } volatile xu = {0}, yu = {0};
xu.f = *x;
yu.f = *y;
return (xu.i ^ extended_sign) <= (yu.i ^ extended_sign);
diff --git a/lib/totalorderl.c b/lib/totalorderl.c
index cc0698febe..6982762a0f 100644
--- a/lib/totalorderl.c
+++ b/lib/totalorderl.c
@@ -18,8 +18,11 @@
#include <config.h>
+/* Specification. */
#include <math.h>
+#include <float.h>
+
#ifndef LDBL_SIGNBIT_WORD
# define LDBL_SIGNBIT_WORD (-1)
#endif
@@ -27,10 +30,6 @@
int
totalorderl (long double const *x, long double const *y)
{
- /* Although the following is not strictly portable, and won't work
- on some obsolete platforms (e.g., PA-RISC, MIPS before alignment
- to IEEE 754-2008), it should be good enough nowadays. */
-
/* If the sign bits of *X and *Y differ, the one with non-zero sign bit
is "smaller" than the one with sign bit == 0. */
int xs = signbit (*x);
@@ -56,6 +55,10 @@ totalorderl (long double const *x, long double const *y)
if (sizeof (long double) <= sizeof (unsigned long long))
{
+#if defined __hppa || defined __mips__
+ /* Invert the most significant bit of the mantissa field. Cf. snan.h.
*/
+ extended_sign ^= (1ULL << 51);
+#endif
union { unsigned long long i; long double f; } volatile
xu = {0}, yu = {0};
xu.f = *x;
@@ -63,6 +66,14 @@ totalorderl (long double const *x, long double const *y)
return (xu.i ^ extended_sign) <= (yu.i ^ extended_sign);
}
+ unsigned long long extended_sign_hi = extended_sign;
+#if defined __hppa || defined __mips__
+ /* Invert the most significant bit of the mantissa field. Cf. snan.h. */
+ extended_sign_hi ^=
+ (1ULL << (LDBL_MANT_DIG == 106
+ ? 51 /* double-double representation */
+ : (LDBL_MANT_DIG - 2) - 64)); /* quad precision representation */
+#endif
union u { unsigned long long i[2]; long double f; } volatile xu, yu;
/* Although it is tempting to initialize with {0}, Solaris cc (Sun C 5.8)
on x86_64 miscompiles {0}: it initializes only the lower 80 bits,
@@ -87,8 +98,8 @@ totalorderl (long double const *x, long double const *y)
bigendian = LDBL_SIGNBIT_WORD < sizeof xu.i[0] / sizeof (unsigned);
unsigned long long
- xhi = xu.i[!bigendian] ^ extended_sign,
- yhi = yu.i[!bigendian] ^ extended_sign,
+ xhi = xu.i[!bigendian] ^ extended_sign_hi,
+ yhi = yu.i[!bigendian] ^ extended_sign_hi,
xlo = xu.i[ bigendian] ^ extended_sign,
ylo = yu.i[ bigendian] ^ extended_sign;
return (xhi < yhi) | ((xhi == yhi) & (xlo <= ylo));
diff --git a/modules/totalorderl b/modules/totalorderl
index ed6f801130..76dafb8ecd 100644
--- a/modules/totalorderl
+++ b/modules/totalorderl
@@ -10,6 +10,7 @@ m4/signbit.m4
Depends-on:
math
extensions
+float [test $HAVE_TOTALORDERL = 0 || test $REPLACE_TOTALORDERL = 1]
stdbool [test $HAVE_TOTALORDERL = 0 || test $REPLACE_TOTALORDERL = 1]
isnanl [test $HAVE_TOTALORDERL = 0 || test $REPLACE_TOTALORDERL = 1]
signbit [test $HAVE_TOTALORDERL = 0 || test $REPLACE_TOTALORDERL = 1]
