Correct isnan on x86; improve isnan on other architectures

7 files changed, 63 insertions, 9 deletions

diff --git a/vec_avx_double4.h b/vec_avx_double4.h
index 9102e67..9ed8739 100644
--- a/vec_avx_double4.h
+++ b/vec_avx_double4.h
@@ -569,7 +569,14 @@ namespace vecmathlib {
     intvec_t ilogb() const { return MF::vml_ilogb(*this); }
     boolvec_t isfinite() const { return MF::vml_isfinite(*this); }
     boolvec_t isinf() const { return MF::vml_isinf(*this); }
-    boolvec_t isnan() const { return _mm256_cmp_pd(v, v, _CMP_UNORD_Q); }
+    boolvec_t isnan() const
+    {
+#ifdef VML_HAVE_NAN
+      return _mm256_cmp_pd(v, v, _CMP_UNORD_Q);
+#else
+      return BV(false);
+#endif
+    }
     boolvec_t isnormal() const { return MF::vml_isnormal(*this); }
     realvec ldexp(int_t n) const { return MF::vml_ldexp(*this, n); }
     realvec ldexp(intvec_t n) const { return MF::vml_ldexp(*this, n); }
diff --git a/vec_avx_float8.h b/vec_avx_float8.h
index 6431227..1e50425 100644
--- a/vec_avx_float8.h
+++ b/vec_avx_float8.h
@@ -561,7 +561,14 @@ namespace vecmathlib {
     intvec_t ilogb() const { return MF::vml_ilogb(*this); }
     boolvec_t isfinite() const { return MF::vml_isfinite(*this); }
     boolvec_t isinf() const { return MF::vml_isinf(*this); }
-    boolvec_t isnan() const { return _mm256_cmp_ps(v, v, _CMP_UNORD_Q); }
+    boolvec_t isnan() const
+    {
+#ifdef VML_HAVE_NAN
+      return _mm256_cmp_ps(v, v, _CMP_UNORD_Q);
+#else
+      return BV(false);
+#endif
+    }
     boolvec_t isnormal() const { return MF::vml_isnormal(*this); }
     realvec ldexp(int_t n) const { return MF::vml_ldexp(*this, n); }
     realvec ldexp(intvec_t n) const { return MF::vml_ldexp(*this, n); }
diff --git a/vec_qpx_double4.h b/vec_qpx_double4.h
index afcbcbb..46cd1d0 100644
--- a/vec_qpx_double4.h
+++ b/vec_qpx_double4.h
@@ -574,7 +574,14 @@ namespace vecmathlib {
     }
     boolvec_t isfinite() const { return MF::vml_isfinite(*this); }
     boolvec_t isinf() const { return MF::vml_isinf(*this); }
-    boolvec_t isnan() const { return vec_tstnan(v, v); }
+    boolvec_t isnan() const
+    {
+#ifdef VML_HAVE_NAN
+      return vec_tstnan(v, v);
+#else
+      return BV(false);
+#endif
+    }
     boolvec_t isnormal() const { return MF::vml_isnormal(*this); }
     realvec ldexp(int_t n) const { return ldexp(intvec_t(n)); }
     realvec ldexp(intvec_t n) const
diff --git a/vec_sse_double1.h b/vec_sse_double1.h
index a30cd2d..e4680d9 100644
--- a/vec_sse_double1.h
+++ b/vec_sse_double1.h
@@ -434,7 +434,14 @@ namespace vecmathlib {
     boolvec_t isinf() const { return std::isinf(v); }
     boolvec_t isnan() const
     {
-      return _mm_ucomineq_sd(from_double(v), from_double(v));
+      // This is wrong:
+      // return _mm_ucomineq_sd(from_double(v), from_double(v));
+      // This works:
+      // char r;
+      // __asm__("ucomisd %[v],%[v]; setp %[r]": [r]"=q"(r): [v]"x"(v));
+      // return boolvec_t::scalar_t(r);
+      // This works as well:
+      return std::isnan(v);
     }
     boolvec_t isnormal() const { return std::isnormal(v); }
     realvec ldexp(int_t n) const { return std::ldexp(v, n); }
diff --git a/vec_sse_double2.h b/vec_sse_double2.h
index ee689cd..47586cb 100644
--- a/vec_sse_double2.h
+++ b/vec_sse_double2.h
@@ -556,7 +556,14 @@ namespace vecmathlib {
     intvec_t ilogb() const { return MF::vml_ilogb(*this); }
     boolvec_t isfinite() const { return MF::vml_isfinite(*this); }
     boolvec_t isinf() const { return MF::vml_isinf(*this); }
-    boolvec_t isnan() const { return _mm_cmpunord_pd(v, v);; }
+    boolvec_t isnan() const
+    {
+#ifdef VML_HAVE_NAN
+      return _mm_cmpunord_pd(v, v);
+#else
+      return BV(false);
+#endif
+    }
     boolvec_t isnormal() const { return MF::vml_isnormal(*this); }
     realvec ldexp(int_t n) const { return MF::vml_ldexp(*this, n); }
     realvec ldexp(intvec_t n) const { return MF::vml_ldexp(*this, n); }
diff --git a/vec_sse_float1.h b/vec_sse_float1.h
index 12ab1ae..0e04555 100644
--- a/vec_sse_float1.h
+++ b/vec_sse_float1.h
@@ -426,15 +426,27 @@ namespace vecmathlib {
     intvec_t ilogb() const
     {
       int_t r = std::ilogb(v);
-      if (r == FP_ILOGB0) r = std::numeric_limits<int_t>::min();
-      else if (r == FP_ILOGBNAN) r = std::numeric_limits<int_t>::max();
+      typedef std::numeric_limits<int_t> NL;
+      if (FP_ILOGB0 == NL::min() and r == FP_ILOGB0) r = NL::min();
+      else if (FP_ILOGBNAN == NL::max() and r == FP_ILOGBNAN) r = NL::max();
       return r;
     }
     boolvec_t isfinite() const { return std::isfinite(v); }
     boolvec_t isinf() const { return std::isinf(v); }
     boolvec_t isnan() const
     {
-      return _mm_ucomineq_ss(from_float(v), from_float(v));
+#if defined VML_HAVE_NAN
+      // This is wrong:
+      // return _mm_ucomineq_ss(from_float(v), from_float(v));
+      // This works:
+      // char r;
+      // __asm__("ucomiss %[v],%[v]; setp %[r]": [r]"=q"(r): [v]"x"(v));
+      // return boolvec_t::scalar_t(r);
+      // This works as well:
+      return std::isnan(v);
+#else
+      return BV(false);
+#endif
     }
     boolvec_t isnormal() const { return std::isnormal(v); }
     realvec ldexp(int_t n) const { return std::ldexp(v, n); }
diff --git a/vec_sse_float4.h b/vec_sse_float4.h
index 2163e09..366c7bb 100644
--- a/vec_sse_float4.h
+++ b/vec_sse_float4.h
@@ -550,7 +550,14 @@ namespace vecmathlib {
     intvec_t ilogb() const { return MF::vml_ilogb(*this); }
     boolvec_t isfinite() const { return MF::vml_isfinite(*this); }
     boolvec_t isinf() const { return MF::vml_isinf(*this); }
-    boolvec_t isnan() const { return _mm_cmpunord_ps(v, v);; }
+    boolvec_t isnan() const
+    {
+#if defined VML_HAVE_NAN
+      return _mm_cmpunord_ps(v, v);
+#else
+      return BV(false);
+#endif
+    }
     boolvec_t isnormal() const { return MF::vml_isnormal(*this); }
     realvec ldexp(int_t n) const { return MF::vml_ldexp(*this, n); }
     realvec ldexp(intvec_t n) const { return MF::vml_ldexp(*this, n); }