Fix ulp distance calculation for numbers with different signs

This is a simplification of the fix proposed in #2152, with the
critical function split out so that it can be tested directly,
without having to go through the ULP matcher.

Closes #2152

src/CMakeLists.txt

@@ -65,6 +65,7 @@ set(INTERNAL_HEADERS
     ${SOURCES_DIR}/internal/catch_errno_guard.hpp
     ${SOURCES_DIR}/internal/catch_exception_translator_registry.hpp
     ${SOURCES_DIR}/internal/catch_fatal_condition_handler.hpp
+    ${SOURCES_DIR}/internal/catch_floating_point_helpers.hpp
     ${SOURCES_DIR}/internal/catch_unique_name.hpp
     ${SOURCES_DIR}/generators/catch_generator_exception.hpp
     ${SOURCES_DIR}/generators/catch_generators.hpp
@@ -160,6 +161,7 @@ set(IMPL_SOURCES
     ${SOURCES_DIR}/internal/catch_enum_values_registry.cpp
     ${SOURCES_DIR}/internal/catch_exception_translator_registry.cpp
     ${SOURCES_DIR}/internal/catch_fatal_condition_handler.cpp
+    ${SOURCES_DIR}/internal/catch_floating_point_helpers.cpp
     ${SOURCES_DIR}/generators/internal/catch_generators_combined_tu.cpp
     ${SOURCES_DIR}/interfaces/catch_interfaces_combined_tu.cpp
     ${SOURCES_DIR}/interfaces/catch_interfaces_reporter.cpp

src/catch2/catch_all.hpp

@@ -66,6 +66,7 @@
 #include <catch2/internal/catch_errno_guard.hpp>
 #include <catch2/internal/catch_exception_translator_registry.hpp>
 #include <catch2/internal/catch_fatal_condition_handler.hpp>
+#include <catch2/internal/catch_floating_point_helpers.hpp>
 #include <catch2/internal/catch_lazy_expr.hpp>
 #include <catch2/internal/catch_leak_detector.hpp>
 #include <catch2/internal/catch_list.hpp>

src/catch2/internal/catch_floating_point_helpers.cpp

@@ -0,0 +1,32 @@
+
+//              Copyright Catch2 Authors
+// Distributed under the Boost Software License, Version 1.0.
+//   (See accompanying file LICENSE_1_0.txt or copy at
+//        https://www.boost.org/LICENSE_1_0.txt)
+
+// SPDX-License-Identifier: BSL-1.0
+
+#include <catch2/internal/catch_floating_point_helpers.hpp>
+
+#include <cstring>
+
+namespace Catch {
+    namespace Detail {
+
+        uint32_t convertToBits(float f) {
+            static_assert(sizeof(float) == sizeof(uint32_t), "Important ULP matcher assumption violated");
+            uint32_t i;
+            std::memcpy(&i, &f, sizeof(f));
+            return i;
+        }
+
+        uint64_t convertToBits(double d) {
+            static_assert(sizeof(double) == sizeof(uint64_t), "Important ULP matcher assumption violated");
+            uint64_t i;
+            std::memcpy(&i, &d, sizeof(d));
+            return i;
+        }
+
+    } // end namespace Detail
+} // end namespace Catch
+

src/catch2/internal/catch_floating_point_helpers.hpp

@@ -0,0 +1,88 @@
+
+//              Copyright Catch2 Authors
+// Distributed under the Boost Software License, Version 1.0.
+//   (See accompanying file LICENSE_1_0.txt or copy at
+//        https://www.boost.org/LICENSE_1_0.txt)
+
+// SPDX-License-Identifier: BSL-1.0
+#ifndef CATCH_FLOATING_POINT_HELPERS_HPP_INCLUDED
+#define CATCH_FLOATING_POINT_HELPERS_HPP_INCLUDED
+
+#include <catch2/internal/catch_polyfills.hpp>
+
+#include <cmath>
+#include <cstdint>
+#include <utility>
+#include <limits>
+
+namespace Catch {
+    namespace Detail {
+
+        uint32_t convertToBits(float f);
+        uint64_t convertToBits(double d);
+
+    } // end namespace Detail
+
+
+    /**
+     * Calculates the ULP distance between two floating point numbers
+     *
+     * The ULP distance of two floating point numbers is the count of
+     * valid floating point numbers representable between them.
+     *
+     * There are some exceptions between how this function counts the
+     * distance, and the interpretation of the standard as implemented.
+     * by e.g. `nextafter`. For this function it always holds that:
+     * * `(x == y) => ulpDistance(x, y) == 0` (so `ulpDistance(-0, 0) == 0`)
+     * * `ulpDistance(maxFinite, INF) == 1`
+     * * `ulpDistance(x, -x) == 2 * ulpDistance(x, 0)`
+     *
+     * \pre `!isnan( lhs )`
+     * \pre `!isnan( rhs )`
+     * \pre floating point numbers are represented in IEEE-754 format
+     */
+    template <typename FP>
+    uint64_t ulpDistance( FP lhs, FP rhs ) {
+        assert( std::numeric_limits<FP>::is_iec559 &&
+            "ulpDistance assumes IEEE-754 format for floating point types" );
+        assert( !Catch::isnan( lhs ) &&
+                "Distance between NaN and number is not meaningful" );
+        assert( !Catch::isnan( rhs ) &&
+                "Distance between NaN and number is not meaningful" );
+
+        // We want X == Y to imply 0 ULP distance even if X and Y aren't
+        // bit-equal (-0 and 0), or X - Y != 0 (same sign infinities).
+        if ( lhs == rhs ) { return 0; }
+
+        // We need a properly typed positive zero for type inference.
+        static constexpr FP positive_zero{};
+
+        // We want to ensure that +/- 0 is always represented as positive zero
+        if ( lhs == positive_zero ) { lhs = positive_zero; }
+        if ( rhs == positive_zero ) { rhs = positive_zero; }
+
+        // If arguments have different signs, we can handle them by summing
+        // how far are they from 0 each.
+        if ( std::signbit( lhs ) != std::signbit( rhs ) ) {
+            return ulpDistance( std::abs( lhs ), positive_zero ) +
+                   ulpDistance( std::abs( rhs ), positive_zero );
+        }
+
+        // When both lhs and rhs are of the same sign, we can just
+        // read the numbers bitwise as integers, and then subtract them
+        // (assuming IEEE).
+        uint64_t lc = Detail::convertToBits( lhs );
+        uint64_t rc = Detail::convertToBits( rhs );
+
+        // The ulp distance between two numbers is symmetric, so to avoid
+        // dealing with overflows we want the bigger converted number on the lhs
+        if ( lc < rc ) {
+            std::swap( lc, rc );
+        }
+
+        return lc - rc;
+    }
+
+} // end namespace Catch
+
+#endif // CATCH_FLOATING_POINT_HELPERS_HPP_INCLUDED

src/catch2/matchers/catch_matchers_floating_point.cpp

@@ -10,12 +10,12 @@
 #include <catch2/internal/catch_polyfills.hpp>
 #include <catch2/internal/catch_to_string.hpp>
 #include <catch2/catch_tostring.hpp>
+#include <catch2/internal/catch_floating_point_helpers.hpp>
 
 #include <algorithm>
 #include <cmath>
 #include <cstdlib>
 #include <cstdint>
-#include <cstring>
 #include <sstream>
 #include <iomanip>
 #include <limits>
@@ -24,20 +24,6 @@
 namespace Catch {
 namespace {
 
-    int32_t convert(float f) {
-        static_assert(sizeof(float) == sizeof(int32_t), "Important ULP matcher assumption violated");
-        int32_t i;
-        std::memcpy(&i, &f, sizeof(f));
-        return i;
-    }
-
-    int64_t convert(double d) {
-        static_assert(sizeof(double) == sizeof(int64_t), "Important ULP matcher assumption violated");
-        int64_t i;
-        std::memcpy(&i, &d, sizeof(d));
-        return i;
-    }
-
     template <typename FP>
     bool almostEqualUlps(FP lhs, FP rhs, uint64_t maxUlpDiff) {
         // Comparison with NaN should always be false.
@@ -46,17 +32,10 @@ namespace {
             return false;
         }
 
-        auto lc = convert(lhs);
-        auto rc = convert(rhs);
+        // This should also handle positive and negative zeros, infinities
+        const auto ulpDist = ulpDistance(lhs, rhs);
 
-        if ((lc < 0) != (rc < 0)) {
-            // Potentially we can have +0 and -0
-            return lhs == rhs;
-        }
-
-        // static cast as a workaround for IBM XLC
-        auto ulpDiff = std::abs(static_cast<FP>(lc - rc));
-        return static_cast<uint64_t>(ulpDiff) <= maxUlpDiff;
+        return ulpDist <= maxUlpDiff;
     }
 
 #if defined(CATCH_CONFIG_GLOBAL_NEXTAFTER)
@@ -131,6 +110,9 @@ namespace Detail {
         CATCH_ENFORCE(m_type == Detail::FloatingPointKind::Double
                    || m_ulps < (std::numeric_limits<uint32_t>::max)(),
             "Provided ULP is impossibly large for a float comparison.");
+        CATCH_ENFORCE( std::numeric_limits<double>::is_iec559,
+                       "WithinUlp matcher only supports platforms with "
+                       "IEEE-754 compatible floating point representation" );
     }
 
 #if defined(__clang__)