mirror/eigen
2
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2025-02-17 18:09:55 +08:00
Code Issues Packages Projects Releases Wiki Activity

Bug #1785: Introduce numext::rint.

Browse Source 
This provides a new op that matches std::rint and previous behavior of
pround. Also adds corresponding unsupported/../Tensor op.
Performance is the same as e. g. floor (tested SSE/AVX).
This commit is contained in:
Ilya Tokar 2019-12-16 16:00:35 -05:00
parent d0ae052da4
commit 19876ced76
12 changed files with 117 additions and 1 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
Eigen/src/Core/GenericPacketMath.h
View File

@ -92,6 +92,7 @@ struct default_packet_traits
HasBetaInc = 0,
HasRound = 0,
HasRint = 0,
HasFloor = 0,
HasCeil = 0,
@ -575,6 +576,11 @@ Packet pround(const Packet& a) { using numext::round; return round(a); }
template<typename Packet> EIGEN_DECLARE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS
Packet pfloor(const Packet& a) { using numext::floor; return floor(a); }
/** \internal \returns the rounded value of \a a (coeff-wise) with current
* rounding mode */
template<typename Packet> EIGEN_DECLARE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS
Packet print(const Packet& a) { using numext::rint; return rint(a); }
/** \internal \returns the ceil of \a a (coeff-wise) */
template<typename Packet> EIGEN_DECLARE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS
Packet pceil(const Packet& a) { using numext::ceil; return ceil(a); }

1
Eigen/src/Core/GlobalFunctions.h
View File

@ -89,6 +89,7 @@ namespace Eigen
EIGEN_ARRAY_DECLARE_GLOBAL_UNARY(rsqrt,scalar_rsqrt_op,reciprocal square root,\sa ArrayBase::rsqrt)
EIGEN_ARRAY_DECLARE_GLOBAL_UNARY(square,scalar_square_op,square (power 2),\sa Eigen::abs2 DOXCOMMA Eigen::pow DOXCOMMA ArrayBase::square)
EIGEN_ARRAY_DECLARE_GLOBAL_UNARY(cube,scalar_cube_op,cube (power 3),\sa Eigen::pow DOXCOMMA ArrayBase::cube)
EIGEN_ARRAY_DECLARE_GLOBAL_UNARY(rint,scalar_rint_op,nearest integer,\sa Eigen::floor DOXCOMMA Eigen::ceil DOXCOMMA ArrayBase::round)
EIGEN_ARRAY_DECLARE_GLOBAL_UNARY(round,scalar_round_op,nearest integer,\sa Eigen::floor DOXCOMMA Eigen::ceil DOXCOMMA ArrayBase::round)
EIGEN_ARRAY_DECLARE_GLOBAL_UNARY(floor,scalar_floor_op,nearest integer not greater than the giben value,\sa Eigen::ceil DOXCOMMA ArrayBase::floor)
EIGEN_ARRAY_DECLARE_GLOBAL_UNARY(ceil,scalar_ceil_op,nearest integer not less than the giben value,\sa Eigen::floor DOXCOMMA ArrayBase::ceil)

39
Eigen/src/Core/MathFunctions.h
View File

@ -429,6 +429,38 @@ struct round_retval
typedef Scalar type;
};
/****************************************************************************
* Implementation of rint *
****************************************************************************/
template<typename Scalar>
struct rint_impl {
static inline Scalar run(const Scalar& x)
{
EIGEN_STATIC_ASSERT((!NumTraits<Scalar>::IsComplex), NUMERIC_TYPE_MUST_BE_REAL)
#if EIGEN_HAS_CXX11_MATH
EIGEN_USING_STD_MATH(rint);
#endif
return rint(x);
}
};
#if !EIGEN_HAS_CXX11_MATH
template<>
struct rint_impl<float> {
static inline float run(const float& x)
{
return rintf(x);
}
};
#endif
template<typename Scalar>
struct rint_retval
{
typedef Scalar type;
};
/****************************************************************************
* Implementation of arg *
****************************************************************************/
@ -1194,6 +1226,13 @@ SYCL_SPECIALIZE_FLOATING_TYPES_UNARY_FUNC_RET_TYPE(isinf, isinf, bool)
SYCL_SPECIALIZE_FLOATING_TYPES_UNARY_FUNC_RET_TYPE(isfinite, isfinite, bool)
#endif
template<typename Scalar>
EIGEN_DEVICE_FUNC
inline EIGEN_MATHFUNC_RETVAL(rint, Scalar) rint(const Scalar& x)
{
return EIGEN_MATHFUNC_IMPL(rint, Scalar)::run(x);
}
template<typename Scalar>
EIGEN_DEVICE_FUNC
inline EIGEN_MATHFUNC_RETVAL(round, Scalar) round(const Scalar& x)

5
Eigen/src/Core/arch/AVX/PacketMath.h
View File

@ -81,7 +81,8 @@ template<> struct packet_traits<float> : default_packet_traits
HasBlend = 1,
HasRound = 1,
HasFloor = 1,
HasCeil = 1
HasCeil = 1,
HasRint = 1
};
};
template<> struct packet_traits<double> : default_packet_traits
@ -316,6 +317,8 @@ template<> EIGEN_STRONG_INLINE Packet8i pcmp_eq(const Packet8i& a, const Packet8
#endif
}
template<> EIGEN_STRONG_INLINE Packet8f print<Packet8f>(const Packet8f& a) { return _mm256_round_ps(a, _MM_FROUND_CUR_DIRECTION); }
template<> EIGEN_STRONG_INLINE Packet4d print<Packet4d>(const Packet4d& a) { return _mm256_round_pd(a, _MM_FROUND_CUR_DIRECTION); }
template<> EIGEN_STRONG_INLINE Packet8f pceil<Packet8f>(const Packet8f& a) { return _mm256_ceil_ps(a); }
template<> EIGEN_STRONG_INLINE Packet4d pceil<Packet4d>(const Packet4d& a) { return _mm256_ceil_pd(a); }

5
Eigen/src/Core/arch/SSE/PacketMath.h
View File

@ -124,6 +124,7 @@ struct packet_traits<float> : default_packet_traits {
#ifdef EIGEN_VECTORIZE_SSE4_1
,
HasRint = 1,
HasRound = 1,
HasCeil = 1
#endif
@ -148,6 +149,7 @@ struct packet_traits<double> : default_packet_traits {
#ifdef EIGEN_VECTORIZE_SSE4_1
,
HasRound = 1,
HasRint = 1,
HasFloor = 1,
HasCeil = 1
#endif
@ -443,6 +445,9 @@ template<> EIGEN_STRONG_INLINE Packet2d pround<Packet2d>(const Packet2d& a)
return _mm_round_pd(padd(por(pand(a, mask), prev0dot5), a), _MM_FROUND_TO_ZERO);
}
template<> EIGEN_STRONG_INLINE Packet4f print<Packet4f>(const Packet4f& a) { return _mm_round_ps(a, _MM_FROUND_CUR_DIRECTION); }
template<> EIGEN_STRONG_INLINE Packet2d print<Packet2d>(const Packet2d& a) { return _mm_round_pd(a, _MM_FROUND_CUR_DIRECTION); }
template<> EIGEN_STRONG_INLINE Packet4f pceil<Packet4f>(const Packet4f& a) { return _mm_ceil_ps(a); }
template<> EIGEN_STRONG_INLINE Packet2d pceil<Packet2d>(const Packet2d& a) { return _mm_ceil_pd(a); }

19
Eigen/src/Core/functors/UnaryFunctors.h
View File

@ -735,6 +735,25 @@ struct functor_traits<scalar_floor_op<Scalar> >
};
};
/** \internal
* \brief Template functor to compute the rounded (with current rounding mode) value of a scalar
* \sa class CwiseUnaryOp, ArrayBase::rint()
*/
template<typename Scalar> struct scalar_rint_op {
EIGEN_EMPTY_STRUCT_CTOR(scalar_rint_op)
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar operator() (const Scalar& a) const { return numext::rint(a); }
template <typename Packet>
EIGEN_DEVICE_FUNC inline Packet packetOp(const Packet& a) const { return internal::print(a); }
};
template<typename Scalar>
struct functor_traits<scalar_rint_op<Scalar> >
{
enum {
Cost = NumTraits<Scalar>::MulCost,
PacketAccess = packet_traits<Scalar>::HasRint
};
};
/** \internal
* \brief Template functor to compute the ceil of a scalar
* \sa class CwiseUnaryOp, ArrayBase::ceil()

15
Eigen/src/plugins/ArrayCwiseUnaryOps.h
View File

@ -32,6 +32,7 @@ typedef CwiseUnaryOp<internal::scalar_cosh_op<Scalar>, const Derived> CoshReturn
typedef CwiseUnaryOp<internal::scalar_square_op<Scalar>, const Derived> SquareReturnType;
typedef CwiseUnaryOp<internal::scalar_cube_op<Scalar>, const Derived> CubeReturnType;
typedef CwiseUnaryOp<internal::scalar_round_op<Scalar>, const Derived> RoundReturnType;
typedef CwiseUnaryOp<internal::scalar_rint_op<Scalar>, const Derived> RintReturnType;
typedef CwiseUnaryOp<internal::scalar_floor_op<Scalar>, const Derived> FloorReturnType;
typedef CwiseUnaryOp<internal::scalar_ceil_op<Scalar>, const Derived> CeilReturnType;
typedef CwiseUnaryOp<internal::scalar_isnan_op<Scalar>, const Derived> IsNaNReturnType;
@ -427,6 +428,20 @@ cube() const
return CubeReturnType(derived());
}
/** \returns an expression of the coefficient-wise rint of *this.
*
* Example: \include Cwise_rint.cpp
* Output: \verbinclude Cwise_rint.out
*
* \sa <a href="group__CoeffwiseMathFunctions.html#cwisetable_rint">Math functions</a>, ceil(), floor()
*/
EIGEN_DEVICE_FUNC
inline const RintReturnType
rint() const
{
return RintReturnType(derived());
}
/** \returns an expression of the coefficient-wise round of *this.
*
* Example: \include Cwise_round.cpp

11
doc/CoeffwiseMathFunctionsTable.dox
View File

@ -394,6 +394,17 @@ This also means that, unless specified, if the function \c std::foo is available
plus \c using <a href="http://en.cppreference.com/w/cpp/numeric/math/round">\c std::round </a>; \cpp11</td>
<td>SSE4,AVX,ZVector (f,d)</td>
</tr>
<tr>
<td class="code">
\anchor cwisetable_rint
a.\link ArrayBase::rint rint\endlink(); \n
\link Eigen::rint rint\endlink(a);
</td>
<td>nearest integer, \n rounding to nearest even in halfway cases</td>
<td>built-in generic implementation using <a href="http://en.cppreference.com/w/cpp/numeric/math/rint">\c std::rint</a>
or <a href="http://en.cppreference.com/w/c/numeric/math/rint">\c rintf</a>; </td>
<td>SSE4,AVX (f,d)</td>
</tr>
<tr>
<th colspan="4">Floating point manipulation functions</th>
</tr>

3
doc/snippets/Cwise_rint.cpp Normal file
View File

@ -0,0 +1,3 @@
ArrayXd v = ArrayXd::LinSpaced(7,-2,2);
cout << v << endl << endl;
cout << rint(v) << endl;

6
test/array_cwise.cpp
View File

@ -296,6 +296,7 @@ template<typename ArrayType> void array_real(const ArrayType& m)
VERIFY_IS_APPROX(m1.arg(), arg(m1));
VERIFY_IS_APPROX(m1.round(), round(m1));
VERIFY_IS_APPROX(m1.rint(), rint(m1));
VERIFY_IS_APPROX(m1.floor(), floor(m1));
VERIFY_IS_APPROX(m1.ceil(), ceil(m1));
VERIFY((m1.isNaN() == (Eigen::isnan)(m1)).all());
@ -331,6 +332,11 @@ template<typename ArrayType> void array_real(const ArrayType& m)
VERIFY_IS_APPROX(logistic(m1), (1.0/(1.0+exp(-m1))));
VERIFY_IS_APPROX(arg(m1), ((m1<0).template cast<Scalar>())*std::acos(-1.0));
VERIFY((round(m1) <= ceil(m1) && round(m1) >= floor(m1)).all());
VERIFY((rint(m1) <= ceil(m1) && rint(m1) >= floor(m1)).all());
VERIFY(((ceil(m1) - round(m1)) <= Scalar(0.5) || (round(m1) - floor(m1)) <= Scalar(0.5)).all());
VERIFY(((ceil(m1) - round(m1)) <= Scalar(1.0) && (round(m1) - floor(m1)) <= Scalar(1.0)).all());
VERIFY(((ceil(m1) - rint(m1)) <= Scalar(0.5) || (rint(m1) - floor(m1)) <= Scalar(0.5)).all());
VERIFY(((ceil(m1) - rint(m1)) <= Scalar(1.0) && (rint(m1) - floor(m1)) <= Scalar(1.0)).all());
VERIFY((Eigen::isnan)((m1*0.0)/0.0).all());
VERIFY((Eigen::isinf)(m4/0.0).all());
VERIFY(((Eigen::isfinite)(m1) && (!(Eigen::isfinite)(m1*0.0/0.0)) && (!(Eigen::isfinite)(m4/0.0))).all());

2
test/packetmath.cpp
View File

@ -518,6 +518,7 @@ template<typename Scalar,typename Packet> void packetmath_real()
CHECK_CWISE1_IF(PacketTraits::HasRound, numext::round, internal::pround);
CHECK_CWISE1_IF(PacketTraits::HasCeil, numext::ceil, internal::pceil);
CHECK_CWISE1_IF(PacketTraits::HasFloor, numext::floor, internal::pfloor);
CHECK_CWISE1_IF(PacketTraits::HasRint, numext::rint, internal::print);
// See bug 1785.
for (int i=0; i<size; ++i)
@ -526,6 +527,7 @@ template<typename Scalar,typename Packet> void packetmath_real()
data2[i] = -1.5 + i;
}
CHECK_CWISE1_IF(PacketTraits::HasRound, numext::round, internal::pround);
CHECK_CWISE1_IF(PacketTraits::HasRint, numext::rint, internal::print);
for (int i=0; i<size; ++i)
{

6
unsupported/Eigen/CXX11/src/Tensor/TensorBase.h
View File

@ -422,6 +422,12 @@ class TensorBase<Derived, ReadOnlyAccessors>
return unaryExpr(internal::scalar_round_op<Scalar>());
}
EIGEN_DEVICE_FUNC
EIGEN_STRONG_INLINE const TensorCwiseUnaryOp<internal::scalar_rint_op<Scalar>, const Derived>
rint() const {
return unaryExpr(internal::scalar_rint_op<Scalar>());
}
EIGEN_DEVICE_FUNC
EIGEN_STRONG_INLINE const TensorCwiseUnaryOp<internal::scalar_ceil_op<Scalar>, const Derived>
ceil() const {