bug #1521: avoid signalling NaN in hypot and make it std::complex<> friendly.

This commit is contained in:
Gael Guennebaud 2018-04-04 13:47:23 +02:00
parent 73729025a4
commit e116f6847e
2 changed files with 20 additions and 25 deletions

View File

@ -347,31 +347,7 @@ struct norm1_retval
* Implementation of hypot *
****************************************************************************/
template<typename Scalar>
struct hypot_impl
{
typedef typename NumTraits<Scalar>::Real RealScalar;
static inline RealScalar run(const Scalar& x, const Scalar& y)
{
EIGEN_USING_STD_MATH(abs);
EIGEN_USING_STD_MATH(sqrt);
RealScalar _x = abs(x);
RealScalar _y = abs(y);
Scalar p, qp;
if(_x>_y)
{
p = _x;
qp = _y / p;
}
else
{
p = _y;
qp = _x / p;
}
if(p==RealScalar(0)) return RealScalar(0);
return p * sqrt(RealScalar(1) + qp*qp);
}
};
template<typename Scalar> struct hypot_impl;
template<typename Scalar>
struct hypot_retval

View File

@ -66,6 +66,25 @@ T generic_fast_tanh_float(const T& a_x)
return pdiv(p, q);
}
template<typename Scalar>
struct hypot_impl
{
typedef typename NumTraits<Scalar>::Real RealScalar;
static inline RealScalar run(const Scalar& x, const Scalar& y)
{
EIGEN_USING_STD_MATH(abs);
EIGEN_USING_STD_MATH(sqrt);
RealScalar _x = abs(x);
RealScalar _y = abs(y);
RealScalar p, qp;
p = numext::maxi(_x,_y);
if(p==RealScalar(0)) return RealScalar(0);
qp = numext::mini(_y,_x) / p;
return p * sqrt(RealScalar(1) + qp*qp);
}
};
} // end namespace internal
} // end namespace Eigen