[svn-r7875] Purpose:

    Omnibus floating-point bug fix changes

Description:
    There are a number of problems in the floating-point conversion code that
were exposed by Ray's recent int<->float checkin:
    - The 'my_isnan' code in test/dtypes.c was broken and would always return
        true.  The meant that the actual values in the float<->float conversion
        tests were _never_ checked, hiding the other bugs included in this
        checkin.
    - A recent change I made to the type conversion code used "FLT_MIN" instead
        of "-FLT_MAX" for the most negative 'float' value for the double->float
        conversion, which meant that any the negative number that was converted
        from a double to a float would have been mapped to zero, essentially.
    - A change that Robb appeared to have made ~2.5 years ago to the "generic"
        float->float conversion routine appears to be incorrect and I've backed
        it out.
    - Floating-point conversions on SGI's which converted denormalized values
        would be mapped to zero instead of being propertly preserved in the new
        type.  This was addressed by an SGI-specific system call to prevent the
        behavior.

Solution:
    Described above, generally.

Platforms tested:
    FreeBSD 4.9 (sleipnir)
    h5committest

Misc. update:
    release_docs/RELEASE update forthcoming...

4 files changed, 210 insertions, 2 deletions

diff --git a/src/H5T.c b/src/H5T.c
index 0565a4c..e1e463c 100644
--- a/src/H5T.c
+++ b/src/H5T.c
@@ -35,6 +35,11 @@
 #include "H5Pprivate.h"		/* Property Lists			  */
 #include "H5Tpkg.h"		/*data-type functions			  */
 
+/* Check for header needed for SGI floating-point code */
+#ifdef H5_HAVE_SYS_FPU_H
+#include <sys/fpu.h>
+#endif /* H5_HAVE_SYS_FPU_H */
+
 /* Interface initialization */
 static int interface_initialize_g = 0;
 #define INTERFACE_INIT H5T_init_interface
@@ -215,6 +220,13 @@ size_t H5T_NATIVE_UINT_LEAST64_ALIGN_g	= 0;
 size_t H5T_NATIVE_INT_FAST64_ALIGN_g	= 0;
 size_t H5T_NATIVE_UINT_FAST64_ALIGN_g	= 0;
 
+/* Useful floating-point values for conversion routines */
+/* (+/- Inf for all floating-point types) */
+float H5T_NATIVE_FLOAT_POS_INF_g        = 0.0;
+float H5T_NATIVE_FLOAT_NEG_INF_g        = 0.0;
+double H5T_NATIVE_DOUBLE_POS_INF_g      = 0.0;
+double H5T_NATIVE_DOUBLE_NEG_INF_g      = 0.0;
+
 
 /*
  * The path database. Each path has a source and destination data type pair
@@ -467,6 +479,165 @@ done:
 }
 
 
+/*-------------------------------------------------------------------------
+ * Function:	H5T_init_inf
+ *
+ * Purpose:	Initialize the +/- Infinity floating-poing values for type
+ *              conversion.
+ *
+ * Return:	Success:	non-negative
+ *
+ *		Failure:	negative
+ *
+ * Programmer:	Quincey Koziol
+ *              Saturday, November 22, 2003
+ *
+ * Modifications:
+ *
+ *-------------------------------------------------------------------------
+ */
+static herr_t
+H5T_init_inf(void)
+{
+    H5T_t	*dst_p;		/* Datatype type operate on */
+    H5T_atomic_t *dst;		/* Datatype's atomic info */
+    uint8_t	*d;             /* Pointer to value to set */
+    size_t	half_size;	/* Half the type size */
+    size_t      u;              /* Local index value */
+    herr_t ret_value=SUCCEED;   /* Return value */
+
+    FUNC_ENTER_NOINIT(H5T_init_inf);
+
+    /* Get the float datatype */
+    if (NULL==(dst_p=H5I_object(H5T_NATIVE_FLOAT_g)))
+        HGOTO_ERROR (H5E_ARGS, H5E_BADTYPE, FAIL, "not a datatype");
+    dst = &dst_p->u.atomic;
+
+    /* Check that we can re-order the bytes correctly */
+    if (H5T_ORDER_LE!=dst->order && H5T_ORDER_BE!=dst->order)
+        HGOTO_ERROR (H5E_DATATYPE, H5E_UNSUPPORTED, FAIL, "unsupported byte order");
+
+    /* +Inf */
+    d=(uint8_t *)&H5T_NATIVE_FLOAT_POS_INF_g;
+    H5T_bit_set (d, dst->u.f.sign, 1, FALSE);
+    H5T_bit_set (d, dst->u.f.epos, dst->u.f.esize, TRUE);
+    H5T_bit_set (d, dst->u.f.mpos, dst->u.f.msize, FALSE);
+
+    /* Swap the bytes if the machine architecture is big-endian */
+    if (H5T_ORDER_BE==dst->order) {
+        half_size = dst_p->size/2;
+        for (u=0; u<half_size; u++) {
+            uint8_t tmp = d[dst_p->size-(u+1)];
+            d[dst_p->size-(u+1)] = d[u];
+            d[u] = tmp;
+        }
+    }
+
+    /* -Inf */
+    d=(uint8_t *)&H5T_NATIVE_FLOAT_NEG_INF_g;
+    H5T_bit_set (d, dst->u.f.sign, 1, TRUE);
+    H5T_bit_set (d, dst->u.f.epos, dst->u.f.esize, TRUE);
+    H5T_bit_set (d, dst->u.f.mpos, dst->u.f.msize, FALSE);
+
+    /* Swap the bytes if the machine architecture is big-endian */
+    if (H5T_ORDER_BE==dst->order) {
+        half_size = dst_p->size/2;
+        for (u=0; u<half_size; u++) {
+            uint8_t tmp = d[dst_p->size-(u+1)];
+            d[dst_p->size-(u+1)] = d[u];
+            d[u] = tmp;
+        }
+    }
+
+    /* Get the double datatype */
+    if (NULL==(dst_p=H5I_object(H5T_NATIVE_DOUBLE_g)))
+        HGOTO_ERROR (H5E_ARGS, H5E_BADTYPE, FAIL, "not a datatype");
+    dst = &dst_p->u.atomic;
+
+    /* Check that we can re-order the bytes correctly */
+    if (H5T_ORDER_LE!=dst->order && H5T_ORDER_BE!=dst->order)
+        HGOTO_ERROR (H5E_DATATYPE, H5E_UNSUPPORTED, FAIL, "unsupported byte order");
+
+    /* +Inf */
+    d=(uint8_t *)&H5T_NATIVE_DOUBLE_POS_INF_g;
+    H5T_bit_set (d, dst->u.f.sign, 1, FALSE);
+    H5T_bit_set (d, dst->u.f.epos, dst->u.f.esize, TRUE);
+    H5T_bit_set (d, dst->u.f.mpos, dst->u.f.msize, FALSE);
+
+    /* Swap the bytes if the machine architecture is big-endian */
+    if (H5T_ORDER_BE==dst->order) {
+        half_size = dst_p->size/2;
+        for (u=0; u<half_size; u++) {
+            uint8_t tmp = d[dst_p->size-(u+1)];
+            d[dst_p->size-(u+1)] = d[u];
+            d[u] = tmp;
+        }
+    }
+
+    /* -Inf */
+    d=(uint8_t *)&H5T_NATIVE_DOUBLE_NEG_INF_g;
+    H5T_bit_set (d, dst->u.f.sign, 1, TRUE);
+    H5T_bit_set (d, dst->u.f.epos, dst->u.f.esize, TRUE);
+    H5T_bit_set (d, dst->u.f.mpos, dst->u.f.msize, FALSE);
+
+    /* Swap the bytes if the machine architecture is big-endian */
+    if (H5T_ORDER_BE==dst->order) {
+        half_size = dst_p->size/2;
+        for (u=0; u<half_size; u++) {
+            uint8_t tmp = d[dst_p->size-(u+1)];
+            d[dst_p->size-(u+1)] = d[u];
+            d[u] = tmp;
+        }
+    }
+
+done:
+    FUNC_LEAVE_NOAPI(ret_value);
+}
+
+
+/*-------------------------------------------------------------------------
+ * Function:	H5T_init_hw
+ *
+ * Purpose:	Perform hardware specific [floating-point] initialization
+ *
+ * Return:	Success:	non-negative
+ *
+ *		Failure:	negative
+ *
+ * Programmer:	Quincey Koziol
+ *              Monday, November 24, 2003
+ *
+ * Modifications:
+ *
+ *-------------------------------------------------------------------------
+ */
+static herr_t
+H5T_init_hw(void)
+{
+#ifdef H5_HAVE_GET_FPC_CSR
+    union fpc_csr csr;          /* Union to hold results of floating-point status register query */
+#endif /* H5_HAVE_GET_FPC_CSR */
+    herr_t ret_value=SUCCEED;   /* Return value */
+
+    FUNC_ENTER_NOINIT(H5T_init_hw);
+
+#ifdef H5_HAVE_GET_FPC_CSR
+    /* [This code is specific to SGI machines] */
+
+    /* Get the floating-point status register */
+    csr.fc_word=get_fpc_csr();
+
+    /* If the "flush denormalized values to zero" flag is set, unset it */
+    if(csr.fc_struct.flush) {
+        csr.fc_struct.flush=0;
+        set_fpc_csr(csr.fc_word);
+    } /* end if */
+#endif /* H5_HAVE_GET_FPC_CSR */
+
+    FUNC_LEAVE_NOAPI(ret_value);
+}
+
+
 /*--------------------------------------------------------------------------
 NAME
    H5T_init_interface -- Initialize interface-specific information
@@ -529,6 +700,10 @@ H5T_init_interface(void)
     /* Only 16 (numbered 0-15) are supported in the current file format */
     assert(H5T_NCLASSES<16);
 
+    /* Perform any necessary hardware initializations */
+    if(H5T_init_hw()<0)
+	HGOTO_ERROR (H5E_DATATYPE, H5E_CANTINIT, FAIL, "unable to initialize interface");
+
     /*
      * Initialize pre-defined native data types from code generated during
      * the library configuration by H5detect.
@@ -993,6 +1168,9 @@ H5T_init_interface(void)
      */
     status |= H5T_register(H5T_PERS_HARD, "no-op", native_int, native_int, H5T_conv_noop, H5AC_dxpl_id);
 
+    /* Initialize the +/- Infinity values for floating-point types */
+    status |= H5T_init_inf();
+
     if (status<0)
 	HGOTO_ERROR(H5E_DATATYPE, H5E_CANTINIT, FAIL, "unable to register conversion function(s)");
     
diff --git a/src/H5Tconv.c b/src/H5Tconv.c
index 8e40b0a..53f0b29 100644
--- a/src/H5Tconv.c
+++ b/src/H5Tconv.c
@@ -307,9 +307,23 @@ H5FL_BLK_DEFINE_STATIC(array_seq);
     H5T_CONV(H5T_CONV_xX, double, STYPE, DTYPE, ST, DT, D_MIN, D_MAX)	      \
 }
 
+/* Same as H5T_CONV_Xx_CORE, except that instead of using D_MAX and D_MIN
+ * when an overflow occurs, use the 'float' infinity values.
+ */
+#define H5T_CONV_Ff_CORE(S,D,STYPE,DTYPE,ST,DT,D_MIN,D_MAX) {		      \
+    if (*((ST*)S) > (DT)(D_MAX)) {					      \
+        if (!H5T_overflow_g || (H5T_overflow_g)(src_id, dst_id, S, D)<0)      \
+            *((DT*)D) = (H5T_NATIVE_FLOAT_POS_INF_g);			      \
+    } else if (*((ST*)S) < (D_MIN)) {					      \
+        if (!H5T_overflow_g || (H5T_overflow_g)(src_id, dst_id, S, D)<0)      \
+            *((DT*)D) = (H5T_NATIVE_FLOAT_NEG_INF_g);			      \
+    } else								      \
+        *((DT*)D) = (DT)(*((ST*)S));					      \
+}
+
 #define H5T_CONV_Ff(STYPE,DTYPE,ST,DT,D_MIN,D_MAX) {			      \
     assert(sizeof(ST)>=sizeof(DT));					      \
-    H5T_CONV(H5T_CONV_Xx, double, STYPE, DTYPE, ST, DT, D_MIN, D_MAX)	      \
+    H5T_CONV(H5T_CONV_Ff, double, STYPE, DTYPE, ST, DT, D_MIN, D_MAX)	      \
 }
 
 #define H5T_CONV_xF(STYPE,DTYPE,ST,DT,D_MIN,D_MAX) {			      \
@@ -3351,6 +3365,8 @@ H5T_conv_f_f (hid_t src_id, hid_t dst_id, H5T_cdata_t *cdata, hsize_t nelmts,
                 }
                     
                 /* Write the exponent */
+#ifdef OLD_WAY  
+/* It appears to be incorrect to increment the exponent when the carry is set -QAK */
                 if (carry) {
                     expo++;
                     if (expo>=expo_max) {
@@ -3379,6 +3395,7 @@ H5T_conv_f_f (hid_t src_id, hid_t dst_id, H5T_cdata_t *cdata, hsize_t nelmts,
                         H5T_bit_set(d, dst.u.f.mpos, dst.u.f.msize, FALSE);
                     }
                 }
+#endif /* OLD_WAY */
                 H5_CHECK_OVERFLOW(expo,hssize_t,hsize_t);
                 H5T_bit_set_d(d, dst.u.f.epos, dst.u.f.esize, (hsize_t)expo);
 
@@ -6674,7 +6691,7 @@ H5T_conv_double_float (hid_t src_id, hid_t dst_id, H5T_cdata_t *cdata,
 
     FUNC_ENTER_NOAPI(H5T_conv_float_double, FAIL);
 
-    H5T_CONV_Ff(DOUBLE, FLOAT, double, float, FLT_MIN, FLT_MAX);
+    H5T_CONV_Ff(DOUBLE, FLOAT, double, float, -FLT_MAX, FLT_MAX);
 
 done:
     FUNC_LEAVE_NOAPI(ret_value);
diff --git a/src/H5Tpkg.h b/src/H5Tpkg.h
index 9058d90..d87cc2f 100644
--- a/src/H5Tpkg.h
+++ b/src/H5Tpkg.h
@@ -305,6 +305,13 @@ H5_DLLVAR size_t	H5T_NATIVE_UINT_LEAST64_ALIGN_g;
 H5_DLLVAR size_t	H5T_NATIVE_INT_FAST64_ALIGN_g;
 H5_DLLVAR size_t	H5T_NATIVE_UINT_FAST64_ALIGN_g;
 
+/* Useful floating-point values for conversion routines */
+/* (+/- Inf for all floating-point types) */
+H5_DLLVAR float H5T_NATIVE_FLOAT_POS_INF_g;
+H5_DLLVAR float H5T_NATIVE_FLOAT_NEG_INF_g;
+H5_DLLVAR double H5T_NATIVE_DOUBLE_POS_INF_g;
+H5_DLLVAR double H5T_NATIVE_DOUBLE_NEG_INF_g;
+
 /* Common functions */
 H5_DLL herr_t H5T_init_interface(void);
 H5_DLL H5T_t *H5T_create(H5T_class_t type, size_t size);
diff --git a/src/H5config.h.in b/src/H5config.h.in
index 20ea728..f2348c4 100644
--- a/src/H5config.h.in
+++ b/src/H5config.h.in
@@ -69,6 +69,9 @@
 /* Define to 1 if you have the `gettimeofday' function. */
 #undef HAVE_GETTIMEOFDAY
 
+/* Define to 1 if you have the `get_fpc_csr' function. */
+#undef HAVE_GET_FPC_CSR
+
 /* Define to 1 if you have the <globus_common.h> header file. */
 #undef HAVE_GLOBUS_COMMON_H
 
@@ -272,6 +275,9 @@
 /* Define to 1 if you have the <sys/filio.h> header file. */
 #undef HAVE_SYS_FILIO_H
 
+/* Define to 1 if you have the <sys/fpu.h> header file. */
+#undef HAVE_SYS_FPU_H
+
 /* Define to 1 if you have the <sys/ioctl.h> header file. */
 #undef HAVE_SYS_IOCTL_H