Correct some frexp/ilogb implementations

3 files changed, 43 insertions, 11 deletions

diff --git a/vec_pseudo.h b/vec_pseudo.h
index 98948b2..1e78fb7 100644
--- a/vec_pseudo.h
+++ b/vec_pseudo.h
@@ -752,9 +752,12 @@ namespace vecmathlib {
       for (int d=0; d<size; ++d) {
         int ir;
         real_t r = std::frexp(v[d], &ir);
-        typedef std::numeric_limits<int_t> NL;
-        if (FP_ILOGB0 != NL::min() and ir == FP_ILOGB0) ir = NL::min();
-        else if (FP_ILOGBNAN != NL::max() and ir == FP_ILOGBNAN) ir = NL::max();
+#if defined VML_HAVE_INF
+        if (std::isinf(v[d])) ir = std::numeric_limits<int_t>::max();
+#endif
+#if defined VML_HAVE_NAN
+        if (std::isnan(v[d])) ir = std::numeric_limits<int_t>::max();
+#endif
         res.v[d] = r;
         ires.v[d] = ir;
       }
@@ -767,8 +770,13 @@ namespace vecmathlib {
       for (int d=0; d<size; ++d) {
         int_t r = std::ilogb(v[d]);
         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();
+        if (FP_ILOGB0 != NL::min() and v[d] == R(0.0)) {
+          r = NL::min();
+#if defined VML_HAVE_NAN
+        } else if (FP_ILOGBNAN != NL::max() and std::isnan(v[d])) {
+          r = NL::max();
+#endif
+        }
         res.v[d] = r;
       }
       return res;
diff --git a/vec_sse_double1.h b/vec_sse_double1.h
index b05168f..76487af 100644
--- a/vec_sse_double1.h
+++ b/vec_sse_double1.h
@@ -419,14 +419,27 @@ namespace vecmathlib {
       return to_double(_mm_min_sd(from_double(v), from_double(y.v)));
     }
     realvec fmod(realvec y) const { return std::fmod(v, y.v); }
-    realvec frexp(intvec_t& r) const { return MF::vml_frexp(*this, r); }
+    realvec frexp(intvec_t& ir) const
+    {
+      int iri;
+      realvec r = std::frexp(v, &iri);
+      ir.v = iri;
+      if (isinf()) ir.v = std::numeric_limits<int_t>::max();
+      if (isnan()) ir.v = std::numeric_limits<int_t>::max();
+      return r;
+    }
     realvec hypot(realvec y) const { return MF::vml_hypot(*this, y); }
     intvec_t ilogb() const
     {
       int_t r = std::ilogb(v);
       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();
+      if (FP_ILOGB0 != NL::min() and v == R(0.0)) {
+        r = NL::min();
+#if defined VML_HAVE_NAN
+      } else if (FP_ILOGBNAN != NL::max() and std::isnan(v)) {
+        r = NL::max();
+#endif
+      }
       return r;
     }
     boolvec_t isfinite() const { return std::isfinite(v); }
diff --git a/vec_sse_float1.h b/vec_sse_float1.h
index 6aa42dc..d3039ad 100644
--- a/vec_sse_float1.h
+++ b/vec_sse_float1.h
@@ -419,14 +419,25 @@ namespace vecmathlib {
       return to_float(_mm_min_ss(from_float(v), from_float(y.v)));
     }
     realvec fmod(realvec y) const { return std::fmod(v, y.v); }
-    realvec frexp(intvec_t& r) const { return MF::vml_frexp(*this, r); }
+    realvec frexp(intvec_t& ir) const
+    {
+      realvec r = std::frexp(v, &ir.v);
+      if (isinf()) ir.v = std::numeric_limits<int_t>::max();
+      if (isnan()) ir.v = std::numeric_limits<int_t>::max();
+      return r;
+    }
     realvec hypot(realvec y) const { return MF::vml_hypot(*this, y); }
     intvec_t ilogb() const
     {
       int_t r = std::ilogb(v);
       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();
+      if (FP_ILOGB0 != NL::min() and *this == RV(R(0.0))) {
+        r = NL::min();
+#if defined VML_HAVE_NAN
+      } else if (FP_ILOGBNAN != NL::max() and isnan()) {
+        r = NL::max();
+#endif
+      }
       return r;
     }
     boolvec_t isfinite() const { return std::isfinite(v); }