[libcxx-commits] [libcxx] [libc++][math] Fix undue overflowing of `std::hypot(x, y, z)` (PR #93350)

via libcxx-commits libcxx-commits at lists.llvm.org
Thu Jul 4 11:55:20 PDT 2024


https://github.com/PaulXiCao updated https://github.com/llvm/llvm-project/pull/93350

>From c650ee99b7f35e8459687cac069887b581d79802 Mon Sep 17 00:00:00 2001
From: Paul <{ID}+{username}@users.noreply.github.com>
Date: Fri, 24 May 2024 15:15:25 +0200
Subject: [PATCH 1/8] make use of C++17 is_same_v

---
 .../test/std/numerics/c.math/cmath.pass.cpp   | 30 +++++++++----------
 1 file changed, 15 insertions(+), 15 deletions(-)

diff --git a/libcxx/test/std/numerics/c.math/cmath.pass.cpp b/libcxx/test/std/numerics/c.math/cmath.pass.cpp
index 9379084499792..ad4b46b3713fc 100644
--- a/libcxx/test/std/numerics/c.math/cmath.pass.cpp
+++ b/libcxx/test/std/numerics/c.math/cmath.pass.cpp
@@ -1136,21 +1136,21 @@ void test_hypot()
     assert(std::hypot(3,4) == 5);
 
 #if TEST_STD_VER > 14
-    static_assert((std::is_same<decltype(std::hypot((float)0, (float)0, (float)0)), float>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (bool)0, (float)0)), double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (unsigned short)0, (double)0)), double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (int)0, (long double)0)), long double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (double)0, (long)0)), double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (long double)0, (unsigned long)0)), long double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (int)0, (long long)0)), double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (int)0, (unsigned long long)0)), double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (double)0, (double)0)), double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (long double)0, (long double)0)), long double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (float)0, (double)0)), double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (float)0, (long double)0)), long double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((float)0, (double)0, (long double)0)), long double>::value), "");
-    static_assert((std::is_same<decltype(std::hypot((int)0, (int)0, (int)0)), double>::value), "");
-    static_assert((std::is_same<decltype(hypot(Ambiguous(), Ambiguous(), Ambiguous())), Ambiguous>::value), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0, (float)0)), float>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (bool)0, (float)0)), double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (unsigned short)0, (double)0)), double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0, (long double)0)), long double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0, (long)0)), double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (long double)0, (unsigned long)0)), long double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0, (long long)0)), double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0, (unsigned long long)0)), double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0, (double)0)), double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (long double)0, (long double)0)), long double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0, (double)0)), double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0, (long double)0)), long double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0, (long double)0)), long double>), "");
+    static_assert((std::is_same_v<decltype(std::hypot((int)0, (int)0, (int)0)), double>), "");
+    static_assert((std::is_same_v<decltype(hypot(Ambiguous(), Ambiguous(), Ambiguous())), Ambiguous>), "");
 
     assert(std::hypot(2,3,6) == 7);
     assert(std::hypot(1,4,8) == 9);

>From ba77ceae6d2c8bb53e071c6cbfc734cf0724e6a2 Mon Sep 17 00:00:00 2001
From: Paul <{ID}+{username}@users.noreply.github.com>
Date: Fri, 24 May 2024 15:20:28 +0200
Subject: [PATCH 2/8] style change: align argument types

---
 .../test/std/numerics/c.math/cmath.pass.cpp   | 32 ++++++++++---------
 1 file changed, 17 insertions(+), 15 deletions(-)

diff --git a/libcxx/test/std/numerics/c.math/cmath.pass.cpp b/libcxx/test/std/numerics/c.math/cmath.pass.cpp
index ad4b46b3713fc..a16aa0724dffc 100644
--- a/libcxx/test/std/numerics/c.math/cmath.pass.cpp
+++ b/libcxx/test/std/numerics/c.math/cmath.pass.cpp
@@ -1136,21 +1136,23 @@ void test_hypot()
     assert(std::hypot(3,4) == 5);
 
 #if TEST_STD_VER > 14
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0, (float)0)), float>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (bool)0, (float)0)), double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (unsigned short)0, (double)0)), double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0, (long double)0)), long double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0, (long)0)), double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (long double)0, (unsigned long)0)), long double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0, (long long)0)), double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0, (unsigned long long)0)), double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0, (double)0)), double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (long double)0, (long double)0)), long double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0, (double)0)), double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0, (long double)0)), long double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0, (long double)0)), long double>), "");
-    static_assert((std::is_same_v<decltype(std::hypot((int)0, (int)0, (int)0)), double>), "");
-    static_assert((std::is_same_v<decltype(hypot(Ambiguous(), Ambiguous(), Ambiguous())), Ambiguous>), "");
+    // clang-format off
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0,          (float)0)),              float>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (bool)0,           (float)0)),              double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (unsigned short)0, (double)0)),             double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0,            (long double)0)),        long double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0,         (long)0)),               double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (long double)0,    (unsigned long)0)),      long double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0,            (long long)0)),          double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (int)0,            (unsigned long long)0)), double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0,         (double)0)),             double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (long double)0,    (long double)0)),        long double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0,          (double)0)),             double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (float)0,          (long double)0)),        long double>));
+    static_assert((std::is_same_v<decltype(std::hypot((float)0, (double)0,         (long double)0)),        long double>));
+    static_assert((std::is_same_v<decltype(std::hypot((int)0,   (int)0,            (int)0)),                double>));
+    static_assert((std::is_same_v<decltype(hypot(Ambiguous(), Ambiguous(), Ambiguous())), Ambiguous>));
+    // clang-format on
 
     assert(std::hypot(2,3,6) == 7);
     assert(std::hypot(1,4,8) == 9);

>From 4f2f6b32f5851ba3a059725b0c941c6775dab1ac Mon Sep 17 00:00:00 2001
From: Paul <{ID}+{username}@users.noreply.github.com>
Date: Fri, 24 May 2024 23:18:03 +0200
Subject: [PATCH 3/8] implementation+testing

---
 libcxx/include/__math/hypot.h                 | 78 +++++++++++++++++++
 libcxx/include/cmath                          | 25 +-----
 .../test/std/numerics/c.math/cmath.pass.cpp   | 45 +++++++++++
 3 files changed, 124 insertions(+), 24 deletions(-)

diff --git a/libcxx/include/__math/hypot.h b/libcxx/include/__math/hypot.h
index 8e2c0f5b6b6d8..6b68e7bfed341 100644
--- a/libcxx/include/__math/hypot.h
+++ b/libcxx/include/__math/hypot.h
@@ -9,11 +9,16 @@
 #ifndef _LIBCPP___MATH_HYPOT_H
 #define _LIBCPP___MATH_HYPOT_H
 
+#include <__algorithm/max.h>
 #include <__config>
+#include <__math/abs.h>
+#include <__math/roots.h>
 #include <__type_traits/enable_if.h>
 #include <__type_traits/is_arithmetic.h>
 #include <__type_traits/is_same.h>
 #include <__type_traits/promote.h>
+#include <array>
+#include <limits>
 
 #if !defined(_LIBCPP_HAS_NO_PRAGMA_SYSTEM_HEADER)
 #  pragma GCC system_header
@@ -41,6 +46,79 @@ inline _LIBCPP_HIDE_FROM_ABI typename __promote<_A1, _A2>::type hypot(_A1 __x, _
   return __math::hypot((__result_type)__x, (__result_type)__y);
 }
 
+#if _LIBCPP_STD_VER >= 17
+template <class _Real>
+struct __hypot_factors {
+  _Real __threshold;
+  _Real __scale_xyz;
+  _Real __scale_M;
+};
+
+// returns [underflow_factors, overflow_factors]
+template <class _Real>
+std::array<__hypot_factors<_Real>, 2> __create_factors() {
+  static_assert(std::numeric_limits<_Real>::is_iec559);
+
+  __hypot_factors<_Real> __underflow, __overflow;
+  if constexpr (std::is_same_v<_Real, float>) {
+    static_assert(-125 == std::numeric_limits<_Real>::min_exponent);
+    static_assert(+128 == std::numeric_limits<_Real>::max_exponent);
+    __underflow = __hypot_factors<_Real>{0x1.0p-62f, 0x1.0p70f, 0x1.0p-70f};
+    __overflow  = __hypot_factors<_Real>{0x1.0p62f, 0x1.0p-70f, 0x1.0p+70f};
+  } else if constexpr (std::is_same_v<_Real, double>) {
+    static_assert(-1021 == std::numeric_limits<_Real>::min_exponent);
+    static_assert(+1024 == std::numeric_limits<_Real>::max_exponent);
+    __underflow = __hypot_factors<_Real>{0x1.0p-510, 0x1.0p600, 0x1.0p-600};
+    __overflow  = __hypot_factors<_Real>{0x1.0p510, 0x1.0p-600, 0x1.0p+600};
+  } else { // long double
+    static_assert(std::is_same_v<_Real, long double>);
+    static_assert(-16'381 == std::numeric_limits<_Real>::min_exponent);
+    static_assert(+16'384 == std::numeric_limits<_Real>::max_exponent);
+    __underflow = __hypot_factors<_Real>{0x1.0p-8'190l, 0x1.0p9'000l, 0x1.0p-9'000l};
+    __overflow  = __hypot_factors<_Real>{0x1.0p8'190l, 0x1.0p-9'000l, 0x1.0p+9'000l};
+  }
+  return {__underflow, __overflow};
+}
+
+template <class _Real>
+_Real __hypot(_Real __x, _Real __y, _Real __z) {
+  const auto [__underflow, __overflow] = __create_factors<_Real>();
+  _Real __M                            = std::max({__math::fabs(__x), __math::fabs(__y), __math::fabs(__z)});
+  if (__M > __overflow.__threshold) { // x*x + y*y + z*z might overflow
+    __x *= __overflow.__scale_xyz;
+    __y *= __overflow.__scale_xyz;
+    __z *= __overflow.__scale_xyz;
+    __M = __overflow.__scale_M;
+  } else if (__M < __underflow.__threshold) { // x*x + y*y + z*z might underflow
+    __x *= __underflow.__scale_xyz;
+    __y *= __underflow.__scale_xyz;
+    __z *= __underflow.__scale_xyz;
+    __M = __underflow.__scale_M;
+  } else
+    __M = 1;
+  return __M * __math::sqrt(__x * __x + __y * __y + __z * __z);
+}
+
+inline _LIBCPP_HIDE_FROM_ABI float hypot(float __x, float __y, float __z) { return __hypot(__x, __y, __z); }
+
+inline _LIBCPP_HIDE_FROM_ABI double hypot(double __x, double __y, double __z) { return __hypot(__x, __y, __z); }
+
+inline _LIBCPP_HIDE_FROM_ABI long double hypot(long double __x, long double __y, long double __z) {
+  return __hypot(__x, __y, __z);
+}
+
+template <class _A1,
+          class _A2,
+          class _A3,
+          std::enable_if_t< is_arithmetic_v<_A1> && is_arithmetic_v<_A2> && is_arithmetic_v<_A3>, int> = 0 >
+_LIBCPP_HIDE_FROM_ABI typename __promote<_A1, _A2, _A3>::type hypot(_A1 __x, _A2 __y, _A3 __z) _NOEXCEPT {
+  using __result_type = typename __promote<_A1, _A2, _A3>::type;
+  static_assert(!(
+      std::is_same_v<_A1, __result_type> && std::is_same_v<_A2, __result_type> && std::is_same_v<_A3, __result_type>));
+  return __hypot(static_cast<__result_type>(__x), static_cast<__result_type>(__y), static_cast<__result_type>(__z));
+}
+#endif
+
 } // namespace __math
 
 _LIBCPP_END_NAMESPACE_STD
diff --git a/libcxx/include/cmath b/libcxx/include/cmath
index dd194bbb55896..49f35caa7b74e 100644
--- a/libcxx/include/cmath
+++ b/libcxx/include/cmath
@@ -305,6 +305,7 @@ constexpr long double lerp(long double a, long double b, long double t) noexcept
 */
 
 #include <__config>
+#include <__math/hypot.h>
 #include <__type_traits/enable_if.h>
 #include <__type_traits/is_arithmetic.h>
 #include <__type_traits/is_constant_evaluated.h>
@@ -544,30 +545,6 @@ using ::scalbnl _LIBCPP_USING_IF_EXISTS;
 using ::tgammal _LIBCPP_USING_IF_EXISTS;
 using ::truncl _LIBCPP_USING_IF_EXISTS;
 
-#if _LIBCPP_STD_VER >= 17
-inline _LIBCPP_HIDE_FROM_ABI float hypot(float __x, float __y, float __z) {
-  return sqrt(__x * __x + __y * __y + __z * __z);
-}
-inline _LIBCPP_HIDE_FROM_ABI double hypot(double __x, double __y, double __z) {
-  return sqrt(__x * __x + __y * __y + __z * __z);
-}
-inline _LIBCPP_HIDE_FROM_ABI long double hypot(long double __x, long double __y, long double __z) {
-  return sqrt(__x * __x + __y * __y + __z * __z);
-}
-
-template <class _A1, class _A2, class _A3>
-inline _LIBCPP_HIDE_FROM_ABI
-    typename enable_if_t< is_arithmetic<_A1>::value && is_arithmetic<_A2>::value && is_arithmetic<_A3>::value,
-                          __promote<_A1, _A2, _A3> >::type
-    hypot(_A1 __lcpp_x, _A2 __lcpp_y, _A3 __lcpp_z) _NOEXCEPT {
-  typedef typename __promote<_A1, _A2, _A3>::type __result_type;
-  static_assert((!(is_same<_A1, __result_type>::value && is_same<_A2, __result_type>::value &&
-                   is_same<_A3, __result_type>::value)),
-                "");
-  return std::hypot((__result_type)__lcpp_x, (__result_type)__lcpp_y, (__result_type)__lcpp_z);
-}
-#endif
-
 template <class _A1, __enable_if_t<is_floating_point<_A1>::value, int> = 0>
 _LIBCPP_HIDE_FROM_ABI _LIBCPP_CONSTEXPR bool __constexpr_isnan(_A1 __lcpp_x) _NOEXCEPT {
 #if __has_builtin(__builtin_isnan)
diff --git a/libcxx/test/std/numerics/c.math/cmath.pass.cpp b/libcxx/test/std/numerics/c.math/cmath.pass.cpp
index a16aa0724dffc..544c197ceb1af 100644
--- a/libcxx/test/std/numerics/c.math/cmath.pass.cpp
+++ b/libcxx/test/std/numerics/c.math/cmath.pass.cpp
@@ -12,14 +12,17 @@
 
 // <cmath>
 
+#include <array>
 #include <cmath>
 #include <limits>
 #include <type_traits>
 #include <cassert>
 
+#include "fp_compare.h"
 #include "test_macros.h"
 #include "hexfloat.h"
 #include "truncate_fp.h"
+#include "type_algorithms.h"
 
 // convertible to int/float/double/etc
 template <class T, int N=0>
@@ -1113,6 +1116,44 @@ void test_fmin()
     assert(std::fmin(1,0) == 0);
 }
 
+struct TestHypot3 {
+    template <class Real>
+    static void operator()() {
+        const auto check = [](Real elem, Real abs_tol) {
+            assert(std::isfinite(std::hypot(elem, Real(0), Real(0))));
+            assert(fptest_close(std::hypot(elem, Real(0), Real(0)), elem, abs_tol));
+            assert(std::isfinite(std::hypot(elem, elem, Real(0))));
+            assert(fptest_close(std::hypot(elem, elem, Real(0)), std::sqrt(Real(2)) * elem, abs_tol));
+            assert(std::isfinite(std::hypot(elem, elem, elem)));
+            assert(fptest_close(std::hypot(elem, elem, elem), std::sqrt(Real(3)) * elem, abs_tol));
+        };
+
+        { // check for overflow
+            const auto [elem, abs_tol] = []() -> std::array<Real, 2> {
+                if constexpr (std::is_same_v<Real, float>)
+                    return {1e20f, 1e16f};
+                else if constexpr (std::is_same_v<Real, double>)
+                    return {1e300, 1e287};
+                else // long double
+                    return {1e4000l, 1e3985l};
+            }();
+            check(elem, abs_tol);
+        }
+
+        { // check for underflow
+            const auto [elem, abs_tol] = []() -> std::array<Real, 2> {
+                if constexpr (std::is_same_v<Real, float>)
+                    return {1e-20f, 1e-24f};
+                else if constexpr (std::is_same_v<Real, double>)
+                    return {1e-287, 1e-300};
+                else // long double
+                    return {1e-3985l, 1e-4000l};
+            }();
+            check(elem, abs_tol);
+        }
+    }
+};
+
 void test_hypot()
 {
     static_assert((std::is_same<decltype(std::hypot((float)0, (float)0)), float>::value), "");
@@ -1156,6 +1197,10 @@ void test_hypot()
 
     assert(std::hypot(2,3,6) == 7);
     assert(std::hypot(1,4,8) == 9);
+
+    // Check for undue over-/underflows of intermediate results.
+    // See discussion at https://github.com/llvm/llvm-project/issues/92782.
+    types::for_each(types::floating_point_types(), TestHypot3());
 #endif
 }
 

>From bc7ef4868665c13f4b5ec14f165a9a3099b2fa6c Mon Sep 17 00:00:00 2001
From: Paul <{ID}+{username}@users.noreply.github.com>
Date: Sat, 25 May 2024 13:08:36 +0200
Subject: [PATCH 4/8] handle case: sizeof(double) == sizeof(long double)

---
 libcxx/include/__math/hypot.h | 12 ++++++++----
 1 file changed, 8 insertions(+), 4 deletions(-)

diff --git a/libcxx/include/__math/hypot.h b/libcxx/include/__math/hypot.h
index 6b68e7bfed341..835bc20f6d1ab 100644
--- a/libcxx/include/__math/hypot.h
+++ b/libcxx/include/__math/hypot.h
@@ -72,10 +72,14 @@ std::array<__hypot_factors<_Real>, 2> __create_factors() {
     __overflow  = __hypot_factors<_Real>{0x1.0p510, 0x1.0p-600, 0x1.0p+600};
   } else { // long double
     static_assert(std::is_same_v<_Real, long double>);
-    static_assert(-16'381 == std::numeric_limits<_Real>::min_exponent);
-    static_assert(+16'384 == std::numeric_limits<_Real>::max_exponent);
-    __underflow = __hypot_factors<_Real>{0x1.0p-8'190l, 0x1.0p9'000l, 0x1.0p-9'000l};
-    __overflow  = __hypot_factors<_Real>{0x1.0p8'190l, 0x1.0p-9'000l, 0x1.0p+9'000l};
+    if constexpr (sizeof(_Real) == sizeof(double))
+      return static_cast<std::array<__hypot_factors<_Real>, 2>>(__create_factors<double>());
+    else {
+      static_assert(-16'381 == std::numeric_limits<_Real>::min_exponent);
+      static_assert(+16'384 == std::numeric_limits<_Real>::max_exponent);
+      __underflow = __hypot_factors<_Real>{0x1.0p-8'190l, 0x1.0p9'000l, 0x1.0p-9'000l};
+      __overflow  = __hypot_factors<_Real>{0x1.0p8'190l, 0x1.0p-9'000l, 0x1.0p+9'000l};
+    }
   }
   return {__underflow, __overflow};
 }

>From ada9edbb84e415b443fb23bc792effcba873762f Mon Sep 17 00:00:00 2001
From: Paul <{ID}+{username}@users.noreply.github.com>
Date: Thu, 4 Jul 2024 20:05:58 +0200
Subject: [PATCH 5/8] add hide_from_abi macro to internal functions

---
 libcxx/include/__math/hypot.h | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/libcxx/include/__math/hypot.h b/libcxx/include/__math/hypot.h
index 835bc20f6d1ab..e8965b957fa31 100644
--- a/libcxx/include/__math/hypot.h
+++ b/libcxx/include/__math/hypot.h
@@ -56,7 +56,7 @@ struct __hypot_factors {
 
 // returns [underflow_factors, overflow_factors]
 template <class _Real>
-std::array<__hypot_factors<_Real>, 2> __create_factors() {
+_LIBCPP_HIDE_FROM_ABI std::array<__hypot_factors<_Real>, 2> __create_factors() {
   static_assert(std::numeric_limits<_Real>::is_iec559);
 
   __hypot_factors<_Real> __underflow, __overflow;
@@ -85,7 +85,7 @@ std::array<__hypot_factors<_Real>, 2> __create_factors() {
 }
 
 template <class _Real>
-_Real __hypot(_Real __x, _Real __y, _Real __z) {
+_LIBCPP_HIDE_FROM_ABI _Real __hypot(_Real __x, _Real __y, _Real __z) {
   const auto [__underflow, __overflow] = __create_factors<_Real>();
   _Real __M                            = std::max({__math::fabs(__x), __math::fabs(__y), __math::fabs(__z)});
   if (__M > __overflow.__threshold) { // x*x + y*y + z*z might overflow

>From ad123177c71e3992d238a33f8935a9e29a0f7222 Mon Sep 17 00:00:00 2001
From: Paul <{ID}+{username}@users.noreply.github.com>
Date: Thu, 4 Jul 2024 20:14:00 +0200
Subject: [PATCH 6/8] fully qualify local function calls

---
 libcxx/include/__math/hypot.h | 28 ++++++++++++++--------------
 1 file changed, 14 insertions(+), 14 deletions(-)

diff --git a/libcxx/include/__math/hypot.h b/libcxx/include/__math/hypot.h
index e8965b957fa31..397f77c6c4387 100644
--- a/libcxx/include/__math/hypot.h
+++ b/libcxx/include/__math/hypot.h
@@ -56,29 +56,29 @@ struct __hypot_factors {
 
 // returns [underflow_factors, overflow_factors]
 template <class _Real>
-_LIBCPP_HIDE_FROM_ABI std::array<__hypot_factors<_Real>, 2> __create_factors() {
+_LIBCPP_HIDE_FROM_ABI std::array<__math::__hypot_factors<_Real>, 2> __create_factors() {
   static_assert(std::numeric_limits<_Real>::is_iec559);
 
-  __hypot_factors<_Real> __underflow, __overflow;
+  __math::__hypot_factors<_Real> __underflow, __overflow;
   if constexpr (std::is_same_v<_Real, float>) {
     static_assert(-125 == std::numeric_limits<_Real>::min_exponent);
     static_assert(+128 == std::numeric_limits<_Real>::max_exponent);
-    __underflow = __hypot_factors<_Real>{0x1.0p-62f, 0x1.0p70f, 0x1.0p-70f};
-    __overflow  = __hypot_factors<_Real>{0x1.0p62f, 0x1.0p-70f, 0x1.0p+70f};
+    __underflow = __math::__hypot_factors<_Real>{0x1.0p-62f, 0x1.0p70f, 0x1.0p-70f};
+    __overflow  = __math::__hypot_factors<_Real>{0x1.0p62f, 0x1.0p-70f, 0x1.0p+70f};
   } else if constexpr (std::is_same_v<_Real, double>) {
     static_assert(-1021 == std::numeric_limits<_Real>::min_exponent);
     static_assert(+1024 == std::numeric_limits<_Real>::max_exponent);
-    __underflow = __hypot_factors<_Real>{0x1.0p-510, 0x1.0p600, 0x1.0p-600};
-    __overflow  = __hypot_factors<_Real>{0x1.0p510, 0x1.0p-600, 0x1.0p+600};
+    __underflow = __math::__hypot_factors<_Real>{0x1.0p-510, 0x1.0p600, 0x1.0p-600};
+    __overflow  = __math::__hypot_factors<_Real>{0x1.0p510, 0x1.0p-600, 0x1.0p+600};
   } else { // long double
     static_assert(std::is_same_v<_Real, long double>);
     if constexpr (sizeof(_Real) == sizeof(double))
-      return static_cast<std::array<__hypot_factors<_Real>, 2>>(__create_factors<double>());
+      return static_cast<std::array<__math::__hypot_factors<_Real>, 2>>(__math::__create_factors<double>());
     else {
       static_assert(-16'381 == std::numeric_limits<_Real>::min_exponent);
       static_assert(+16'384 == std::numeric_limits<_Real>::max_exponent);
-      __underflow = __hypot_factors<_Real>{0x1.0p-8'190l, 0x1.0p9'000l, 0x1.0p-9'000l};
-      __overflow  = __hypot_factors<_Real>{0x1.0p8'190l, 0x1.0p-9'000l, 0x1.0p+9'000l};
+      __underflow = __math::__hypot_factors<_Real>{0x1.0p-8'190l, 0x1.0p9'000l, 0x1.0p-9'000l};
+      __overflow  = __math::__hypot_factors<_Real>{0x1.0p8'190l, 0x1.0p-9'000l, 0x1.0p+9'000l};
     }
   }
   return {__underflow, __overflow};
@@ -86,7 +86,7 @@ _LIBCPP_HIDE_FROM_ABI std::array<__hypot_factors<_Real>, 2> __create_factors() {
 
 template <class _Real>
 _LIBCPP_HIDE_FROM_ABI _Real __hypot(_Real __x, _Real __y, _Real __z) {
-  const auto [__underflow, __overflow] = __create_factors<_Real>();
+  const auto [__underflow, __overflow] = __math::__create_factors<_Real>();
   _Real __M                            = std::max({__math::fabs(__x), __math::fabs(__y), __math::fabs(__z)});
   if (__M > __overflow.__threshold) { // x*x + y*y + z*z might overflow
     __x *= __overflow.__scale_xyz;
@@ -103,12 +103,12 @@ _LIBCPP_HIDE_FROM_ABI _Real __hypot(_Real __x, _Real __y, _Real __z) {
   return __M * __math::sqrt(__x * __x + __y * __y + __z * __z);
 }
 
-inline _LIBCPP_HIDE_FROM_ABI float hypot(float __x, float __y, float __z) { return __hypot(__x, __y, __z); }
+inline _LIBCPP_HIDE_FROM_ABI float hypot(float __x, float __y, float __z) { return __math::__hypot(__x, __y, __z); }
 
-inline _LIBCPP_HIDE_FROM_ABI double hypot(double __x, double __y, double __z) { return __hypot(__x, __y, __z); }
+inline _LIBCPP_HIDE_FROM_ABI double hypot(double __x, double __y, double __z) { return __math::__hypot(__x, __y, __z); }
 
 inline _LIBCPP_HIDE_FROM_ABI long double hypot(long double __x, long double __y, long double __z) {
-  return __hypot(__x, __y, __z);
+  return __math::__hypot(__x, __y, __z);
 }
 
 template <class _A1,
@@ -119,7 +119,7 @@ _LIBCPP_HIDE_FROM_ABI typename __promote<_A1, _A2, _A3>::type hypot(_A1 __x, _A2
   using __result_type = typename __promote<_A1, _A2, _A3>::type;
   static_assert(!(
       std::is_same_v<_A1, __result_type> && std::is_same_v<_A2, __result_type> && std::is_same_v<_A3, __result_type>));
-  return __hypot(static_cast<__result_type>(__x), static_cast<__result_type>(__y), static_cast<__result_type>(__z));
+  return __math::__hypot(static_cast<__result_type>(__x), static_cast<__result_type>(__y), static_cast<__result_type>(__z));
 }
 #endif
 

>From a80b11e49e67a816c6c7cb116ea155ad9957764b Mon Sep 17 00:00:00 2001
From: Paul <{ID}+{username}@users.noreply.github.com>
Date: Thu, 4 Jul 2024 20:33:56 +0200
Subject: [PATCH 7/8] remove for operator(): friendlier for older compilers

---
 libcxx/test/std/numerics/c.math/cmath.pass.cpp | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/libcxx/test/std/numerics/c.math/cmath.pass.cpp b/libcxx/test/std/numerics/c.math/cmath.pass.cpp
index 544c197ceb1af..6ed37c272192e 100644
--- a/libcxx/test/std/numerics/c.math/cmath.pass.cpp
+++ b/libcxx/test/std/numerics/c.math/cmath.pass.cpp
@@ -1118,7 +1118,7 @@ void test_fmin()
 
 struct TestHypot3 {
     template <class Real>
-    static void operator()() {
+    void operator()() const {
         const auto check = [](Real elem, Real abs_tol) {
             assert(std::isfinite(std::hypot(elem, Real(0), Real(0))));
             assert(fptest_close(std::hypot(elem, Real(0), Real(0)), elem, abs_tol));

>From d5d3ec524baaa2d07216f3b218352255cff4e07b Mon Sep 17 00:00:00 2001
From: Paul <{ID}+{username}@users.noreply.github.com>
Date: Thu, 4 Jul 2024 20:55:03 +0200
Subject: [PATCH 8/8] document non-naive implementation

---
 libcxx/include/__math/hypot.h | 8 +++++++-
 1 file changed, 7 insertions(+), 1 deletion(-)

diff --git a/libcxx/include/__math/hypot.h b/libcxx/include/__math/hypot.h
index 397f77c6c4387..960d6919d5941 100644
--- a/libcxx/include/__math/hypot.h
+++ b/libcxx/include/__math/hypot.h
@@ -47,6 +47,7 @@ inline _LIBCPP_HIDE_FROM_ABI typename __promote<_A1, _A2>::type hypot(_A1 __x, _
 }
 
 #if _LIBCPP_STD_VER >= 17
+// Factors needed to determine if over-/underflow might happen for `std::hypot(x,y,z)`.
 template <class _Real>
 struct __hypot_factors {
   _Real __threshold;
@@ -54,7 +55,8 @@ struct __hypot_factors {
   _Real __scale_M;
 };
 
-// returns [underflow_factors, overflow_factors]
+// Computes `__hypot_factors` needed to determine if over-/underflow might happen for `std::hypot(x,y,z)`.
+// Returns: [underflow_factors, overflow_factors]
 template <class _Real>
 _LIBCPP_HIDE_FROM_ABI std::array<__math::__hypot_factors<_Real>, 2> __create_factors() {
   static_assert(std::numeric_limits<_Real>::is_iec559);
@@ -84,6 +86,10 @@ _LIBCPP_HIDE_FROM_ABI std::array<__math::__hypot_factors<_Real>, 2> __create_fac
   return {__underflow, __overflow};
 }
 
+// Computes the three-dimensional hypotenuse: `std::hypot(x,y,z)`.
+// The naive implementation might over-/underflow which is why this implementation is more involved:
+//    If the square of an argument might run into issues, we scale the arguments appropriately.
+// See https://github.com/llvm/llvm-project/issues/92782 for a detailed discussion and summary.
 template <class _Real>
 _LIBCPP_HIDE_FROM_ABI _Real __hypot(_Real __x, _Real __y, _Real __z) {
   const auto [__underflow, __overflow] = __math::__create_factors<_Real>();



More information about the libcxx-commits mailing list