diff --git a/base/base_tests/math_test.cpp b/base/base_tests/math_test.cpp index 12f22b96ea..16eb1e5f9c 100644 --- a/base/base_tests/math_test.cpp +++ b/base/base_tests/math_test.cpp @@ -64,6 +64,22 @@ UNIT_TEST(AlmostEqualULPs_Smoke) TEST(!my::AlmostEqualULPs(0.0, eps), ()); } +UNIT_TEST(AlmostEqual_Smoke) +{ + double const small = 1e-18; + double const eps = 1e-10; + + TEST(my::AlmostEqualAbs(0.0, 0.0 + small, eps), ()); + TEST(!my::AlmostEqualRel(0.0, 0.0 + small, eps), ()); + TEST(!my::AlmostEqualULPs(0.0, 0.0 + small), ()); + + TEST(my::AlmostEqualAbs(1.0, 1.0 + small, eps), ()); + TEST(my::AlmostEqualRel(1.0, 1.0 + small, eps), ()); + TEST(my::AlmostEqualULPs(1.0, 1.0 + small), ()); + + TEST(my::AlmostEqualRel(123456789.0, 123456780.0, 1e-7), ()); +} + namespace { diff --git a/base/math.hpp b/base/math.hpp index b5cec68fd8..d8ecd2418c 100644 --- a/base/math.hpp +++ b/base/math.hpp @@ -21,19 +21,21 @@ template inline T Abs(T x) // maxULPs - number of closest floating point values that are considered equal. // Infinity is treated as almost equal to the largest possible floating point values. // NaN produces undefined result. -// See http://www.cygnus-software.com/papers/comparingfloats/comparingfloats.htm for details. -template bool AlmostEqualULPs(FloatT x, FloatT y, unsigned int maxULPs = 256) +// See https://randomascii.wordpress.com/2012/02/25/comparing-floating-point-numbers-2012-edition/ +// for details. +template +bool AlmostEqualULPs(TFloat x, TFloat y, unsigned int maxULPs = 256) { - STATIC_ASSERT(is_floating_point::value); - STATIC_ASSERT(numeric_limits::is_iec559); - STATIC_ASSERT(!numeric_limits::is_exact); - STATIC_ASSERT(!numeric_limits::is_integer); + STATIC_ASSERT(is_floating_point::value); + STATIC_ASSERT(numeric_limits::is_iec559); + STATIC_ASSERT(!numeric_limits::is_exact); + STATIC_ASSERT(!numeric_limits::is_integer); // Make sure maxUlps is non-negative and small enough that the - // default NAN won't compare as equal to anything. + // default NaN won't compare as equal to anything. ASSERT_LESS(maxULPs, 4 * 1024 * 1024, ()); - int const bits = 8 * sizeof(FloatT); + int const bits = 8 * sizeof(TFloat); typedef typename boost::int_t::exact IntType; typedef typename boost::uint_t::exact UIntType; @@ -49,7 +51,27 @@ template bool AlmostEqualULPs(FloatT x, FloatT y, unsigned int UIntType const diff = Abs(xInt - yInt); - return (diff <= maxULPs); + return diff <= maxULPs; +} + +// Returns true if x and y are equal up to the absolute difference eps. +// Does not produce a sensible result if any of the arguments is NaN or infinity. +// The default value for eps is deliberately not provided: the intended usage +// is for the client to choose the precision according to the problem domain, +// explicitly define the precision constant and call this function. +template +inline bool AlmostEqualAbs(TFloat x, TFloat y, TFloat eps) +{ + return fabs(x - y) < eps; +} + +// Returns true if x and y are equal up to the relative difference eps. +// Does not produce a sensible result if any of the arguments is NaN, infinity or zero. +// The same considerations as in AlmostEqualAbs apply. +template +inline bool AlmostEqualRel(TFloat x, TFloat y, TFloat eps) +{ + return fabs(x - y) < eps * max(fabs(x), fabs(y)); } template inline TFloat DegToRad(TFloat deg)