isnanl: recognize non-IEEE numbers

[Bruno Haible]

As noticed in the thread with Andreas Schwab a few days ago, isnanl()
should better recognize numbers that are not operating according to what
IEEE 754 mandates. On IA-64, some of these operate like signalling NaNs,
but there is also IA-32 and AMD64, which have the same 'long double' format.

Note that on HP-UX for IA-64, the processor is in big-endian mode.


2007-06-05  Bruno Haible  <address@hidden>

        Fix isnanl so that it recognizes non-IEEE numbers on i386, x86_64, ia64.
        * m4/isnanl.c (gl_FUNC_ISNANL_WORKS): Require AC_C_BIGENDIAN. Tested
        non-IEEE numbers on i386, x86_64, ia64.
        (gl_LONG_DOUBLE_EXPONENT_LOCATION): Require AC_C_BIGENDIAN.
        * lib/isnan.c (FUNC): Add special code for i386, x86_64, ia64.
        * tests/test-isnanl.h: Include float.h.
        (main): Check also non-IEEE numbers on i386, x86_64, ia64.

*** m4/isnanl.m4        6 Apr 2007 14:36:56 -0000       1.6
--- m4/isnanl.m4        6 Jun 2007 01:46:26 -0000
***************
*** 1,4 ****
! # isnanl.m4 serial 3
  dnl Copyright (C) 2007 Free Software Foundation, Inc.
  dnl This file is free software; the Free Software Foundation
  dnl gives unlimited permission to copy and/or distribute it,
--- 1,4 ----
! # isnanl.m4 serial 4
  dnl Copyright (C) 2007 Free Software Foundation, Inc.
  dnl This file is free software; the Free Software Foundation
  dnl gives unlimited permission to copy and/or distribute it,
***************
*** 97,110 ****
  ])
  
  dnl Test whether isnanl() recognizes all numbers which are neither finite nor
! dnl infinite. This test fails e.g. on NetBSD/i386.
  AC_DEFUN([gl_FUNC_ISNANL_WORKS],
  [
    AC_REQUIRE([AC_PROG_CC])
    AC_REQUIRE([AC_CANONICAL_HOST]) dnl for cross-compiles
    AC_CACHE_CHECK([whether isnanl works], [gl_cv_func_isnanl_works],
      [
        AC_TRY_RUN([
  #include <limits.h>
  #include <math.h>
  #ifdef isnan
--- 97,112 ----
  ])
  
  dnl Test whether isnanl() recognizes all numbers which are neither finite nor
! dnl infinite. This test fails e.g. on NetBSD/i386 and on glibc/ia64.
  AC_DEFUN([gl_FUNC_ISNANL_WORKS],
  [
    AC_REQUIRE([AC_PROG_CC])
+   AC_REQUIRE([AC_C_BIGENDIAN])
    AC_REQUIRE([AC_CANONICAL_HOST]) dnl for cross-compiles
    AC_CACHE_CHECK([whether isnanl works], [gl_cv_func_isnanl_works],
      [
        AC_TRY_RUN([
+ #include <float.h>
  #include <limits.h>
  #include <math.h>
  #ifdef isnan
***************
*** 113,119 ****
  #endif
  #define NWORDS \
    ((sizeof (long double) + sizeof (unsigned int) - 1) / sizeof (unsigned int))
! typedef union { long double value; unsigned int word[NWORDS]; }
          memory_long_double;
  int main ()
  {
--- 115,121 ----
  #endif
  #define NWORDS \
    ((sizeof (long double) + sizeof (unsigned int) - 1) / sizeof (unsigned int))
! typedef union { unsigned int word[NWORDS]; long double value; }
          memory_long_double;
  int main ()
  {
***************
*** 131,136 ****
--- 133,203 ----
    if (!isnanl (m.value))
      return 1;
  
+ #if ((defined __ia64 && LDBL_MANT_DIG == 64) || (defined __x86_64__ || 
defined __amd64__) || (defined __i386 || defined __i386__ || defined _I386 || 
defined _M_IX86 || defined _X86_))
+ /* Representation of an 80-bit 'long double' as an initializer for a sequence
+    of 'unsigned int' words.  */
+ # ifdef WORDS_BIGENDIAN
+ #  define LDBL80_WORDS(exponent,manthi,mantlo) \
+      { ((unsigned int) (exponent) << 16) | ((unsigned int) (manthi) >> 16), \
+        ((unsigned int) (manthi) << 16) | (unsigned int) (mantlo) >> 16),    \
+        (unsigned int) (mantlo) << 16                                        \
+      }
+ # else
+ #  define LDBL80_WORDS(exponent,manthi,mantlo) \
+      { mantlo, manthi, exponent }
+ # endif
+   { /* Quiet NaN.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0xFFFF, 0xC3333333, 0x00000000) };
+     if (!isnanl (x.value))
+       return 1;
+   }
+   {
+     /* Signalling NaN.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0xFFFF, 0x83333333, 0x00000000) };
+     if (!isnanl (x.value))
+       return 1;
+   }
+   /* The isnanl function should recognize Pseudo-NaNs, Pseudo-Infinities,
+      Pseudo-Zeroes, Unnormalized Numbers, and Pseudo-Denormals, as defined in
+        Intel IA-64 Architecture Software Developer's Manual, Volume 1:
+        Application Architecture.
+        Table 5-2 "Floating-Point Register Encodings"
+        Figure 5-6 "Memory to Floating-Point Register Data Translation"
+    */
+   { /* Pseudo-NaN.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0xFFFF, 0x40000001, 0x00000000) };
+     if (!isnanl (x.value))
+       return 1;
+   }
+   { /* Pseudo-Infinity.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0xFFFF, 0x00000000, 0x00000000) };
+     if (!isnanl (x.value))
+       return 1;
+   }
+   { /* Pseudo-Zero.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0x4004, 0x00000000, 0x00000000) };
+     if (!isnanl (x.value))
+       return 1;
+   }
+   { /* Unnormalized number.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0x4000, 0x63333333, 0x00000000) };
+     if (!isnanl (x.value))
+       return 1;
+   }
+   { /* Pseudo-Denormal.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0x0000, 0x83333333, 0x00000000) };
+     if (!isnanl (x.value))
+       return 1;
+   }
+ #endif
+ 
    return 0;
  }], [gl_cv_func_isnanl_works=yes], [gl_cv_func_isnanl_works=no],
        [case "$host_os" in
***************
*** 143,148 ****
--- 210,216 ----
  
  AC_DEFUN([gl_LONG_DOUBLE_EXPONENT_LOCATION],
  [
+   AC_REQUIRE([AC_C_BIGENDIAN])
    AC_CACHE_CHECK([where to find the exponent in a 'long double'],
      [gl_cv_cc_long_double_expbit0],
      [
*** lib/isnan.c 19 May 2007 23:54:48 -0000      1.7
--- lib/isnan.c 6 Jun 2007 01:46:26 -0000
***************
*** 72,80 ****
  FUNC (DOUBLE x)
  {
  #ifdef KNOWN_EXPBIT0_LOCATION
    /* Be careful to not do any floating-point operation on x, such as x == x,
       because x may be a signaling NaN.  */
! # if defined __SUNPRO_C || defined __DECC || (defined __sgi && !defined 
__GNUC__)
    /* The Sun C 5.0 compilers and the Compaq (ex-DEC) 6.4 compilers don't
       recognize the initializers as constant expressions.  The latter compiler
       also fails when constant-folding 0.0 / 0.0 even when constant-folding is
--- 72,112 ----
  FUNC (DOUBLE x)
  {
  #ifdef KNOWN_EXPBIT0_LOCATION
+ # if defined USE_LONG_DOUBLE && ((defined __ia64 && LDBL_MANT_DIG == 64) || 
(defined __x86_64__ || defined __amd64__) || (defined __i386 || defined 
__i386__ || defined _I386 || defined _M_IX86 || defined _X86_))
+   /* Special CPU dependent code is needed to treat bit patterns outside the
+      IEEE 754 specification (such as Pseudo-NaNs, Pseudo-Infinities,
+      Pseudo-Zeroes, Unnormalized Numbers, and Pseudo-Denormals) as NaNs.
+      These bit patterns are:
+        - exponent = 0x0001..0x7FFF, mantissa bit 63 = 0,
+        - exponent = 0x0000, mantissa bit 63 = 1.
+      The NaN bit pattern is:
+        - exponent = 0x7FFF, mantissa >= 0x8000000000000001.  */
+   memory_double m;
+   unsigned int exponent;
+ 
+   m.value = x;
+   exponent = (m.word[EXPBIT0_WORD] >> EXPBIT0_BIT) & EXP_MASK;
+ #  ifdef WORDS_BIGENDIAN
+   /* Big endian: EXPBIT0_WORD = 0, EXPBIT0_BIT = 16.  */
+   if (exponent == 0)
+     return 1 & (m.word[0] >> 15);
+   else if (exponent == EXP_MASK)
+     return (((m.word[0] ^ 0x8000U) << 16) | m.word[1] | (m.word[2] >> 16)) != 
0;
+   else
+     return 1 & ~(m.word[0] >> 15);
+ #  else
+   /* Little endian: EXPBIT0_WORD = 2, EXPBIT0_BIT = 0.  */
+   if (exponent == 0)
+     return (m.word[1] >> 31);
+   else if (exponent == EXP_MASK)
+     return ((m.word[1] ^ 0x80000000U) | m.word[0]) != 0;
+   else
+     return (m.word[1] >> 31) ^ 1;
+ #  endif
+ # else
    /* Be careful to not do any floating-point operation on x, such as x == x,
       because x may be a signaling NaN.  */
! #  if defined __SUNPRO_C || defined __DECC || (defined __sgi && !defined 
__GNUC__)
    /* The Sun C 5.0 compilers and the Compaq (ex-DEC) 6.4 compilers don't
       recognize the initializers as constant expressions.  The latter compiler
       also fails when constant-folding 0.0 / 0.0 even when constant-folding is
***************
*** 85,95 ****
    DOUBLE plus_inf = L_(1.0) / L_(0.0);
    DOUBLE minus_inf = -L_(1.0) / L_(0.0);
    nan.value = zero / zero;
! # else
    static memory_double nan = { L_(0.0) / L_(0.0) };
    static DOUBLE plus_inf = L_(1.0) / L_(0.0);
    static DOUBLE minus_inf = -L_(1.0) / L_(0.0);
! # endif
    {
      memory_double m;
  
--- 117,127 ----
    DOUBLE plus_inf = L_(1.0) / L_(0.0);
    DOUBLE minus_inf = -L_(1.0) / L_(0.0);
    nan.value = zero / zero;
! #  else
    static memory_double nan = { L_(0.0) / L_(0.0) };
    static DOUBLE plus_inf = L_(1.0) / L_(0.0);
    static DOUBLE minus_inf = -L_(1.0) / L_(0.0);
! #  endif
    {
      memory_double m;
  
***************
*** 104,109 ****
--- 136,142 ----
      else
        return 0;
    }
+ # endif
  #else
    /* The configuration did not find sufficient information.  Give up about
       the signaling NaNs, handle only the quiet NaNs.  */
*** tests/test-isnanl.h 29 Apr 2007 09:15:13 -0000      1.2
--- tests/test-isnanl.h 6 Jun 2007 01:46:27 -0000
***************
*** 17,22 ****
--- 17,23 ----
  
  /* Written by Bruno Haible <address@hidden>, 2007.  */
  
+ #include <float.h>
  #include <limits.h>
  #include <stdio.h>
  #include <stdlib.h>
***************
*** 35,40 ****
--- 36,46 ----
  int
  main ()
  {
+   #define NWORDS \
+     ((sizeof (long double) + sizeof (unsigned int) - 1) / sizeof (unsigned 
int))
+   typedef union { unsigned int word[NWORDS]; long double value; }
+           memory_long_double;
+ 
    /* Finite values.  */
    ASSERT (!isnanl (3.141L));
    ASSERT (!isnanl (3.141e30L));
***************
*** 47,59 ****
    ASSERT (!isnanl (-1.0L / 0.0L));
    /* Quiet NaN.  */
    ASSERT (isnanl (0.0L / 0.0L));
  #if defined LDBL_EXPBIT0_WORD && defined LDBL_EXPBIT0_BIT
!   /* Signalling NaN.  */
    {
-     #define NWORDS \
-       ((sizeof (long double) + sizeof (unsigned int) - 1) / sizeof (unsigned 
int))
-     typedef union { long double value; unsigned int word[NWORDS]; }
-             memory_long_double;
      memory_long_double m;
      m.value = 0.0L / 0.0L;
  # if LDBL_EXPBIT0_BIT > 0
--- 53,63 ----
    ASSERT (!isnanl (-1.0L / 0.0L));
    /* Quiet NaN.  */
    ASSERT (isnanl (0.0L / 0.0L));
+ 
  #if defined LDBL_EXPBIT0_WORD && defined LDBL_EXPBIT0_BIT
!   /* A bit pattern that is different from a Quiet NaN.  With a bit of luck,
!      it's a Signalling NaN.  */
    {
      memory_long_double m;
      m.value = 0.0L / 0.0L;
  # if LDBL_EXPBIT0_BIT > 0
***************
*** 67,71 ****
--- 71,134 ----
      ASSERT (isnanl (m.value));
    }
  #endif
+ 
+ #if ((defined __ia64 && LDBL_MANT_DIG == 64) || (defined __x86_64__ || 
defined __amd64__) || (defined __i386 || defined __i386__ || defined _I386 || 
defined _M_IX86 || defined _X86_))
+ /* Representation of an 80-bit 'long double' as an initializer for a sequence
+    of 'unsigned int' words.  */
+ # ifdef WORDS_BIGENDIAN
+ #  define LDBL80_WORDS(exponent,manthi,mantlo) \
+      { ((unsigned int) (exponent) << 16) | ((unsigned int) (manthi) >> 16), \
+        ((unsigned int) (manthi) << 16) | (unsigned int) (mantlo) >> 16),    \
+        (unsigned int) (mantlo) << 16                                        \
+      }
+ # else
+ #  define LDBL80_WORDS(exponent,manthi,mantlo) \
+      { mantlo, manthi, exponent }
+ # endif
+   { /* Quiet NaN.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0xFFFF, 0xC3333333, 0x00000000) };
+     ASSERT (isnanl (x.value));
+   }
+   {
+     /* Signalling NaN.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0xFFFF, 0x83333333, 0x00000000) };
+     ASSERT (isnanl (x.value));
+   }
+   /* The isnanl function should recognize Pseudo-NaNs, Pseudo-Infinities,
+      Pseudo-Zeroes, Unnormalized Numbers, and Pseudo-Denormals, as defined in
+        Intel IA-64 Architecture Software Developer's Manual, Volume 1:
+        Application Architecture.
+        Table 5-2 "Floating-Point Register Encodings"
+        Figure 5-6 "Memory to Floating-Point Register Data Translation"
+    */
+   { /* Pseudo-NaN.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0xFFFF, 0x40000001, 0x00000000) };
+     ASSERT (isnanl (x.value));
+   }
+   { /* Pseudo-Infinity.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0xFFFF, 0x00000000, 0x00000000) };
+     ASSERT (isnanl (x.value));
+   }
+   { /* Pseudo-Zero.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0x4004, 0x00000000, 0x00000000) };
+     ASSERT (isnanl (x.value));
+   }
+   { /* Unnormalized number.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0x4000, 0x63333333, 0x00000000) };
+     ASSERT (isnanl (x.value));
+   }
+   { /* Pseudo-Denormal.  */
+     static memory_long_double x =
+       { LDBL80_WORDS (0x0000, 0x83333333, 0x00000000) };
+     ASSERT (isnanl (x.value));
+   }
+ #endif
+ 
    return 0;
  }