diff --git a/Eigen/src/Core/MathFunctions.h b/Eigen/src/Core/MathFunctions.h index 5df2d8bca..d415bc719 100644 --- a/Eigen/src/Core/MathFunctions.h +++ b/Eigen/src/Core/MathFunctions.h @@ -332,37 +332,45 @@ inline NewType cast(const OldType& x) * Implementation of atanh2 * ****************************************************************************/ -template -struct atanh2_default_impl +template +struct atanh2_impl { - typedef Scalar retval; - typedef typename NumTraits::Real RealScalar; - static inline Scalar run(const Scalar& x, const Scalar& y) + static inline Scalar run(const Scalar& x, const Scalar& r) + { + EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar) + #if __cplusplus >= 201103L + using std::log1p; + return log1p(2 * x / (r - x)) / 2; + #else + using std::abs; + using std::log; + using std::sqrt; + Scalar z = x / r; + if (r == 0 || abs(z) > sqrt(NumTraits::epsilon())) + return log((r + x) / (r - x)) / 2; + else + return z + z*z*z / 3; + #endif + } +}; + +template +struct atanh2_impl > +{ + typedef std::complex Scalar; + static inline Scalar run(const Scalar& x, const Scalar& r) { - using std::abs; using std::log; + using std::norm; using std::sqrt; - Scalar z = x / y; - if (y == Scalar(0) || abs(z) > sqrt(NumTraits::epsilon())) - return RealScalar(0.5) * log((y + x) / (y - x)); + Scalar z = x / r; + if (r == Scalar(0) || norm(z) > NumTraits::epsilon()) + return RealScalar(0.5) * log((r + x) / (r - x)); else return z + z*z*z / RealScalar(3); } }; -template -struct atanh2_default_impl -{ - static inline Scalar run(const Scalar&, const Scalar&) - { - EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar) - return Scalar(0); - } -}; - -template -struct atanh2_impl : atanh2_default_impl::IsInteger> {}; - template struct atanh2_retval {