Skip to content

Add single source of truth for absolute and relative differences to Mantid::Kernel::FloatingPointComparison #37931

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 3 commits into from
Sep 17, 2024
Merged

Add single source of truth for absolute and relative differences to Mantid::Kernel::FloatingPointComparison #37931

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
03fb669
add absolute and relative difference operators plus tests
rboston628 Sep 6, 2024
c615c2b
handle NaNs according to PEP 485
rboston628 Sep 13, 2024
65e39a1
move release notes
rboston628 Sep 13, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 11 additions & 3 deletions Framework/Kernel/inc/MantidKernel/FloatingPointComparison.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -14,10 +14,18 @@
namespace Mantid {
namespace Kernel {
/// Test for equality of doubles using compiler-defined precision
template <typename T> MANTID_KERNEL_DLL bool equals(const T x, const T y);
template <typename T> MANTID_KERNEL_DLL bool equals(T const x, T const y);
/// Test whether x<=y within machine precision
template <typename T> MANTID_KERNEL_DLL bool ltEquals(const T x, const T y);
template <typename T> MANTID_KERNEL_DLL bool ltEquals(T const x, T const y);
/// Test whether x>=y within machine precision
template <typename T> MANTID_KERNEL_DLL bool gtEquals(const T x, const T y);
template <typename T> MANTID_KERNEL_DLL bool gtEquals(T const x, T const y);
/// Calculate absolute difference between x, y
template <typename T> MANTID_KERNEL_DLL T absoluteDifference(T const x, T const y);
/// Calculate relative difference between x, y
template <typename T> MANTID_KERNEL_DLL T relativeDifference(T const x, T const y);
/// Test whether x, y are within absolute tolerance tol
template <typename T> MANTID_KERNEL_DLL bool withinAbsoluteDifference(T const x, T const y, T const tolerance);
/// Test whether x, y are within relative tolerance tol
template <typename T> MANTID_KERNEL_DLL bool withinRelativeDifference(T const x, T const y, T const tolerance);
} // namespace Kernel
} // namespace Mantid
110 changes: 98 additions & 12 deletions Framework/Kernel/src/FloatingPointComparison.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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 <typename TYPE> bool equals(const TYPE x, const TYPE y) {
return !(std::fabs(x - y) > std::numeric_limits<TYPE>::epsilon());
}
template <typename T> bool equals(T const x, T const y) { return std::abs(x - y) <= std::numeric_limits<T>::epsilon(); }

/**
* Compare two floating-point numbers as to whether they satisfy x<=y within
Expand All @@ -35,7 +33,7 @@ template <typename TYPE> bool equals(const TYPE x, const TYPE y) {
* @returns True if the numbers are considered <= within the machine tolerance,
* false otherwise
*/
template <typename T> MANTID_KERNEL_DLL bool ltEquals(const T x, const T y) { return (equals(x, y) || x < y); }
template <typename T> 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
Expand All @@ -45,15 +43,103 @@ template <typename T> 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 <typename T> MANTID_KERNEL_DLL bool gtEquals(const T x, const T y) { return (equals(x, y) || x > y); }
template <typename T> 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 <typename T> 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 <typename T> T relativeDifference(T const x, T const y) {
// calculate numerator |x-y|
T const num = absoluteDifference<T>(x, y);
if (num <= std::numeric_limits<T>::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<T>(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 <typename T> bool withinAbsoluteDifference(T const x, T const y, T const tolerance) {
return ltEquals(absoluteDifference<T>(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 <typename T> 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<T>(x, y);
if (!(num > std::numeric_limits<T>::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<T>(0.5 * (std::abs(x) + std::abs(y)));
Copy link
Member

@peterfpeterson peterfpeterson Sep 16, 2024 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In the case of the two values being on opposite sides of zero, this gives a very different value than the average of the two

    T const denom = static_cast<T>(0.5 * std::abs(x + y));

Is that intentional?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It does lead to different values. I'm pretty convinced that the average of the absolute values is the correct norm to use. I am following what is used, for instance, in TableColumn,h. It's also safer.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Another option is using the max of the absolute values. That is what is used in python's math.isclose(). In the case of a positive and negative, the max and average absolute value will remain similar, while the absolute value of the average could be zero.

E.g. x=-2, y=3.

method rel. diff.
max(|x|,|y|) 5/3 = 2.66
avg(|x|, |y|) 5/2.5 = 2
|avg(x, y)| 5/1 = 5

E.g. x=-10, y=10.

method rel. diff.
max(|x|,|y|) 20/10 = 2
avg(|x|, |y|) 20/10 = 2
|avg(x, y)| 20/0 = +inf

I'm willing to try the max instead, but that might change some expected behavior

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I feel sufficiently schooled.

// 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<double>(const double, const double);
template DLLExport bool equals<float>(const float, const float);
template DLLExport bool ltEquals<double>(const double, const double);
template DLLExport bool ltEquals<float>(const float, const float);
template DLLExport bool gtEquals<double>(const double, const double);
template DLLExport bool gtEquals<float>(const float, const float);
template DLLExport bool equals<double>(double const, double const);
template DLLExport bool equals<float>(float const, float const);
template DLLExport bool ltEquals<double>(double const, double const);
template DLLExport bool ltEquals<float>(float const, float const);
template DLLExport bool gtEquals<double>(double const, double const);
template DLLExport bool gtEquals<float>(float const, float const);
//
template DLLExport double absoluteDifference<double>(double const, double const);
template DLLExport float absoluteDifference<float>(float const, float const);
template DLLExport double relativeDifference<double>(double const, double const);
template DLLExport float relativeDifference<float>(float const, float const);
//
template DLLExport bool withinAbsoluteDifference<double>(double const, double const, double const);
template DLLExport bool withinAbsoluteDifference<float>(float const, float const, float const);
template DLLExport bool withinRelativeDifference<double>(double const, double const, double const);
template DLLExport bool withinRelativeDifference<float>(float const, float const, float const);
///@endcond

} // namespace Mantid::Kernel
114 changes: 114 additions & 0 deletions Framework/Kernel/test/FloatingPointComparisonTest.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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<double>::quiet_NaN();
constexpr double bnan = std::numeric_limits<double>::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);
Expand Down Expand Up @@ -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<double>::quiet_NaN();
constexpr double bnan = std::numeric_limits<double>::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<double>::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<double>::quiet_NaN();
constexpr double bnan = std::numeric_limits<double>::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<double>::quiet_NaN();
constexpr double bnan = std::numeric_limits<double>::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<double>::quiet_NaN();
constexpr double bnan = std::numeric_limits<double>::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);
}
};
1 change: 1 addition & 0 deletions docs/source/release/v6.12.0/Framework/Algorithms/New_features/37931.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
- add new functions for efficiently calculating absolute and relative differences
Loading