From c06840fe4e1caa5e40c8f8e8d236b9720d03b7f2 Mon Sep 17 00:00:00 2001 From: Erik Lindahl Date: Sun, 31 Jul 2016 18:34:36 +0200 Subject: [PATCH] Fix proper relative tolerance checking in testutils The previous relative tolerances were convereted to ULP checks, but since the relative difference for 1 ULP depends on whether the exponent changes there is no strict 1-to-1 correspondence between ULP and relative differences. Fixed by storing the original magnitude of the numbers in FloatingPointDifference and using a traditional relative tolerance check for the case when a relative tolerance is set from floating point values; other checks still use ULP. Change-Id: I95331b46da00cd0b4e369f087991c553cd05f18f --- src/testutils/testasserts.cpp | 67 +++++++++++++++++++++++++++++++--------- src/testutils/testasserts.h | 71 ++++++++++++++++++++++++++++++++----------- 2 files changed, 107 insertions(+), 31 deletions(-) diff --git a/src/testutils/testasserts.cpp b/src/testutils/testasserts.cpp index 4d9b499e9e..a2c515120e 100644 --- a/src/testutils/testasserts.cpp +++ b/src/testutils/testasserts.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2014,2015, by the GROMACS development team, led by + * Copyright (c) 2014,2015,2016, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -193,16 +193,18 @@ gmx_uint64_t relativeToleranceToUlp(FloatType tolerance) * FloatingPointDifference */ -FloatingPointDifference::FloatingPointDifference(float value1, float value2) +FloatingPointDifference::FloatingPointDifference(float ref, float value) + : termMagnitude_(std::abs(ref)) { - initDifference(value1, value2, + initDifference(ref, value, &absoluteDifference_, &ulpDifference_, &bSignDifference_); bDouble_ = false; } -FloatingPointDifference::FloatingPointDifference(double value1, double value2) +FloatingPointDifference::FloatingPointDifference(double ref, double value) + : termMagnitude_(std::abs(ref)) { - initDifference(value1, value2, + initDifference(ref, value, &absoluteDifference_, &ulpDifference_, &bSignDifference_); bDouble_ = true; } @@ -214,11 +216,29 @@ bool FloatingPointDifference::isNaN() const std::string FloatingPointDifference::toString() const { - const double eps = isDouble() ? GMX_DOUBLE_EPS : GMX_FLOAT_EPS; - return formatString("%g (%" GMX_PRIu64 " %s-prec. ULPs, rel. %.3g)%s", + std::string relDiffStr; + + if (termMagnitude_ > 0) + { + // If the reference value is finite we calculate the proper quotient + relDiffStr = formatString("%.3g", std::abs(absoluteDifference_/termMagnitude_)); + } + else if (absoluteDifference_ == 0.0) + { + // If the numbers are identical the quotient is strictly NaN here, but + // there no reason to worry when we have a perfect match. + relDiffStr = formatString("%.3g", 0.0); + } + else + { + // If the reference value is zero and numbers are non-identical, relative difference is infinite. + relDiffStr = formatString("Inf"); + } + + return formatString("%g (%" GMX_PRIu64 " %s-prec. ULPs, rel. %s)%s", absoluteDifference_, ulpDifference_, isDouble() ? "double" : "single", - ulpDifference_ * eps, + relDiffStr.c_str(), bSignDifference_ ? ", signs differ" : ""); } @@ -246,12 +266,23 @@ bool FloatingPointTolerance::isWithin( return true; } + // By using smaller-than-or-equal below, we allow the test to pass if + // the numbers are identical, even if the term magnitude is 0, which seems + // a reasonable thing to do... + const double relativeTolerance + = difference.isDouble() ? doubleRelativeTolerance_ : singleRelativeTolerance_; + if (difference.asAbsolute() <= relativeTolerance * difference.termMagnitude()) + { + return true; + } + const gmx_uint64_t ulpTolerance = difference.isDouble() ? doubleUlpTolerance_ : singleUlpTolerance_; if (ulpTolerance < GMX_UINT64_MAX && difference.asUlps() <= ulpTolerance) { return true; } + return false; } @@ -260,22 +291,30 @@ std::string FloatingPointTolerance::toString(const FloatingPointDifference &diff std::string result; const double absoluteTolerance = difference.isDouble() ? doubleAbsoluteTolerance_ : singleAbsoluteTolerance_; + const double relativeTolerance + = difference.isDouble() ? doubleRelativeTolerance_ : singleRelativeTolerance_; const gmx_uint64_t ulpTolerance = difference.isDouble() ? doubleUlpTolerance_ : singleUlpTolerance_; - const double eps - = difference.isDouble() ? GMX_DOUBLE_EPS : GMX_FLOAT_EPS; if (absoluteTolerance > 0.0) { result.append(formatString("abs. %g", absoluteTolerance)); } + if (relativeTolerance > 0.0) + { + if (!result.empty()) + { + result.append(", "); + } + result.append(formatString("rel. %.3g", relativeTolerance)); + } if (ulpTolerance < GMX_UINT64_MAX) { if (!result.empty()) { result.append(", "); } - result.append(formatString("%" GMX_PRIu64 " ULPs (rel. %.3g)", ulpTolerance, ulpTolerance * eps)); + result.append(formatString("%" GMX_PRIu64 " ULPs", ulpTolerance)); } if (bSignMustMatch_) { @@ -293,10 +332,10 @@ std::string FloatingPointTolerance::toString(const FloatingPointDifference &diff FloatingPointTolerance relativeToleranceAsFloatingPoint(double magnitude, double tolerance) { - const double absoluteTolerance = magnitude * tolerance; + const double absoluteTolerance = std::abs(magnitude) * tolerance; return FloatingPointTolerance(absoluteTolerance, absoluteTolerance, - relativeToleranceToUlp(tolerance), - relativeToleranceToUlp(tolerance), + tolerance, tolerance, + GMX_UINT64_MAX, GMX_UINT64_MAX, false); } //! \endcond diff --git a/src/testutils/testasserts.h b/src/testutils/testasserts.h index dc4b68bba2..dcdaf17db4 100644 --- a/src/testutils/testasserts.h +++ b/src/testutils/testasserts.h @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2013,2014,2015, by the GROMACS development team, led by + * Copyright (c) 2013,2014,2015,2016, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -209,10 +209,30 @@ void processExpectedException(const std::exception &ex); class FloatingPointDifference { public: - //! Initializes a single-precision difference. - FloatingPointDifference(float value1, float value2); - //! Initializes a double-precision difference. - FloatingPointDifference(double value1, double value2); + + /*! \brief Initializes a single-precision difference. + * + * \param ref First term in difference + * \param value Second term in difference + * + * For absolute and ULP differences the two parameters are equivalent, + * since the difference is symmetric. For relative differences + * the first term is interpreted as the reference value, from which + * we extract the magnitude to compare with. + */ + FloatingPointDifference(float ref, float value); + + /*! \brief Initializes a double-precision difference. + * + * \param ref First term in difference + * \param value Second term in difference + * + * For absolute and ULP differences the two parameters are equivalent, + * since the difference is symmetric. For relative differences + * the first term is interpreted as the reference value, from which + * we extract the magnitude to compare with. + */ + FloatingPointDifference(double ref, double value); /*! \brief * Whether one or both of the compared values were NaN. @@ -244,7 +264,13 @@ class FloatingPointDifference //! Formats the difference as a string for assertion failure messages. std::string toString() const; + //! Returns the magnitude of the original second term of the difference. + double termMagnitude() const { return termMagnitude_; } + private: + + //! Save the magnitude of the reference value for relative (i.e., not ULP) tolerance + double termMagnitude_; //! Stores the absolute difference, or NaN if one or both values were NaN. double absoluteDifference_; gmx_uint64_t ulpDifference_; @@ -270,6 +296,9 @@ class FloatingPointDifference * the given tolerance for the check to pass. * Setting the absolute tolerance to zero disables the absolute tolerance * check. + * - _relative tolerance_: the absolute difference between the numbers must + * be smaller than the tolerance multiplied by the first number. Setting + * the relative tolerance to zero disables this check. * - _ULP tolerance_: ULP (units of least precision) difference between the * values must be smaller than the given tolerance for the check to pass. * Setting the ULP tolerance to zero requires exact match. @@ -279,13 +308,13 @@ class FloatingPointDifference * check (note that this also applies to `0.0` and `-0.0`: a value with a * different sign than the zero will fail the check). * - * Either an absolute or a ULP tolerance must always be specified. - * If both are specified, then the check passes if either of the tolerances is - * satisfied. + * Either an absolute, relative, or ULP tolerance must always be specified. + * If several of them are specified, then the check passes if either of the + * tolerances is satisfied. * - * Any combination of absolute and ULP tolerance can be combined with the sign - * check. In this case, the sign check must succeed for the check to pass, - * even if other tolerances are satisfied. + * Any combination of absolute, relative, and ULP tolerance can be combined with + * the sign check. In this case, the sign check must succeed for the check to + * pass, even if other tolerances are satisfied. * * The tolerances can be specified separately for single and double precision * comparison. Different initialization functions have different semantics on @@ -312,6 +341,10 @@ class FloatingPointTolerance * Allowed absolute difference in a single-precision number. * \param[in] doubleAbsoluteTolerance * Allowed absolute difference in a double-precision number. + * \param[in] singleRelativeTolerance + * Allowed relative difference in a single-precision number. + * \param[in] doubleRelativeTolerance + * Allowed relative difference in a double-precision number. * \param[in] singleUlpTolerance * Allowed ULP difference in a single-precision number. * \param[in] doubleUlpTolerance @@ -321,11 +354,15 @@ class FloatingPointTolerance */ FloatingPointTolerance(float singleAbsoluteTolerance, double doubleAbsoluteTolerance, + float singleRelativeTolerance, + double doubleRelativeTolerance, gmx_uint64_t singleUlpTolerance, gmx_uint64_t doubleUlpTolerance, bool bSignMustMatch) : singleAbsoluteTolerance_(singleAbsoluteTolerance), doubleAbsoluteTolerance_(doubleAbsoluteTolerance), + singleRelativeTolerance_(singleRelativeTolerance), + doubleRelativeTolerance_(doubleRelativeTolerance), singleUlpTolerance_(singleUlpTolerance), doubleUlpTolerance_(doubleUlpTolerance), bSignMustMatch_(bSignMustMatch) @@ -345,6 +382,8 @@ class FloatingPointTolerance private: float singleAbsoluteTolerance_; double doubleAbsoluteTolerance_; + float singleRelativeTolerance_; + double doubleRelativeTolerance_; gmx_uint64_t singleUlpTolerance_; gmx_uint64_t doubleUlpTolerance_; bool bSignMustMatch_; @@ -361,7 +400,7 @@ class FloatingPointTolerance static inline FloatingPointTolerance ulpTolerance(gmx_uint64_t ulpDiff) { - return FloatingPointTolerance(0.0, 0.0, ulpDiff, ulpDiff, false); + return FloatingPointTolerance(0.0, 0.0, 0.0, 0.0, ulpDiff, ulpDiff, false); } /*! \brief @@ -371,7 +410,7 @@ ulpTolerance(gmx_uint64_t ulpDiff) * \param[in] magnitude Magnitude of the numbers the computation operates in. * \param[in] tolerance Relative tolerance permitted (e.g. 1e-4). * - * In addition to setting an ULP tolerance equivalent to \p tolerance for both + * In addition to setting an relative tolerance for both * precisions, this sets the absolute tolerance such that values close to zero * (in general, smaller than \p magnitude) do not fail the check if they * differ by less than \p tolerance evaluated at \p magnitude. This accounts @@ -379,9 +418,6 @@ ulpTolerance(gmx_uint64_t ulpDiff) * accuracy of values much less than \p magnitude do not matter for * correctness. * - * The ULP tolerance for different precisions will be different to make them - * both match \p tolerance. - * * \related FloatingPointTolerance */ FloatingPointTolerance @@ -412,6 +448,7 @@ relativeToleranceAsPrecisionDependentUlp(double magnitude, { return FloatingPointTolerance(magnitude*singleUlpDiff*GMX_FLOAT_EPS, magnitude*doubleUlpDiff*GMX_DOUBLE_EPS, + 0.0, 0.0, singleUlpDiff, doubleUlpDiff, false); } @@ -423,7 +460,7 @@ relativeToleranceAsPrecisionDependentUlp(double magnitude, static inline FloatingPointTolerance absoluteTolerance(double tolerance) { - return FloatingPointTolerance(tolerance, tolerance, + return FloatingPointTolerance(tolerance, tolerance, 0.0, 0.0, GMX_UINT64_MAX, GMX_UINT64_MAX, false); } -- 2.11.4.GIT