diff --git a/Framework/Kernel/inc/MantidKernel/FloatingPointComparison.h b/Framework/Kernel/inc/MantidKernel/FloatingPointComparison.h index 72810f09bb0a..eb4ee09791fc 100644 --- a/Framework/Kernel/inc/MantidKernel/FloatingPointComparison.h +++ b/Framework/Kernel/inc/MantidKernel/FloatingPointComparison.h @@ -14,10 +14,18 @@ namespace Mantid { namespace Kernel { /// Test for equality of doubles using compiler-defined precision -template MANTID_KERNEL_DLL bool equals(const T x, const T y); +template MANTID_KERNEL_DLL bool equals(T const x, T const y); /// Test whether x<=y within machine precision -template MANTID_KERNEL_DLL bool ltEquals(const T x, const T y); +template MANTID_KERNEL_DLL bool ltEquals(T const x, T const y); /// Test whether x>=y within machine precision -template MANTID_KERNEL_DLL bool gtEquals(const T x, const T y); +template MANTID_KERNEL_DLL bool gtEquals(T const x, T const y); +/// Calculate absolute difference between x, y +template MANTID_KERNEL_DLL T absoluteDifference(T const x, T const y); +/// Calculate relative difference between x, y +template MANTID_KERNEL_DLL T relativeDifference(T const x, T const y); +/// Test whether x, y are within absolute tolerance tol +template MANTID_KERNEL_DLL bool withinAbsoluteDifference(T const x, T const y, T const tolerance); +/// Test whether x, y are within relative tolerance tol +template MANTID_KERNEL_DLL bool withinRelativeDifference(T const x, T const y, T const tolerance); } // namespace Kernel } // namespace Mantid diff --git a/Framework/Kernel/src/FloatingPointComparison.cpp b/Framework/Kernel/src/FloatingPointComparison.cpp index f3cb9888b97a..f2c15a5492e1 100644 --- a/Framework/Kernel/src/FloatingPointComparison.cpp +++ b/Framework/Kernel/src/FloatingPointComparison.cpp @@ -21,11 +21,9 @@ namespace Mantid::Kernel { * @param x :: LHS comparator * @param y :: RHS comparator * @returns True if the numbers are considered equal within the given tolerance, - * false otherwise + * false otherwise. False if any value is NaN. */ -template bool equals(const TYPE x, const TYPE y) { - return !(std::fabs(x - y) > std::numeric_limits::epsilon()); -} +template bool equals(T const x, T const y) { return std::abs(x - y) <= std::numeric_limits::epsilon(); } /** * Compare two floating-point numbers as to whether they satisfy x<=y within @@ -35,7 +33,7 @@ template bool equals(const TYPE x, const TYPE y) { * @returns True if the numbers are considered <= within the machine tolerance, * false otherwise */ -template MANTID_KERNEL_DLL bool ltEquals(const T x, const T y) { return (equals(x, y) || x < y); } +template MANTID_KERNEL_DLL bool ltEquals(T const x, T const y) { return (equals(x, y) || x < y); } /** * Compare two floating-point numbers as to whether they satisfy x>=y within @@ -45,15 +43,103 @@ template MANTID_KERNEL_DLL bool ltEquals(const T x, const T y) { re * @returns True if the numbers are considered <= within the machine tolerance, * false otherwise */ -template MANTID_KERNEL_DLL bool gtEquals(const T x, const T y) { return (equals(x, y) || x > y); } +template MANTID_KERNEL_DLL bool gtEquals(T const x, T const y) { return (equals(x, y) || x > y); } + +/// ABSOLUTE AND RELATIVE DIFFERENCE + +/** + * Calculate the absolute difference of two floating-point numbers + * @param x :: first value + * @param y :: second value + * @returns the value of the absolute difference + */ +template T absoluteDifference(T const x, T const y) { return std::abs(x - y); } + +/** + * Calculate the relative difference of two floating-point numbers + * @param x :: first value + * @param y :: second value + * @returns the value of the relative difference. Do NOT use this + * to compare the result to a tolerance; for that, use withinRelativeDifference + * instead, as it will be more efficient at the comparison. + */ +template T relativeDifference(T const x, T const y) { + // calculate numerator |x-y| + T const num = absoluteDifference(x, y); + if (num <= std::numeric_limits::epsilon()) { + // if |x-y| == 0.0 (within machine tolerance), relative difference is zero + return 0.0; + } else { + // otherwise we have to calculate the denominator + T const denom = static_cast(0.5 * (std::abs(x) + std::abs(y))); + // NOTE if we made it this far, at least one of x or y is nonzero, so denom will be nonzero + return num / denom; + } +} + +/** + * Compare floating point numbers for absolute difference to + * within the given tolerance + * @param x :: first value + * @param y :: second value + * @param tolerance :: the tolerance + * @returns True if the numbers are considered equal within the given tolerance, + * false otherwise. False if either value is NaN. + */ +template bool withinAbsoluteDifference(T const x, T const y, T const tolerance) { + return ltEquals(absoluteDifference(x, y), tolerance); +} + +/** + * Compare floating point numbers for relative difference to + * within the given tolerance + * @param x :: first value + * @param y :: second value + * @param tolerance :: the tolerance + * @returns True if the numbers are considered equal within the given tolerance, + * false otherwise. False if either value is NaN. + */ +template bool withinRelativeDifference(T const x, T const y, T const tolerance) { + // handle the case of NaNs + if (std::isnan(x) || std::isnan(y)) { + return false; + } + T const num = absoluteDifference(x, y); + if (!(num > std::numeric_limits::epsilon())) { + // if |x-y| == 0.0 (within machine tolerance), this test passes + return true; + } else { + // otherwise we have to calculate the denominator + T const denom = static_cast(0.5 * (std::abs(x) + std::abs(y))); + // if denom <= 1, then |x-y| > tol implies |x-y|/denom > tol, can return early + // NOTE can only return early if BOTH denom > tol AND |x-y| > tol. + if (denom <= 1. && !ltEquals(num, tolerance)) { + return false; + } else { + // avoid division for improved performance + return ltEquals(num, denom * tolerance); + } + } +} ///@cond // Concrete instantiations -template DLLExport bool equals(const double, const double); -template DLLExport bool equals(const float, const float); -template DLLExport bool ltEquals(const double, const double); -template DLLExport bool ltEquals(const float, const float); -template DLLExport bool gtEquals(const double, const double); -template DLLExport bool gtEquals(const float, const float); +template DLLExport bool equals(double const, double const); +template DLLExport bool equals(float const, float const); +template DLLExport bool ltEquals(double const, double const); +template DLLExport bool ltEquals(float const, float const); +template DLLExport bool gtEquals(double const, double const); +template DLLExport bool gtEquals(float const, float const); +// +template DLLExport double absoluteDifference(double const, double const); +template DLLExport float absoluteDifference(float const, float const); +template DLLExport double relativeDifference(double const, double const); +template DLLExport float relativeDifference(float const, float const); +// +template DLLExport bool withinAbsoluteDifference(double const, double const, double const); +template DLLExport bool withinAbsoluteDifference(float const, float const, float const); +template DLLExport bool withinRelativeDifference(double const, double const, double const); +template DLLExport bool withinRelativeDifference(float const, float const, float const); ///@endcond + } // namespace Mantid::Kernel diff --git a/Framework/Kernel/test/FloatingPointComparisonTest.h b/Framework/Kernel/test/FloatingPointComparisonTest.h index f6072eea7c61..9fff197fb606 100644 --- a/Framework/Kernel/test/FloatingPointComparisonTest.h +++ b/Framework/Kernel/test/FloatingPointComparisonTest.h @@ -30,6 +30,25 @@ class FloatingPointComparisonTest : public CxxTest::TestSuite { TS_ASSERT_EQUALS(Mantid::Kernel::equals(0.1, 1.0001 * tol), false); } + void test_with_NaN() { + // everything compares false with an NaN + constexpr double anan = std::numeric_limits::quiet_NaN(); + constexpr double bnan = std::numeric_limits::quiet_NaN(); + constexpr double real = 3.0; + // equals + TS_ASSERT_EQUALS(Mantid::Kernel::equals(anan, real), false); + TS_ASSERT_EQUALS(Mantid::Kernel::equals(real, anan), false); + TS_ASSERT_EQUALS(Mantid::Kernel::equals(anan, bnan), false); + // ltEquals + TS_ASSERT_EQUALS(Mantid::Kernel::ltEquals(anan, real), false); + TS_ASSERT_EQUALS(Mantid::Kernel::ltEquals(real, anan), false); + TS_ASSERT_EQUALS(Mantid::Kernel::ltEquals(anan, bnan), false); + // gtEquals + TS_ASSERT_EQUALS(Mantid::Kernel::gtEquals(anan, real), false); + TS_ASSERT_EQUALS(Mantid::Kernel::gtEquals(real, anan), false); + TS_ASSERT_EQUALS(Mantid::Kernel::gtEquals(anan, bnan), false); + } + void test_LtEquals_With_X_Equal_To_Y_Produces_True() { TS_ASSERT_EQUALS(Mantid::Kernel::ltEquals(0.1, 0.1), true); TS_ASSERT_EQUALS(Mantid::Kernel::ltEquals(-0.1, -0.1), true); @@ -63,4 +82,99 @@ class FloatingPointComparisonTest : public CxxTest::TestSuite { TS_ASSERT_EQUALS(Mantid::Kernel::gtEquals(-5.56, 0.23), false); TS_ASSERT_EQUALS(Mantid::Kernel::gtEquals(-0.00002, -0.00001), false); } + + void test_absoluteDifference() { + constexpr double left = 1.1, right = 1.0; + // test value + TS_ASSERT_EQUALS(Mantid::Kernel::absoluteDifference(left, right), std::abs(left - right)); + // test positive-definiteness + TS_ASSERT_LESS_THAN(0.0, Mantid::Kernel::absoluteDifference(left, -right)); + TS_ASSERT_LESS_THAN(0.0, Mantid::Kernel::absoluteDifference(-left, right)); + TS_ASSERT_LESS_THAN(0.0, Mantid::Kernel::absoluteDifference(-left, -right)); + // test symmetry + TS_ASSERT_EQUALS(Mantid::Kernel::absoluteDifference(left, right), Mantid::Kernel::absoluteDifference(right, left)); + // absolute difference with NaN is NaN + constexpr double anan = std::numeric_limits::quiet_NaN(); + constexpr double bnan = std::numeric_limits::quiet_NaN(); + TS_ASSERT(std::isnan(Mantid::Kernel::absoluteDifference(left, anan))); + TS_ASSERT(std::isnan(Mantid::Kernel::absoluteDifference(bnan, anan))); + } + + void test_relativeDifference() { + constexpr double point3 = 0.3, notquitepoint3 = 0.2 + 0.1; + TS_ASSERT_EQUALS(Mantid::Kernel::relativeDifference(point3, notquitepoint3), 0.0); + TS_ASSERT_EQUALS(Mantid::Kernel::relativeDifference(2.3, 2.3), 0.0); + TS_ASSERT_EQUALS(Mantid::Kernel::relativeDifference(2.3e208, 2.3e208), 0.0); + // check no errors using zero + TS_ASSERT_THROWS_NOTHING(Mantid::Kernel::relativeDifference(0.0, 0.0)) + TS_ASSERT(!std::isnan(Mantid::Kernel::relativeDifference(0.0, 0.0))); + TS_ASSERT_EQUALS(Mantid::Kernel::relativeDifference(0.0, 0.0), 0.0); + // check no errors using machine epsilon + constexpr double realsmall = std::numeric_limits::epsilon(); + TS_ASSERT_THROWS_NOTHING(Mantid::Kernel::relativeDifference(0.0, realsmall)) + TS_ASSERT(!std::isnan(Mantid::Kernel::relativeDifference(0.0, realsmall))); + TS_ASSERT_EQUALS(Mantid::Kernel::relativeDifference(0.0, realsmall), 0.0); + // check we get correct values for normal situations + const double left = 2.6, right = 2.7; + const double reldiff = 2.0 * std::abs(left - right) / (left + right); + TS_ASSERT_EQUALS(Mantid::Kernel::relativeDifference(left, right), reldiff); + TS_ASSERT_EQUALS(Mantid::Kernel::relativeDifference(right, left), reldiff); + // relative difference with NaN is NaN + constexpr double anan = std::numeric_limits::quiet_NaN(); + constexpr double bnan = std::numeric_limits::quiet_NaN(); + TS_ASSERT(std::isnan(Mantid::Kernel::relativeDifference(left, anan))); + TS_ASSERT(std::isnan(Mantid::Kernel::absoluteDifference(bnan, anan))); + } + + void test_withinAbsoluteDifference() { + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(0.3, 0.2, 0.1), true); + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(0.3, 0.1, 0.1), false); + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(0.01, 0.011, 0.01), true); + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(0.01, -0.011, 0.01), false); + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(100.1, 100.15, 0.1), true); + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(12345678923456.789, 12345679023456.788, 0.0001), false); + // case of NaNs -- nothing is close to an NaN + constexpr double anan = std::numeric_limits::quiet_NaN(); + constexpr double bnan = std::numeric_limits::quiet_NaN(); + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(anan, 0.3, 0.1), false); + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(anan, bnan, 0.1), false); + } + + void test_withinRelativeDifference() { + // things difference at machine epsilon are equal + const double point3 = 0.3, notquitepoint3 = 0.2 + 0.1; + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(point3, notquitepoint3, 1.e-307), true); + // some cases with zero difference + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(2.3, 2.3, 1.e-307), true); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(2.3e208, 2.3e208, 1.e-307), true); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(2.3e-208, 2.3e-208, 0.0), true); + // case of large magnitude values -- even though abs diff would always fail, rel diff can still pass + // - passing + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(2.31e208, 2.32e208, 0.01), false); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(2.31e208, 2.32e208, 0.01), true); + // - failing + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(2.3e208, 2.4e208, 0.01), false); + // case of small magnitude values -- even though abs diff would always pass, rel diff still can fail + // - passing + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(2.31e-10, 2.32e-10, 0.01), true); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(2.31e-10, 2.32e-10, 0.01), true); + // - failing + TS_ASSERT_EQUALS(Mantid::Kernel::withinAbsoluteDifference(2.3e-10, 2.4e-10, 0.01), true); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(2.3e-10, 2.4e-10, 0.01), false); + // case of normal-sized values + const double left = 2.6, right = 2.7, far = 3.0; + const double reldiff = 2.0 * std::abs(left - right) / (left + right); + const double tolerance = 1.01 * reldiff; + // - passing + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(left, right, tolerance), true); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(right, left, tolerance), true); + // - failing + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(left, far, tolerance), false); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(far, left, tolerance), false); + // case of NaNs -- nothing is close to an NaN + constexpr double anan = std::numeric_limits::quiet_NaN(); + constexpr double bnan = std::numeric_limits::quiet_NaN(); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(anan, 0.3, 0.1), false); + TS_ASSERT_EQUALS(Mantid::Kernel::withinRelativeDifference(anan, bnan, 0.1), false); + } }; diff --git a/docs/source/release/v6.12.0/Framework/Algorithms/New_features/37931.rst b/docs/source/release/v6.12.0/Framework/Algorithms/New_features/37931.rst new file mode 100644 index 000000000000..fbc9a00c48c4 --- /dev/null +++ b/docs/source/release/v6.12.0/Framework/Algorithms/New_features/37931.rst @@ -0,0 +1 @@ +- add new functions for efficiently calculating absolute and relative differences \ No newline at end of file