Bug #1777: make the scalar and packet path consistent for the logistic function + respective unit test

2 files changed, 19 insertions, 14 deletions

diff --git a/Eigen/src/Core/functors/UnaryFunctors.h b/Eigen/src/Core/functors/UnaryFunctors.h
index a55d7b74e..9a81fa365 100644
--- a/Eigen/src/Core/functors/UnaryFunctors.h
+++ b/Eigen/src/Core/functors/UnaryFunctors.h
@@ -951,8 +951,7 @@ template <typename T>
 struct scalar_logistic_op {
   EIGEN_EMPTY_STRUCT_CTOR(scalar_logistic_op)
   EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE T operator()(const T& x) const {
-    const T one = T(1);
-    return one / (one + numext::exp(-x));
+    return packetOp(x);
   }
 
   template <typename Packet> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
@@ -978,17 +977,7 @@ template <>
 struct scalar_logistic_op<float> {
   EIGEN_EMPTY_STRUCT_CTOR(scalar_logistic_op)
   EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE float operator()(const float& x) const {
-    // The upper cut-off is the smallest x for which the rational approximation evaluates to 1.
-    // Choosing this value saves us a few instructions clamping the results at the end.
-#ifdef EIGEN_VECTORIZE_FMA
-    const float cutoff_upper = 15.7243833541870117f;
-#else
-    const float cutoff_upper = 15.6437711715698242f;
-#endif
-    const float cutoff_lower = -9.f;
-    if (x > cutoff_upper) return 1.0f;
-    else if (x < cutoff_lower) return numext::exp(x);
-    else return 1.0f / (1.0f + numext::exp(-x));
+    return packetOp(x);
   }
 
   template <typename Packet> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
@@ -997,7 +986,8 @@ struct scalar_logistic_op<float> {
     const Packet lt_mask = pcmp_lt<Packet>(_x, cutoff_lower);
     const bool any_small = predux(lt_mask);
 
-    // Clamp the input to be at most 'cutoff_upper'.
+    // The upper cut-off is the smallest x for which the rational approximation evaluates to 1.
+    // Choosing this value saves us a few instructions clamping the results at the end.
 #ifdef EIGEN_VECTORIZE_FMA
     const Packet cutoff_upper = pset1<Packet>(15.7243833541870117f);
 #else
diff --git a/test/packetmath.cpp b/test/packetmath.cpp
index 7f2e69a50..032c0e4c8 100644
--- a/test/packetmath.cpp
+++ b/test/packetmath.cpp
@@ -447,6 +447,21 @@ template<typename Scalar,typename Packet> void packetmath_real()
     VERIFY((numext::isnan)(data2[0]));
   }
 
+  {
+    internal::scalar_logistic_op<Scalar> logistic;
+    for (int i=0; i<size; ++i)
+    {
+      data1[i] = internal::random<Scalar>(-20,20);
+    }
+    internal::pstore(data2, logistic.packetOp(internal::pload<Packet>(data1)));
+    for (int i=0; i<PacketSize; ++i) {
+      VERIFY_IS_APPROX(data2[i],logistic(data1[i]));
+      #ifdef EIGEN_VECTORIZE // don't check for exactness when using the i387 FPU
+      VERIFY_IS_EQUAL(data2[i],logistic(data1[i]));
+      #endif
+    }
+  }
+
 #if EIGEN_HAS_C99_MATH && (__cplusplus > 199711L)
   data1[0] = std::numeric_limits<Scalar>::infinity();
   data1[1] = Scalar(-1);