17
17
#include < __type_traits/is_arithmetic.h>
18
18
#include < __type_traits/is_same.h>
19
19
#include < __type_traits/promote.h>
20
- #include < array >
20
+ #include < __utility/pair.h >
21
21
#include < limits>
22
22
23
23
#if !defined(_LIBCPP_HAS_NO_PRAGMA_SYSTEM_HEADER)
@@ -48,42 +48,29 @@ inline _LIBCPP_HIDE_FROM_ABI typename __promote<_A1, _A2>::type hypot(_A1 __x, _
48
48
49
49
#if _LIBCPP_STD_VER >= 17
50
50
// Factors needed to determine if over-/underflow might happen for `std::hypot(x,y,z)`.
51
+ // returns [overflow_threshold, overflow_scale]
51
52
template <class _Real >
52
- struct __hypot_factors {
53
- _Real __threshold;
54
- _Real __scale_xyz;
55
- _Real __scale_M;
56
- };
57
-
58
- // Computes `__hypot_factors` needed to determine if over-/underflow might happen for `std::hypot(x,y,z)`.
59
- // Returns: [underflow_factors, overflow_factors]
60
- template <class _Real >
61
- _LIBCPP_HIDE_FROM_ABI array<__math::__hypot_factors<_Real>, 2 > __create_factors () {
53
+ _LIBCPP_HIDE_FROM_ABI std::pair<_Real, _Real> __hypot_factors () {
62
54
static_assert (std::numeric_limits<_Real>::is_iec559);
63
55
64
- __math::__hypot_factors<_Real> __underflow, __overflow;
65
56
if constexpr (std::is_same_v<_Real, float >) {
66
57
static_assert (-125 == std::numeric_limits<_Real>::min_exponent);
67
58
static_assert (+128 == std::numeric_limits<_Real>::max_exponent);
68
- __underflow = __math::__hypot_factors<_Real>{0x1 .0p-62f , 0x1 .0p70f, 0x1 .0p-70f };
69
- __overflow = __math::__hypot_factors<_Real>{0x1 .0p62f, 0x1 .0p-70f , 0x1 .0p+70f };
59
+ return {0x1 .0p+62f , 0x1 .0p-70f };
70
60
} else if constexpr (std::is_same_v<_Real, double >) {
71
61
static_assert (-1021 == std::numeric_limits<_Real>::min_exponent);
72
62
static_assert (+1024 == std::numeric_limits<_Real>::max_exponent);
73
- __underflow = __math::__hypot_factors<_Real>{0x1 .0p-510 , 0x1 .0p600, 0x1 .0p-600 };
74
- __overflow = __math::__hypot_factors<_Real>{0x1 .0p510, 0x1 .0p-600 , 0x1 .0p+600 };
63
+ return {0x1 .0p+510 , 0x1 .0p-600 };
75
64
} else { // long double
76
65
static_assert (std::is_same_v<_Real, long double >);
77
66
if constexpr (sizeof (_Real) == sizeof (double ))
78
- return static_cast <std::array<__math::__hypot_factors< _Real>, 2 >>(__math::__create_factors <double >());
67
+ return static_cast <std::pair< _Real, _Real >>(__math::__hypot_factors <double >());
79
68
else {
80
69
static_assert (-16'381 == std::numeric_limits<_Real>::min_exponent);
81
70
static_assert (+16'384 == std::numeric_limits<_Real>::max_exponent);
82
- __underflow = __math::__hypot_factors<_Real>{0x1 .0p-8'190l , 0x1 .0p9' 000l, 0x1.0p-9' 000l };
83
- __overflow = __math::__hypot_factors<_Real>{0x1 .0p8' 190l, 0x1.0p-9' 000l , 0x1 .0p+9'000l };
71
+ return {0x1 .0p+8'190l , 0x1 .0p-9'000l };
84
72
}
85
73
}
86
- return {__underflow, __overflow};
87
74
}
88
75
89
76
// Computes the three-dimensional hypotenuse: `std::hypot(x,y,z)`.
@@ -92,21 +79,22 @@ _LIBCPP_HIDE_FROM_ABI array<__math::__hypot_factors<_Real>, 2> __create_factors(
92
79
// See https://github.com/llvm/llvm-project/issues/92782 for a detailed discussion and summary.
93
80
template <class _Real >
94
81
_LIBCPP_HIDE_FROM_ABI _Real __hypot (_Real __x, _Real __y, _Real __z) {
95
- const auto [__underflow, __overflow] = __math::__create_factors<_Real>();
96
- _Real __M = std::max ({__math::fabs (__x), __math::fabs (__y), __math::fabs (__z)});
97
- if (__M > __overflow.__threshold ) { // x*x + y*y + z*z might overflow
98
- __x *= __overflow.__scale_xyz ;
99
- __y *= __overflow.__scale_xyz ;
100
- __z *= __overflow.__scale_xyz ;
101
- __M = __overflow.__scale_M ;
102
- } else if (__M < __underflow.__threshold ) { // x*x + y*y + z*z might underflow
103
- __x *= __underflow.__scale_xyz ;
104
- __y *= __underflow.__scale_xyz ;
105
- __z *= __underflow.__scale_xyz ;
106
- __M = __underflow.__scale_M ;
82
+ const _Real __max_abs = std::max ({__math::fabs (__x), __math::fabs (__y), __math::fabs (__z)});
83
+ const auto [__overflow_threshold, __overflow_scale] = __math::__hypot_factors<_Real>();
84
+ _Real __scale;
85
+ if (__max_abs > __overflow_threshold) { // x*x + y*y + z*z might overflow
86
+ __scale = __overflow_scale;
87
+ __x *= __scale;
88
+ __y *= __scale;
89
+ __z *= __scale;
90
+ } else if (__max_abs < 1 / __overflow_threshold) { // x*x + y*y + z*z might underflow
91
+ __scale = 1 / __overflow_scale;
92
+ __x *= __scale;
93
+ __y *= __scale;
94
+ __z *= __scale;
107
95
} else
108
- __M = 1 ;
109
- return __M * __math::sqrt (__x * __x + __y * __y + __z * __z);
96
+ __scale = 1 ;
97
+ return __math::sqrt (__x * __x + __y * __y + __z * __z) / __scale ;
110
98
}
111
99
112
100
inline _LIBCPP_HIDE_FROM_ABI float hypot (float __x, float __y, float __z) { return __math::__hypot (__x, __y, __z); }
@@ -125,7 +113,8 @@ _LIBCPP_HIDE_FROM_ABI typename __promote<_A1, _A2, _A3>::type hypot(_A1 __x, _A2
125
113
using __result_type = typename __promote<_A1, _A2, _A3>::type;
126
114
static_assert (!(
127
115
std::is_same_v<_A1, __result_type> && std::is_same_v<_A2, __result_type> && std::is_same_v<_A3, __result_type>));
128
- return __math::__hypot (static_cast <__result_type>(__x), static_cast <__result_type>(__y), static_cast <__result_type>(__z));
116
+ return __math::__hypot (
117
+ static_cast <__result_type>(__x), static_cast <__result_type>(__y), static_cast <__result_type>(__z));
129
118
}
130
119
#endif
131
120
0 commit comments