boost/math
2
0
Fork 0
mirror of https://github.com/boostorg/math.git synced 2026-01-19 04:22:09 +00:00
Code Issues Projects Releases Wiki Activity

Merge pull request #1318 from boostorg/1114

Browse Source 
Add three arg `hypot` to special functions
This commit is contained in:
Matt Borland
2025-09-09 11:37:35 +02:00
committed by GitHub
parent ab265d474f 8c2c5beeaf
commit e97efb104f
7 changed files with 165 additions and 4 deletions
Download Patch File Download Diff File Expand all files Collapse all files

13
doc/sf/powers.qbk
View File

@@ -310,11 +310,17 @@ calculated using NTL::RR at 1000-bit precision.
template <class T1, class T2, class ``__Policy``>
``__sf_result`` hypot(T1 x, T2 y, const ``__Policy``&);
template <class T1, class T2, class T3>
``__sf_result`` hypot(T1 x, T2 y, T3 z);
template <class T1, class T2, class T3, class ``__Policy``>
``__sf_result`` hypot(T1 x, T2 y, T3 z, const ``__Policy``&);
__effects computes [equation hypot]
in such a way as to avoid undue underflow and overflow.
The return type of this function is computed using the __arg_promotion_rules
when T1 and T2 are of different types.
when T1 and T2 (and or T3) are of different types.
[optional_policy]
@@ -328,6 +334,11 @@ The function is even and symmetric in /x/ and /y/, so first take assume
Then if ['x * [epsilon] >= y] we can simply return /x/.
In the 3-arg case the approach is to scale by the largest argument, and then calculate
const auto a = max(max(x, y), z);
return a == 0 ? 0 : a * sqrt(x/a * x/a + y/a * y/a + z/a * z/a);
Otherwise the result is given by:
[equation hypot2]

5
include/boost/math/policies/policy.hpp
View File

@@ -942,7 +942,7 @@ template <typename P>
class is_policy_imp
{
public:
static constexpr bool value = (sizeof(::boost::math::policies::detail::test_is_policy(static_cast<P*>(nullptr))) == sizeof(char));
static constexpr bool value = (sizeof(detail::test_is_policy(static_cast<boost::math::remove_reference_t<boost::math::remove_cv_t<P>>*>(nullptr))) == sizeof(char));
};
}
@@ -955,6 +955,9 @@ public:
using type = boost::math::integral_constant<bool, value>;
};
template <typename T>
BOOST_MATH_INLINE_CONSTEXPR bool is_policy_v = is_policy<T>::value;
//
// Helper traits class for distribution error handling:
//

68
include/boost/math/special_functions/hypot.hpp
View File

@@ -16,6 +16,8 @@
#include <boost/math/tools/type_traits.hpp>
#include <boost/math/policies/error_handling.hpp>
#include <boost/math/special_functions/math_fwd.hpp>
#include <boost/math/tools/utility.hpp>
#include <boost/math/tools/numeric_limits.hpp>
namespace boost{ namespace math{ namespace detail{
@@ -53,6 +55,48 @@ BOOST_MATH_GPU_ENABLED T hypot_imp(T x, T y, const Policy& pol)
return x * sqrt(1 + rat*rat);
} // template <class T> T hypot(T x, T y)
template <class T, class Policy>
BOOST_MATH_GPU_ENABLED T hypot_imp(T x, T y, T z, const Policy& pol)
{
BOOST_MATH_STD_USING
x = fabs(x);
y = fabs(y);
z = fabs(z);
#ifdef _MSC_VER
#pragma warning(push)
#pragma warning(disable: 4127)
#endif
// special case, see C99 Annex F:
BOOST_MATH_IF_CONSTEXPR (boost::math::numeric_limits<T>::has_infinity)
{
if(((x == boost::math::numeric_limits<T>::infinity())
|| (y == boost::math::numeric_limits<T>::infinity())
|| (z == boost::math::numeric_limits<T>::infinity())))
return policies::raise_overflow_error<T>("boost::math::hypot<%1%>(%1%,%1%,%1%)", nullptr, pol);
}
#ifdef _MSC_VER
#pragma warning(pop)
#endif
const T a {(max)((max)(x, y), z)};
if (a == T(0))
{
return a;
}
const T x_div_a {x / a};
const T y_div_a {y / a};
const T z_div_a {z / a};
return a * sqrt(x_div_a * x_div_a
+ y_div_a * y_div_a
+ z_div_a * z_div_a);
}
}
template <class T1, class T2>
@@ -64,7 +108,7 @@ BOOST_MATH_GPU_ENABLED inline typename tools::promote_args<T1, T2>::type
static_cast<result_type>(x), static_cast<result_type>(y), policies::policy<>());
}
template <class T1, class T2, class Policy>
template <class T1, class T2, class Policy, boost::math::enable_if_t<policies::is_policy_v<Policy>, bool>>
BOOST_MATH_GPU_ENABLED inline typename tools::promote_args<T1, T2>::type
hypot(T1 x, T2 y, const Policy& pol)
{
@@ -73,6 +117,28 @@ BOOST_MATH_GPU_ENABLED inline typename tools::promote_args<T1, T2>::type
static_cast<result_type>(x), static_cast<result_type>(y), pol);
}
template <class T1, class T2, class T3, boost::math::enable_if_t<!policies::is_policy_v<T3>, bool>>
BOOST_MATH_GPU_ENABLED inline tools::promote_args_t<T1, T2, T3>
hypot(T1 x, T2 y, T3 z)
{
using result_type = tools::promote_args_t<T1, T2, T3>;
return detail::hypot_imp(static_cast<result_type>(x),
static_cast<result_type>(y),
static_cast<result_type>(z),
policies::policy<>());
}
template <class T1, class T2, class T3, class Policy>
BOOST_MATH_GPU_ENABLED inline tools::promote_args_t<T1, T2, T3>
hypot(T1 x, T2 y, T3 z, const Policy& pol)
{
using result_type = tools::promote_args_t<T1, T2, T3>;
return detail::hypot_imp(static_cast<result_type>(x),
static_cast<result_type>(y),
static_cast<result_type>(z),
pol);
}
} // namespace math
} // namespace boost

10
include/boost/math/special_functions/math_fwd.hpp
View File

@@ -636,10 +636,18 @@ namespace boost
BOOST_MATH_GPU_ENABLED tools::promote_args_t<T1, T2>
hypot(T1 x, T2 y);
template <class T1, class T2, class Policy>
template <class T1, class T2, class Policy, boost::math::enable_if_t<policies::is_policy_v<Policy>, bool> = true>
BOOST_MATH_GPU_ENABLED tools::promote_args_t<T1, T2>
hypot(T1 x, T2 y, const Policy&);
template <class T1, class T2, class T3, boost::math::enable_if_t<!policies::is_policy_v<T3>, bool> = true>
BOOST_MATH_GPU_ENABLED tools::promote_args_t<T1, T2, T3>
hypot(T1 x, T2 y, T3 z);
template <class T1, class T2, class T3, class Policy>
BOOST_MATH_GPU_ENABLED tools::promote_args_t<T1, T2, T3>
hypot(T1 x, T2 y, T3 z, const Policy& pol);
// cbrt - cube root.
template <class RT>
BOOST_MATH_GPU_ENABLED tools::promote_args_t<RT> cbrt(RT z);

1
test/compile_test/instantiate.hpp
View File

@@ -369,6 +369,7 @@ void instantiate(RealType)
boost::math::jacobi_theta4m1(v1, v2);
boost::math::jacobi_theta4m1tau(v1, v2);
boost::math::hypot(v1, v2);
boost::math::hypot(v1, v2, v3);
boost::math::sinc_pi(v1);
boost::math::sinhc_pi(v1);
boost::math::asinh(v1);

3
test/compile_test/sf_hypot_incl_test.cpp
View File

@@ -16,8 +16,11 @@
void compile_and_link_test()
{
check_result<float>(boost::math::hypot<float>(f, f));
check_result<float>(boost::math::hypot(f, f, f));
check_result<double>(boost::math::hypot<double>(d, d));
check_result<double>(boost::math::hypot<double>(d, d,d));
#ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
check_result<long double>(boost::math::hypot<long double>(l, l));
check_result<long double>(boost::math::hypot<long double>(l, l, l));
#endif
}

69
test/hypot_test.cpp
View File

@@ -10,6 +10,7 @@
#include <boost/test/unit_test.hpp>
#include <boost/test/tools/floating_point_comparison.hpp>
#include <boost/math/special_functions/math_fwd.hpp>
#include <boost/math/special_functions/hypot.hpp>
#include <cmath>
@@ -41,6 +42,22 @@ const float boundaries[] = {
std::sqrt((std::numeric_limits<float>::min)()) * 2,
};
void do_test_boundaries(float x, float y, float z)
{
float expected = static_cast<float>((boost::math::hypot)(
#ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
static_cast<long double>(x),
static_cast<long double>(y),
static_cast<long double>(z)));
#else
static_cast<double>(x),
static_cast<double>(y),
static_cast<double>(z));
#endif
float found = (boost::math::hypot)(x, y, z);
BOOST_CHECK_CLOSE(expected, found, tolerance);
}
void do_test_boundaries(float x, float y)
{
float expected = static_cast<float>((boost::math::hypot)(
@@ -55,12 +72,31 @@ void do_test_boundaries(float x, float y)
BOOST_CHECK_CLOSE(expected, found, tolerance);
}
void test_boundaries(float x, float y, float z)
{
do_test_boundaries(x, y, z);
do_test_boundaries(-x, y, z);
do_test_boundaries(x, -y, z);
do_test_boundaries(x, y, -z);
do_test_boundaries(-x, -y, z);
do_test_boundaries(-x, y, -z);
do_test_boundaries(x, -y, -z);
do_test_boundaries(-x, -y, -z);
}
void test_boundaries(float x, float y)
{
do_test_boundaries(x, y);
do_test_boundaries(-x, y);
do_test_boundaries(-x, -y);
do_test_boundaries(x, -y);
for(unsigned i = 0; i < sizeof(boundaries)/sizeof(float); ++i)
{
test_boundaries(x, y, boundaries[i]);
test_boundaries(x, y, boundaries[i] + std::numeric_limits<float>::epsilon()*boundaries[i]);
test_boundaries(x, y, boundaries[i] - std::numeric_limits<float>::epsilon()*boundaries[i]);
}
}
void test_boundaries(float x)
@@ -92,12 +128,26 @@ void test_spots()
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[i], zero), std::fabs(-boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[i], -zero), std::fabs(boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[i], -zero), std::fabs(-boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[i], zero, zero), std::fabs(boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[i], zero, zero), std::fabs(-boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[i], -zero, zero), std::fabs(boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[i], -zero, zero), std::fabs(-boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(zero, boundaries[i], zero), std::fabs(boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(zero, -boundaries[i], zero), std::fabs(-boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(zero, boundaries[i], -zero), std::fabs(boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(zero, -boundaries[i], -zero), std::fabs(-boundaries[i]));
for(unsigned j = 0; j < sizeof(boundaries)/sizeof(float); ++j)
{
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[i], boundaries[j]), boost::math::hypot(boundaries[j], boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[i], boundaries[j]), boost::math::hypot(boundaries[i], -boundaries[j]));
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[i], -boundaries[j]), boost::math::hypot(-boundaries[j], -boundaries[i]));
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[i], -boundaries[j]), boost::math::hypot(-boundaries[i], boundaries[j]));
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[i], boundaries[j], boundaries[j]), boost::math::hypot(boundaries[j], boundaries[i], boundaries[j]));
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[i], boundaries[j], boundaries[j]), boost::math::hypot(boundaries[i], boundaries[j], -boundaries[j]));
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[i], -boundaries[j], boundaries[j]), boost::math::hypot(-boundaries[j], -boundaries[i], boundaries[j]));
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[i], -boundaries[j], boundaries[j]), boost::math::hypot(-boundaries[i], boundaries[j], boundaries[j]));
}
}
if((std::numeric_limits<float>::has_infinity) && (std::numeric_limits<float>::has_quiet_NaN))
@@ -108,6 +158,16 @@ void test_spots()
BOOST_CHECK_EQUAL(boost::math::hypot(-inf, nan), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(nan, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(nan, -inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(inf, nan, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(-inf, nan, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(nan, inf, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(nan, -inf, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(nan, nan, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(-inf, nan, nan), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(nan, inf, nan), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(nan, -inf, nan), inf);
for(unsigned j = 0; j < sizeof(boundaries)/sizeof(float); ++j)
{
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[j], inf), inf);
@@ -118,6 +178,15 @@ void test_spots()
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[j], -inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(-inf, boundaries[j]), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(-inf, -boundaries[j]), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[j], inf, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[j], inf, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(inf, boundaries[j], inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(inf, -boundaries[j], inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(boundaries[j], -inf, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(-boundaries[j], -inf, inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(-inf, boundaries[j], inf), inf);
BOOST_CHECK_EQUAL(boost::math::hypot(-inf, -boundaries[j], inf), inf);
}
}
}