From 8dfe1029a5012774d50bfc39569e49e34d3152e0 Mon Sep 17 00:00:00 2001 From: Antonio Sanchez Date: Tue, 16 Mar 2021 20:12:46 -0700 Subject: Augment NumTraits with min/max_exponent() again. Replace usage of `std::numeric_limits<...>::min/max_exponent` in codebase where possible. Also replaced some other `numeric_limits` usages in affected tests with the `NumTraits` equivalent. The previous MR !443 failed for c++03 due to lack of `constexpr`. Because of this, we need to keep around the `std::numeric_limits` version in enum expressions until the switch to c++11. Fixes #2148 --- Eigen/src/Core/NumTraits.h | 22 +++++++- Eigen/src/Core/StableNorm.h | 6 +- .../Core/arch/Default/GenericPacketMathFunctions.h | 6 +- bench/bench_norm.cpp | 12 ++-- test/array_cwise.cpp | 10 ++-- test/packetmath.cpp | 64 +++++++++++----------- 6 files changed, 70 insertions(+), 50 deletions(-) diff --git a/Eigen/src/Core/NumTraits.h b/Eigen/src/Core/NumTraits.h index 609e11402..fdd4d4f51 100644 --- a/Eigen/src/Core/NumTraits.h +++ b/Eigen/src/Core/NumTraits.h @@ -135,9 +135,18 @@ EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC Tgt bit_cast(const Src& src) { * \li A dummy_precision() function returning a weak epsilon value. It is mainly used as a default * value by the fuzzy comparison operators. * \li highest() and lowest() functions returning the highest and lowest possible values respectively. + * \li digits() function returning the number of radix digits (non-sign digits for integers, mantissa for floating-point). This is + * the analogue of std::numeric_limits::digits + * which is used as the default implementation if specialized. * \li digits10() function returning the number of decimal digits that can be represented without change. This is * the analogue of std::numeric_limits::digits10 * which is used as the default implementation if specialized. + * \li min_exponent() and max_exponent() functions returning the highest and lowest possible values, respectively, + * such that the radix raised to the power exponent-1 is a normalized floating-point number. These are equivalent to + * std::numeric_limits::min_exponent/ + * std::numeric_limits::max_exponent. + * \li infinity() function returning a representation of positive infinity, if available. + * \li quiet_NaN function returning a non-signaling "not-a-number", if available. */ template struct GenericNumTraits @@ -179,6 +188,18 @@ template struct GenericNumTraits return internal::default_digits_impl::run(); } + EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR + static inline int min_exponent() + { + return numext::numeric_limits::min_exponent; + } + + EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR + static inline int max_exponent() + { + return numext::numeric_limits::max_exponent; + } + EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR static inline Real dummy_precision() { @@ -186,7 +207,6 @@ template struct GenericNumTraits return Real(0); } - EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR static inline T highest() { return (numext::numeric_limits::max)(); diff --git a/Eigen/src/Core/StableNorm.h b/Eigen/src/Core/StableNorm.h index 04fe0b371..4a3f0cca8 100644 --- a/Eigen/src/Core/StableNorm.h +++ b/Eigen/src/Core/StableNorm.h @@ -134,9 +134,9 @@ blueNorm_impl(const EigenBase& _vec) // statements can be replaced static const int ibeta = std::numeric_limits::radix; // base for floating-point numbers static const int it = NumTraits::digits(); // number of base-beta digits in mantissa - static const int iemin = std::numeric_limits::min_exponent; // minimum exponent - static const int iemax = std::numeric_limits::max_exponent; // maximum exponent - static const RealScalar rbig = (std::numeric_limits::max)(); // largest floating-point number + static const int iemin = NumTraits::min_exponent(); // minimum exponent + static const int iemax = NumTraits::max_exponent(); // maximum exponent + static const RealScalar rbig = NumTraits::highest(); // largest floating-point number static const RealScalar b1 = RealScalar(pow(RealScalar(ibeta),RealScalar(-((1-iemin)/2)))); // lower boundary of midrange static const RealScalar b2 = RealScalar(pow(RealScalar(ibeta),RealScalar((iemax + 1 - it)/2))); // upper boundary of midrange static const RealScalar s1m = RealScalar(pow(RealScalar(ibeta),RealScalar((2-iemin)/2))); // scaling factor for lower range diff --git a/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h b/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h index 411640ee8..87e8c2703 100644 --- a/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h +++ b/Eigen/src/Core/arch/Default/GenericPacketMathFunctions.h @@ -135,7 +135,7 @@ Packet pldexp_generic(const Packet& a, const Packet& exponent) { // Explicitly multiplies // a * (2^e) // clamping e to the range -// [numeric_limits::min_exponent-2, numeric_limits::max_exponent] +// [NumTraits::min_exponent()-2, NumTraits::max_exponent()] // // This is approx 7x faster than pldexp_impl, but will prematurely over/underflow // if 2^e doesn't fit into a normal floating-point Scalar. @@ -1480,8 +1480,8 @@ Packet generic_pow(const Packet& x, const Packet& y) { const Packet y_is_nan = pandnot(ptrue(y), pcmp_eq(y, y)); const Packet abs_y_is_inf = pcmp_eq(pabs(y), cst_pos_inf); EIGEN_CONSTEXPR Scalar huge_exponent = - (std::numeric_limits::max_exponent * Scalar(EIGEN_LN2)) / - std::numeric_limits::epsilon(); + (NumTraits::max_exponent() * Scalar(EIGEN_LN2)) / + NumTraits::epsilon(); const Packet abs_y_is_huge = pcmp_le(pset1(huge_exponent), pabs(y)); // Predicates for whether y is integer and/or even. diff --git a/bench/bench_norm.cpp b/bench/bench_norm.cpp index a86153945..592f25d66 100644 --- a/bench/bench_norm.cpp +++ b/bench/bench_norm.cpp @@ -111,12 +111,12 @@ EIGEN_DONT_INLINE typename T::Scalar pblueNorm(const T& v) int nbig, ibeta, it, iemin, iemax, iexp; Scalar abig, eps; - nbig = std::numeric_limits::max(); // largest integer - ibeta = std::numeric_limits::radix; //NumTraits::Base; // base for floating-point numbers - it = std::numeric_limits::digits; //NumTraits::Mantissa; // number of base-beta digits in mantissa - iemin = std::numeric_limits::min_exponent; // minimum exponent - iemax = std::numeric_limits::max_exponent; // maximum exponent - rbig = std::numeric_limits::max(); // largest floating-point number + nbig = NumTraits::highest(); // largest integer + ibeta = std::numeric_limits::radix; // NumTraits::Base; // base for floating-point numbers + it = NumTraits::digits(); // NumTraits::Mantissa; // number of base-beta digits in mantissa + iemin = NumTraits::min_exponent(); // minimum exponent + iemax = NumTraits::max_exponent(); // maximum exponent + rbig = NumTraits::highest(); // largest floating-point number // Check the basic machine-dependent constants. if(iemin > 1 - 2*it || 1+it>iemax || (it==2 && ibeta<5) diff --git a/test/array_cwise.cpp b/test/array_cwise.cpp index 6a88e0e27..1bc8e19f9 100644 --- a/test/array_cwise.cpp +++ b/test/array_cwise.cpp @@ -14,18 +14,18 @@ template void pow_test() { const Scalar zero = Scalar(0); - const Scalar eps = std::numeric_limits::epsilon(); + const Scalar eps = Eigen::NumTraits::epsilon(); const Scalar one = Scalar(1); const Scalar two = Scalar(2); const Scalar three = Scalar(3); const Scalar sqrt_half = Scalar(std::sqrt(0.5)); const Scalar sqrt2 = Scalar(std::sqrt(2)); - const Scalar inf = std::numeric_limits::infinity(); - const Scalar nan = std::numeric_limits::quiet_NaN(); + const Scalar inf = Eigen::NumTraits::infinity(); + const Scalar nan = Eigen::NumTraits::quiet_NaN(); const Scalar denorm_min = std::numeric_limits::denorm_min(); const Scalar min = (std::numeric_limits::min)(); const Scalar max = (std::numeric_limits::max)(); - const Scalar max_exp = (static_cast(int(std::numeric_limits::max_exponent)) * Scalar(EIGEN_LN2)) / eps; + const Scalar max_exp = (static_cast(int(Eigen::NumTraits::max_exponent())) * Scalar(EIGEN_LN2)) / eps; const static Scalar abs_vals[] = {zero, denorm_min, @@ -613,7 +613,7 @@ template void min_max(const ArrayType& m) // min/max with various NaN propagation options. if (m1.size() > 1 && !NumTraits::IsInteger) { - m1(0,0) = std::numeric_limits::quiet_NaN(); + m1(0,0) = NumTraits::quiet_NaN(); maxM1 = m1.template maxCoeff(); minM1 = m1.template minCoeff(); VERIFY((numext::isnan)(maxM1)); diff --git a/test/packetmath.cpp b/test/packetmath.cpp index f5cce6436..67d329a67 100644 --- a/test/packetmath.cpp +++ b/test/packetmath.cpp @@ -273,7 +273,7 @@ void packetmath_boolean_mask_ops() { //Test NaN for (int i = 0; i < PacketSize; ++i) { - data1[i] = std::numeric_limits::quiet_NaN(); + data1[i] = NumTraits::quiet_NaN(); data1[i + PacketSize] = internal::random() ? data1[i] : Scalar(0); } CHECK_CWISE2_IF(true, internal::pcmp_eq, internal::pcmp_eq); @@ -634,7 +634,7 @@ void packetmath_real() { if (PacketTraits::HasExp) { // Check denormals: for (int j=0; j<3; ++j) { - data1[0] = Scalar(std::ldexp(1, std::numeric_limits::min_exponent-j)); + data1[0] = Scalar(std::ldexp(1, NumTraits::min_exponent()-j)); CHECK_CWISE1_BYREF1_IF(PacketTraits::HasExp, REF_FREXP, internal::pfrexp); data1[0] = -data1[0]; CHECK_CWISE1_BYREF1_IF(PacketTraits::HasExp, REF_FREXP, internal::pfrexp); @@ -671,10 +671,10 @@ void packetmath_real() { if (PacketTraits::HasExp) { data1[0] = Scalar(-1); // underflow to zero - data1[PacketSize] = Scalar(std::numeric_limits::min_exponent-55); + data1[PacketSize] = Scalar(NumTraits::min_exponent()-55); CHECK_CWISE2_IF(PacketTraits::HasExp, REF_LDEXP, internal::pldexp); // overflow to inf - data1[PacketSize] = Scalar(std::numeric_limits::max_exponent+10); + data1[PacketSize] = Scalar(NumTraits::max_exponent()+10); CHECK_CWISE2_IF(PacketTraits::HasExp, REF_LDEXP, internal::pldexp); // NaN stays NaN data1[0] = NumTraits::quiet_NaN(); @@ -682,21 +682,21 @@ void packetmath_real() { VERIFY((numext::isnan)(data2[0])); // inf stays inf data1[0] = NumTraits::infinity(); - data1[PacketSize] = Scalar(std::numeric_limits::min_exponent-10); + data1[PacketSize] = Scalar(NumTraits::min_exponent()-10); CHECK_CWISE2_IF(PacketTraits::HasExp, REF_LDEXP, internal::pldexp); // zero stays zero data1[0] = Scalar(0); - data1[PacketSize] = Scalar(std::numeric_limits::max_exponent+10); + data1[PacketSize] = Scalar(NumTraits::max_exponent()+10); CHECK_CWISE2_IF(PacketTraits::HasExp, REF_LDEXP, internal::pldexp); // Small number big exponent. - data1[0] = Scalar(std::ldexp(Scalar(1.0), std::numeric_limits::min_exponent-1)); - data1[PacketSize] = Scalar(-std::numeric_limits::min_exponent - +std::numeric_limits::max_exponent); + data1[0] = Scalar(std::ldexp(Scalar(1.0), NumTraits::min_exponent()-1)); + data1[PacketSize] = Scalar(-NumTraits::min_exponent() + +NumTraits::max_exponent()); CHECK_CWISE2_IF(PacketTraits::HasExp, REF_LDEXP, internal::pldexp); // Big number small exponent. - data1[0] = Scalar(std::ldexp(Scalar(1.0), std::numeric_limits::max_exponent-1)); - data1[PacketSize] = Scalar(+std::numeric_limits::min_exponent - -std::numeric_limits::max_exponent); + data1[0] = Scalar(std::ldexp(Scalar(1.0), NumTraits::max_exponent()-1)); + data1[PacketSize] = Scalar(+NumTraits::min_exponent() + -NumTraits::max_exponent()); CHECK_CWISE2_IF(PacketTraits::HasExp, REF_LDEXP, internal::pldexp); } @@ -707,8 +707,8 @@ void packetmath_real() { data1[0] = Scalar(1e-20); CHECK_CWISE1_IF(PacketTraits::HasTanh, std::tanh, internal::ptanh); if (PacketTraits::HasExp && PacketSize >= 2) { - const Scalar small = std::numeric_limits::epsilon(); - data1[0] = std::numeric_limits::quiet_NaN(); + const Scalar small = NumTraits::epsilon(); + data1[0] = NumTraits::quiet_NaN(); data1[1] = small; test::packet_helper h; h.store(data2, internal::pexp(h.load(data1))); @@ -742,7 +742,7 @@ void packetmath_real() { if (PacketTraits::HasTanh) { // NOTE this test migh fail with GCC prior to 6.3, see MathFunctionsImpl.h for details. - data1[0] = std::numeric_limits::quiet_NaN(); + data1[0] = NumTraits::quiet_NaN(); test::packet_helper::HasTanh, Packet> h; h.store(data2, internal::ptanh(h.load(data1))); VERIFY((numext::isnan)(data2[0])); @@ -762,17 +762,17 @@ void packetmath_real() { } #if EIGEN_HAS_C99_MATH && (EIGEN_COMP_CXXVER >= 11) - data1[0] = std::numeric_limits::infinity(); + data1[0] = NumTraits::infinity(); data1[1] = Scalar(-1); CHECK_CWISE1_IF(PacketTraits::HasLog1p, std::log1p, internal::plog1p); - data1[0] = std::numeric_limits::infinity(); - data1[1] = -std::numeric_limits::infinity(); + data1[0] = NumTraits::infinity(); + data1[1] = -NumTraits::infinity(); CHECK_CWISE1_IF(PacketTraits::HasExpm1, std::expm1, internal::pexpm1); #endif if (PacketSize >= 2) { - data1[0] = std::numeric_limits::quiet_NaN(); - data1[1] = std::numeric_limits::epsilon(); + data1[0] = NumTraits::quiet_NaN(); + data1[1] = NumTraits::epsilon(); if (PacketTraits::HasLog) { test::packet_helper h; h.store(data2, internal::plog(h.load(data1))); @@ -782,7 +782,7 @@ void packetmath_real() { VERIFY_IS_APPROX(std::log(data1[1]), data2[1]); } - data1[0] = -std::numeric_limits::epsilon(); + data1[0] = -NumTraits::epsilon(); data1[1] = Scalar(0); h.store(data2, internal::plog(h.load(data1))); VERIFY((numext::isnan)(data2[0])); @@ -813,14 +813,14 @@ void packetmath_real() { h.store(data2, internal::plog(h.load(data1))); VERIFY((numext::isnan)(data2[0])); - data1[0] = std::numeric_limits::infinity(); + data1[0] = NumTraits::infinity(); h.store(data2, internal::plog(h.load(data1))); VERIFY((numext::isinf)(data2[0])); } if (PacketTraits::HasLog1p) { test::packet_helper h; data1[0] = Scalar(-2); - data1[1] = -std::numeric_limits::infinity(); + data1[1] = -NumTraits::infinity(); h.store(data2, internal::plog1p(h.load(data1))); VERIFY((numext::isnan)(data2[0])); VERIFY((numext::isnan)(data2[1])); @@ -831,7 +831,7 @@ void packetmath_real() { if (std::numeric_limits::has_denorm == std::denorm_present) { data1[1] = -std::numeric_limits::denorm_min(); } else { - data1[1] = -std::numeric_limits::epsilon(); + data1[1] = -NumTraits::epsilon(); } h.store(data2, internal::psqrt(h.load(data1))); VERIFY((numext::isnan)(data2[0])); @@ -842,7 +842,7 @@ void packetmath_real() { && !internal::is_same::value && !internal::is_same::value) { test::packet_helper h; - for (Scalar k = Scalar(1); k < Scalar(10000) / std::numeric_limits::epsilon(); k *= Scalar(2)) { + for (Scalar k = Scalar(1); k < Scalar(10000) / NumTraits::epsilon(); k *= Scalar(2)) { for (int k1 = 0; k1 <= 1; ++k1) { data1[0] = Scalar((2 * double(k) + k1) * double(EIGEN_PI) / 2 * internal::random(0.8, 1.2)); data1[1] = Scalar((2 * double(k) + 2 + k1) * double(EIGEN_PI) / 2 * internal::random(0.8, 1.2)); @@ -863,8 +863,8 @@ void packetmath_real() { } } - data1[0] = std::numeric_limits::infinity(); - data1[1] = -std::numeric_limits::infinity(); + data1[0] = NumTraits::infinity(); + data1[1] = -NumTraits::infinity(); h.store(data2, internal::psin(h.load(data1))); VERIFY((numext::isnan)(data2[0])); VERIFY((numext::isnan)(data2[1])); @@ -873,7 +873,7 @@ void packetmath_real() { VERIFY((numext::isnan)(data2[0])); VERIFY((numext::isnan)(data2[1])); - data1[0] = std::numeric_limits::quiet_NaN(); + data1[0] = NumTraits::quiet_NaN(); h.store(data2, internal::psin(h.load(data1))); VERIFY((numext::isnan)(data2[0])); h.store(data2, internal::pcos(h.load(data1))); @@ -997,13 +997,13 @@ void packetmath_notcomplex() { VERIFY(internal::isApprox(ref[0], internal::predux_max(internal::pload(data1))) && "internal::predux_max"); // A single NaN. const size_t index = std::numeric_limits::quiet_NaN() % PacketSize; - data1[index] = std::numeric_limits::quiet_NaN(); + data1[index] = NumTraits::quiet_NaN(); VERIFY(PacketSize==1 || !(numext::isnan)(internal::predux_min(internal::pload(data1)))); VERIFY((numext::isnan)(internal::predux_min(internal::pload(data1)))); VERIFY(PacketSize==1 || !(numext::isnan)(internal::predux_max(internal::pload(data1)))); VERIFY((numext::isnan)(internal::predux_max(internal::pload(data1)))); // All NaNs. - for (int i = 0; i < 4 * PacketSize; ++i) data1[i] = std::numeric_limits::quiet_NaN(); + for (int i = 0; i < 4 * PacketSize; ++i) data1[i] = NumTraits::quiet_NaN(); VERIFY((numext::isnan)(internal::predux_min(internal::pload(data1)))); VERIFY((numext::isnan)(internal::predux_min(internal::pload(data1)))); VERIFY((numext::isnan)(internal::predux_max(internal::pload(data1)))); @@ -1011,8 +1011,8 @@ void packetmath_notcomplex() { // Test NaN propagation for coefficient-wise min and max. for (int i = 0; i < PacketSize; ++i) { - data1[i] = internal::random() ? std::numeric_limits::quiet_NaN() : Scalar(0); - data1[i + PacketSize] = internal::random() ? std::numeric_limits::quiet_NaN() : Scalar(0); + data1[i] = internal::random() ? NumTraits::quiet_NaN() : Scalar(0); + data1[i + PacketSize] = internal::random() ? NumTraits::quiet_NaN() : Scalar(0); } // Note: NaN propagation is implementation defined for pmin/pmax, so we do not test it here. CHECK_CWISE2_IF(PacketTraits::HasMin, propagate_number_min, (internal::pmin)); -- cgit v1.2.3