diff options
author | robot-contrib <robot-contrib@yandex-team.com> | 2023-12-16 23:40:19 +0300 |
---|---|---|
committer | robot-contrib <robot-contrib@yandex-team.com> | 2023-12-17 00:19:38 +0300 |
commit | d69bee97bd68a2689ea0ea66169e7e7472d8a6da (patch) | |
tree | 516262cb2a9131e333072839d2dad97fe0f204ae | |
parent | 993e88069c540e49edba979e9ce8df53cffbdb5b (diff) | |
download | ydb-d69bee97bd68a2689ea0ea66169e7e7472d8a6da.tar.gz |
Update contrib/restricted/boost/math to 1.84.0
17 files changed, 262 insertions, 195 deletions
diff --git a/contrib/restricted/boost/math/include/boost/math/ccmath/abs.hpp b/contrib/restricted/boost/math/include/boost/math/ccmath/abs.hpp index 84e01fba01..db59502ab6 100644 --- a/contrib/restricted/boost/math/include/boost/math/ccmath/abs.hpp +++ b/contrib/restricted/boost/math/include/boost/math/ccmath/abs.hpp @@ -8,22 +8,16 @@ #ifndef BOOST_MATH_CCMATH_ABS #define BOOST_MATH_CCMATH_ABS -#include <cmath> -#include <type_traits> -#include <limits> -#include <boost/math/tools/is_constant_evaluated.hpp> +#include <boost/math/ccmath/detail/config.hpp> + +#ifdef BOOST_MATH_NO_CCMATH +#error "The header <boost/math/abs.hpp> can only be used in C++17 and later." +#endif + #include <boost/math/tools/assert.hpp> #include <boost/math/ccmath/isnan.hpp> #include <boost/math/ccmath/isinf.hpp> -#include <boost/math/tools/is_standalone.hpp> -#ifndef BOOST_MATH_STANDALONE -#include <boost/config.hpp> -#ifdef BOOST_NO_CXX17_IF_CONSTEXPR -#error "The header <boost/math/norms.hpp> can only be used in C++17 and later." -#endif -#endif - namespace boost::math::ccmath { namespace detail { diff --git a/contrib/restricted/boost/math/include/boost/math/ccmath/detail/config.hpp b/contrib/restricted/boost/math/include/boost/math/ccmath/detail/config.hpp new file mode 100644 index 0000000000..61405e781a --- /dev/null +++ b/contrib/restricted/boost/math/include/boost/math/ccmath/detail/config.hpp @@ -0,0 +1,52 @@ +// (C) Copyright John Maddock 2023. +// Use, modification and distribution are subject to the +// Boost Software License, Version 1.0. (See accompanying file +// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt) +// +// Core configuration for ccmath functions, basically will they work or not? + +#ifndef BOOST_MATH_CCMATH_DETAIL_CONFIG +#define BOOST_MATH_CCMATH_DETAIL_CONFIG + +#include <cmath> +#include <type_traits> +#include <limits> +#include <boost/math/tools/is_constant_evaluated.hpp> +#include <boost/math/tools/is_standalone.hpp> + +#ifndef BOOST_MATH_STANDALONE + +#include <boost/config.hpp> +#ifdef BOOST_NO_CXX17_IF_CONSTEXPR +# define BOOST_MATH_NO_CCMATH +#endif + +#else // BOOST_MATH_STANDALONE + +#if defined(_MSC_VER) + +#if defined(_MSVC_LANG) && (_MSVC_LANG < 201703) +# define BOOST_MATH_NO_CCMATH +#endif + +#else // _MSC_VER + +#if (__cplusplus < 201703) +# define BOOST_MATH_NO_CCMATH +#endif + +#endif + +#endif + +#ifndef _MSC_VER +// +// Don't check here for msvc as they didn't get std lib configuration macros at the same time as C++17 <type_traits> +// +#if (__cpp_lib_bool_constant < 201505L) && !defined(BOOST_MATH_NO_CCMATH) +# define BOOST_MATH_NO_CCMATH +#endif +#endif + + +#endif diff --git a/contrib/restricted/boost/math/include/boost/math/ccmath/isinf.hpp b/contrib/restricted/boost/math/include/boost/math/ccmath/isinf.hpp index 7b942dcf72..dfbd1f00af 100644 --- a/contrib/restricted/boost/math/include/boost/math/ccmath/isinf.hpp +++ b/contrib/restricted/boost/math/include/boost/math/ccmath/isinf.hpp @@ -6,18 +6,11 @@ #ifndef BOOST_MATH_CCMATH_ISINF #define BOOST_MATH_CCMATH_ISINF -#include <cmath> -#include <limits> -#include <type_traits> -#include <boost/math/tools/is_constant_evaluated.hpp> #include <boost/math/special_functions/fpclassify.hpp> +#include <boost/math/ccmath/detail/config.hpp> -#include <boost/math/tools/is_standalone.hpp> -#ifndef BOOST_MATH_STANDALONE -#include <boost/config.hpp> -#ifdef BOOST_NO_CXX17_IF_CONSTEXPR -#error "The header <boost/math/norms.hpp> can only be used in C++17 and later." -#endif +#ifdef BOOST_MATH_NO_CCMATH +#error "The header <boost/math/isinf.hpp> can only be used in C++17 and later." #endif namespace boost::math::ccmath { diff --git a/contrib/restricted/boost/math/include/boost/math/ccmath/isnan.hpp b/contrib/restricted/boost/math/include/boost/math/ccmath/isnan.hpp index a5e121a97f..504fa90d2c 100644 --- a/contrib/restricted/boost/math/include/boost/math/ccmath/isnan.hpp +++ b/contrib/restricted/boost/math/include/boost/math/ccmath/isnan.hpp @@ -6,17 +6,11 @@ #ifndef BOOST_MATH_CCMATH_ISNAN #define BOOST_MATH_CCMATH_ISNAN -#include <cmath> -#include <type_traits> -#include <boost/math/tools/is_constant_evaluated.hpp> #include <boost/math/special_functions/fpclassify.hpp> +#include <boost/math/ccmath/detail/config.hpp> -#include <boost/math/tools/is_standalone.hpp> -#ifndef BOOST_MATH_STANDALONE -#include <boost/config.hpp> -#ifdef BOOST_NO_CXX17_IF_CONSTEXPR -#error "The header <boost/math/norms.hpp> can only be used in C++17 and later." -#endif +#ifdef BOOST_MATH_NO_CCMATH +#error "The header <boost/math/isnan.hpp> can only be used in C++17 and later." #endif namespace boost::math::ccmath { diff --git a/contrib/restricted/boost/math/include/boost/math/ccmath/ldexp.hpp b/contrib/restricted/boost/math/include/boost/math/ccmath/ldexp.hpp index 31730b265f..3e2cd3610f 100644 --- a/contrib/restricted/boost/math/include/boost/math/ccmath/ldexp.hpp +++ b/contrib/restricted/boost/math/include/boost/math/ccmath/ldexp.hpp @@ -7,11 +7,13 @@ #ifndef BOOST_MATH_CCMATH_LDEXP_HPP #define BOOST_MATH_CCMATH_LDEXP_HPP -#include <cmath> -#include <limits> -#include <type_traits> +#include <boost/math/ccmath/detail/config.hpp> + +#ifdef BOOST_MATH_NO_CCMATH +#error "The header <boost/math/ldexp.hpp> can only be used in C++17 and later." +#endif + #include <stdexcept> -#include <boost/math/tools/is_constant_evaluated.hpp> #include <boost/math/ccmath/abs.hpp> #include <boost/math/ccmath/isinf.hpp> #include <boost/math/ccmath/isnan.hpp> diff --git a/contrib/restricted/boost/math/include/boost/math/cstdfloat/cstdfloat_cmath.hpp b/contrib/restricted/boost/math/include/boost/math/cstdfloat/cstdfloat_cmath.hpp index 049ac4091e..515f3052a2 100644 --- a/contrib/restricted/boost/math/include/boost/math/cstdfloat/cstdfloat_cmath.hpp +++ b/contrib/restricted/boost/math/include/boost/math/cstdfloat/cstdfloat_cmath.hpp @@ -918,62 +918,6 @@ namespace std using boost::math::cstdfloat::detail::isunordered; // end more functions - // - // Very basic iostream operator: - // - inline std::ostream& operator << (std::ostream& os, __float128 m_value) - { - std::streamsize digits = os.precision(); - std::ios_base::fmtflags f = os.flags(); - std::string s; - - char buf[100]; - std::unique_ptr<char[]> buf2; - std::string format = "%"; - if (f & std::ios_base::showpos) - format += "+"; - if (f & std::ios_base::showpoint) - format += "#"; - format += ".*"; - if (digits == 0) - digits = 36; - format += "Q"; - if (f & std::ios_base::scientific) - format += "e"; - else if (f & std::ios_base::fixed) - format += "f"; - else - format += "g"; - - int v = quadmath_snprintf(buf, 100, format.c_str(), digits, m_value); - - if ((v < 0) || (v >= 99)) - { - int v_max = v; - buf2.reset(new char[v + 3]); - v = quadmath_snprintf(&buf2[0], v_max + 3, format.c_str(), digits, m_value); - if (v >= v_max + 3) - { - BOOST_MATH_THROW_EXCEPTION(std::runtime_error("Formatting of float128_type failed.")); - } - s = &buf2[0]; - } - else - s = buf; - std::streamsize ss = os.width(); - if (ss > static_cast<std::streamsize>(s.size())) - { - char fill = os.fill(); - if ((os.flags() & std::ios_base::left) == std::ios_base::left) - s.append(static_cast<std::string::size_type>(ss - s.size()), fill); - else - s.insert(static_cast<std::string::size_type>(0), static_cast<std::string::size_type>(ss - s.size()), fill); - } - - return os << s; - } - - } // namespace std // We will now remove the preprocessor symbols representing quadruple-precision <cmath> diff --git a/contrib/restricted/boost/math/include/boost/math/cstdfloat/cstdfloat_iostream.hpp b/contrib/restricted/boost/math/include/boost/math/cstdfloat/cstdfloat_iostream.hpp index c01c236f71..7c2c396776 100644 --- a/contrib/restricted/boost/math/include/boost/math/cstdfloat/cstdfloat_iostream.hpp +++ b/contrib/restricted/boost/math/include/boost/math/cstdfloat/cstdfloat_iostream.hpp @@ -32,8 +32,7 @@ #include <boost/math/tools/nothrow.hpp> #include <boost/math/tools/throw_exception.hpp> -// #if (0) - #if defined(__GNUC__) +#if defined(__GNUC__) && !defined(BOOST_MATH_TEST_IO_AS_INTEL_QUAD) // Forward declarations of quadruple-precision string functions. extern "C" int quadmath_snprintf(char *str, size_t size, const char *format, ...) BOOST_MATH_NOTHROW; @@ -96,7 +95,7 @@ // So we have to use dynamic memory allocation for the output // string buffer. - char* my_buffer2 = static_cast<char*>(0U); + char* my_buffer2 = nullptr; #ifndef BOOST_NO_EXCEPTIONS try @@ -160,8 +159,7 @@ } } -// #elif defined(__GNUC__) - #elif defined(__INTEL_COMPILER) +#elif defined(__INTEL_COMPILER) || defined(BOOST_MATH_TEST_IO_AS_INTEL_QUAD) // The section for I/O stream support for the ICC compiler is particularly // long, because these functions must be painstakingly synthesized from @@ -172,6 +170,7 @@ // used in Boost.Multiprecision by John Maddock and Christopher Kormanyos. // This methodology has been slightly modified here for boost::float128_t. + #include <cstring> #include <cctype> @@ -266,7 +265,7 @@ { // Pad out the end with zero's if we need to. - int chars = static_cast<int>(str.size()); + std::ptrdiff_t chars = static_cast<std::ptrdiff_t>(str.size()); chars = digits - chars; if(scientific) @@ -442,7 +441,7 @@ if(isneg) { x = -x; } float_type t; - float_type ten = 10; + constexpr float_type ten = 10; eval_log10(t, x); eval_floor(t, t); @@ -507,6 +506,8 @@ eval_subtract(t, digit); eval_multiply(t, ten); } + if (result.size() == 0) + result = "0"; // Possibly round the result. if(digits >= 0) @@ -522,11 +523,13 @@ if((static_cast<int>(*result.rbegin() - '0') & 1) != 0) { round_string_up_at(result, static_cast<int>(result.size() - 1U), expon); + if (digits == 0) digits = 1; } } else if(cdigit >= 5) { - round_string_up_at(result, static_cast<int>(result.size() - 1), expon); + round_string_up_at(result, static_cast<int>(result.size() - 1u), expon); + if (digits == 0) digits = 1; } } } @@ -569,9 +572,9 @@ { value = 0; - if((p == static_cast<const char*>(0U)) || (*p == static_cast<char>(0))) + if((p == nullptr) || (*p == '\0')) { - return; + return false; } bool is_neg = false; @@ -584,11 +587,11 @@ constexpr int max_digits = std::numeric_limits<float_type>::max_digits10 + 1; - if(*p == static_cast<char>('+')) + if(*p == '+') { ++p; } - else if(*p == static_cast<char>('-')) + else if(*p == '-') { is_neg = true; ++p; @@ -632,7 +635,7 @@ ++digits_seen; } - if(*p == static_cast<char>('.')) + if(*p == '.') { // Grab everything after the point, stop when we've seen // enough digits, even if there are actually more available. @@ -659,15 +662,15 @@ } // Parse the exponent. - if((*p == static_cast<char>('e')) || (*p == static_cast<char>('E'))) + if((*p == 'e') || (*p == 'E')) { ++p; - if(*p == static_cast<char>('+')) + if(*p == '+') { ++p; } - else if(*p == static_cast<char>('-')) + else if(*p == '-') { is_neg_expon = true; ++p; @@ -718,7 +721,7 @@ value = -value; } - return (*p == static_cast<char>(0)); + return (*p == '\0'); } } } } } // boost::math::cstdfloat::detail diff --git a/contrib/restricted/boost/math/include/boost/math/special_functions/beta.hpp b/contrib/restricted/boost/math/include/boost/math/special_functions/beta.hpp index ec824daf72..f15e157324 100644 --- a/contrib/restricted/boost/math/include/boost/math/special_functions/beta.hpp +++ b/contrib/restricted/boost/math/include/boost/math/special_functions/beta.hpp @@ -230,7 +230,10 @@ T ibeta_power_terms(T a, T agh = static_cast<T>(a + Lanczos::g() - 0.5f); T bgh = static_cast<T>(b + Lanczos::g() - 0.5f); T cgh = static_cast<T>(c + Lanczos::g() - 0.5f); - result = Lanczos::lanczos_sum_expG_scaled(c) / (Lanczos::lanczos_sum_expG_scaled(a) * Lanczos::lanczos_sum_expG_scaled(b)); + if ((a < tools::min_value<T>()) || (b < tools::min_value<T>())) + result = 0; // denominator overflows in this case + else + result = Lanczos::lanczos_sum_expG_scaled(c) / (Lanczos::lanczos_sum_expG_scaled(a) * Lanczos::lanczos_sum_expG_scaled(b)); result *= prefix; // combine with the leftover terms from the Lanczos approximation: result *= sqrt(bgh / boost::math::constants::e<T>()); @@ -389,14 +392,15 @@ T ibeta_power_terms(T a, } else { - T p1 = pow(b1, a / b); + // This protects against spurious overflow in a/b: + T p1 = (b1 < 1) && (b < 1) && (tools::max_value<T>() * b < a) ? static_cast<T>(0) : static_cast<T>(pow(b1, a / b)); T l3 = (p1 != 0) && (b2 != 0) ? (log(p1) + log(b2)) * b : tools::max_value<T>(); // arbitrary large value if the logs would fail! if((l3 < tools::log_max_value<T>()) && (l3 > tools::log_min_value<T>())) { result *= pow(p1 * b2, b); } - else + else if(result != 0) // we can elude the calculation below if we're already going to be zero { l2 += l1 + log(result); if(l2 >= tools::log_max_value<T>()) @@ -473,20 +477,108 @@ T ibeta_power_terms(T a, if ((shift_a == 0) && (shift_b == 0)) { T power1, power2; + bool need_logs = false; if (a < b) { - power1 = pow((x * y * c * c) / (a * b), a); - power2 = pow((y * c) / b, b - a); + BOOST_IF_CONSTEXPR(std::numeric_limits<T>::has_infinity) + { + power1 = pow((x * y * c * c) / (a * b), a); + power2 = pow((y * c) / b, b - a); + } + else + { + // We calculate these logs purely so we can check for overflow in the power functions + T l1 = log((x * y * c * c) / (a * b)); + T l2 = log((y * c) / b); + if ((l1 * a > tools::log_min_value<T>()) && (l1 * a < tools::log_max_value<T>()) && (l2 * (b - a) < tools::log_max_value<T>()) && (l2 * (b - a) > tools::log_min_value<T>())) + { + power1 = pow((x * y * c * c) / (a * b), a); + power2 = pow((y * c) / b, b - a); + } + else + { + need_logs = true; + } + } } else { - power1 = pow((x * y * c * c) / (a * b), b); - power2 = pow((x * c) / a, a - b); + BOOST_IF_CONSTEXPR(std::numeric_limits<T>::has_infinity) + { + power1 = pow((x * y * c * c) / (a * b), b); + power2 = pow((x * c) / a, a - b); + } + else + { + // We calculate these logs purely so we can check for overflow in the power functions + T l1 = log((x * y * c * c) / (a * b)) * b; + T l2 = log((x * c) / a) * (a - b); + if ((l1 * a > tools::log_min_value<T>()) && (l1 * a < tools::log_max_value<T>()) && (l2 * (b - a) < tools::log_max_value<T>()) && (l2 * (b - a) > tools::log_min_value<T>())) + { + power1 = pow((x * y * c * c) / (a * b), b); + power2 = pow((x * c) / a, a - b); + } + else + need_logs = true; + } + } + BOOST_IF_CONSTEXPR(std::numeric_limits<T>::has_infinity) + { + if (!(boost::math::isnormal)(power1) || !(boost::math::isnormal)(power2)) + { + need_logs = true; + } } - if (!(boost::math::isnormal)(power1) || !(boost::math::isnormal)(power2)) + if (need_logs) { - // We have to use logs :( - return prefix * exp(a * log(x * c / a) + b * log(y * c / b)) * scaled_tgamma_no_lanczos(c, pol) / (scaled_tgamma_no_lanczos(a, pol) * scaled_tgamma_no_lanczos(b, pol)); + // + // We want: + // + // (xc / a)^a (yc / b)^b + // + // But we know that one or other term will over / underflow and combining the logs will be next to useless as that will cause significant cancellation. + // If we assume b > a and express z ^ b as(z ^ b / a) ^ a with z = (yc / b) then we can move one power term inside the other : + // + // ((xc / a) * (yc / b)^(b / a))^a + // + // However, we're not quite there yet, as the term being exponentiated is quite likely to be close to unity, so let: + // + // xc / a = 1 + (xb - ya) / a + // + // analogously let : + // + // 1 + p = (yc / b) ^ (b / a) = 1 + expm1((b / a) * log1p((ya - xb) / b)) + // + // so putting the two together we have : + // + // exp(a * log1p((xb - ya) / a + p + p(xb - ya) / a)) + // + // Analogously, when a > b we can just swap all the terms around. + // + // Finally, there are a few cases (x or y is unity) when the above logic can't be used + // or where there is no logarithmic cancellation and accuracy is better just using + // the regular formula: + // + T xc_a = x * c / a; + T yc_b = y * c / b; + if ((x == 1) || (y == 1) || (fabs(xc_a - 1) > 0.25) || (fabs(yc_b - 1) > 0.25)) + { + // The above logic fails, the result is almost certainly zero: + power1 = exp(log(xc_a) * a + log(yc_b) * b); + power2 = 1; + } + else if (b > a) + { + T p = boost::math::expm1((b / a) * boost::math::log1p((y * a - x * b) / b)); + power1 = exp(a * boost::math::log1p((x * b - y * a) / a + p * (x * c / a))); + power2 = 1; + } + else + { + T p = boost::math::expm1((a / b) * boost::math::log1p((x * b - y * a) / a)); + power1 = exp(b * boost::math::log1p((y * a - x * b) / b + p * (y * c / b))); + power2 = 1; + } } return prefix * power1 * power2 * scaled_tgamma_no_lanczos(c, pol) / (scaled_tgamma_no_lanczos(a, pol) * scaled_tgamma_no_lanczos(b, pol)); } @@ -568,7 +660,10 @@ T ibeta_series(T a, T b, T x, T s0, const Lanczos&, bool normalised, T* p_deriva T agh = static_cast<T>(a + Lanczos::g() - 0.5f); T bgh = static_cast<T>(b + Lanczos::g() - 0.5f); T cgh = static_cast<T>(c + Lanczos::g() - 0.5f); - result = Lanczos::lanczos_sum_expG_scaled(c) / (Lanczos::lanczos_sum_expG_scaled(a) * Lanczos::lanczos_sum_expG_scaled(b)); + if ((a < tools::min_value<T>()) || (b < tools::min_value<T>())) + result = 0; // denorms cause overflow in the Lanzos series, result will be zero anyway + else + result = Lanczos::lanczos_sum_expG_scaled(c) / (Lanczos::lanczos_sum_expG_scaled(a) * Lanczos::lanczos_sum_expG_scaled(b)); if (!(boost::math::isfinite)(result)) result = 0; @@ -601,10 +696,13 @@ T ibeta_series(T a, T b, T x, T s0, const Lanczos&, bool normalised, T* p_deriva // // Oh dear, we need logs, and this *will* cancel: // - result = log(result) + l1 + l2 + (log(agh) - 1) / 2; - if(p_derivative) - *p_derivative = exp(result + b * log(y)); - result = exp(result); + if (result != 0) // elude calculation when result will be zero. + { + result = log(result) + l1 + l2 + (log(agh) - 1) / 2; + if (p_derivative) + *p_derivative = exp(result + b * log(y)); + result = exp(result); + } } } else @@ -1021,18 +1119,15 @@ T ibeta_imp(T a, T b, T x, const Policy& pol, bool inv, bool normalised, T* p_de BOOST_MATH_ASSERT((p_derivative == 0) || normalised); if(!(boost::math::isfinite)(a)) - return policies::raise_domain_error<T>(function, "The argument a to the incomplete beta function must be >= zero (got a=%1%).", a, pol); + return policies::raise_domain_error<T>(function, "The argument a to the incomplete beta function must be finite (got a=%1%).", a, pol); if(!(boost::math::isfinite)(b)) - return policies::raise_domain_error<T>(function, "The argument b to the incomplete beta function must be >= zero (got b=%1%).", b, pol); - if(!(boost::math::isfinite)(x)) + return policies::raise_domain_error<T>(function, "The argument b to the incomplete beta function must be finite (got b=%1%).", b, pol); + if (!(0 <= x && x <= 1)) return policies::raise_domain_error<T>(function, "The argument x to the incomplete beta function must be in [0,1] (got x=%1%).", x, pol); if(p_derivative) *p_derivative = -1; // value not set. - if((x < 0) || (x > 1)) - return policies::raise_domain_error<T>(function, "Parameter x outside the range [0,1] in the incomplete beta function (got x=%1%).", x, pol); - if(normalised) { if(a < 0) @@ -1422,18 +1517,16 @@ T ibeta_derivative_imp(T a, T b, T x, const Policy& pol) // start with the usual error checks: // if (!(boost::math::isfinite)(a)) - return policies::raise_domain_error<T>(function, "The argument a to the incomplete beta function must be >= zero (got a=%1%).", a, pol); + return policies::raise_domain_error<T>(function, "The argument a to the incomplete beta function must be finite (got a=%1%).", a, pol); if (!(boost::math::isfinite)(b)) - return policies::raise_domain_error<T>(function, "The argument b to the incomplete beta function must be >= zero (got b=%1%).", b, pol); - if (!(boost::math::isfinite)(x)) + return policies::raise_domain_error<T>(function, "The argument b to the incomplete beta function must be finite (got b=%1%).", b, pol); + if (!(0 <= x && x <= 1)) return policies::raise_domain_error<T>(function, "The argument x to the incomplete beta function must be in [0,1] (got x=%1%).", x, pol); if(a <= 0) return policies::raise_domain_error<T>(function, "The argument a to the incomplete beta function must be greater than zero (got a=%1%).", a, pol); if(b <= 0) return policies::raise_domain_error<T>(function, "The argument b to the incomplete beta function must be greater than zero (got b=%1%).", b, pol); - if((x < 0) || (x > 1)) - return policies::raise_domain_error<T>(function, "Parameter x outside the range [0,1] in the incomplete beta function (got x=%1%).", x, pol); // // Now the corner cases: // diff --git a/contrib/restricted/boost/math/include/boost/math/special_functions/detail/ibeta_inverse.hpp b/contrib/restricted/boost/math/include/boost/math/special_functions/detail/ibeta_inverse.hpp index be4bd95399..98cc1d7485 100644 --- a/contrib/restricted/boost/math/include/boost/math/special_functions/detail/ibeta_inverse.hpp +++ b/contrib/restricted/boost/math/include/boost/math/special_functions/detail/ibeta_inverse.hpp @@ -25,23 +25,16 @@ namespace boost{ namespace math{ namespace detail{ template <class T> struct temme_root_finder { - temme_root_finder(const T t_, const T a_) : t(t_), a(a_) {} + temme_root_finder(const T t_, const T a_) : t(t_), a(a_) { + const T x_extrema = 1 / (1 + a); + BOOST_MATH_ASSERT(0 < x_extrema && x_extrema < 1); + } boost::math::tuple<T, T> operator()(T x) { BOOST_MATH_STD_USING // ADL of std names T y = 1 - x; - if(y == 0) - { - T big = tools::max_value<T>() / 4; - return boost::math::make_tuple(static_cast<T>(-big), static_cast<T>(-big)); - } - if(x == 0) - { - T big = tools::max_value<T>() / 4; - return boost::math::make_tuple(static_cast<T>(-big), big); - } T f = log(x) + a * log(y) + t; T f1 = (1 / x) - (a / (y)); return boost::math::make_tuple(f, f1); @@ -410,6 +403,10 @@ T temme_method_3_ibeta_inverse(T a, T b, T p, T q, const Policy& pol) T lower = eta < mu ? cross : 0; T upper = eta < mu ? 1 : cross; T x = (lower + upper) / 2; + + // Early exit for cases with numerical precision issues. + if (cross == 0 || cross == 1) { return cross; } + x = tools::newton_raphson_iterate( temme_root_finder<T>(u, mu), x, lower, upper, policies::digits<T, Policy>() / 2); #ifdef BOOST_INSTRUMENT @@ -828,7 +825,21 @@ T ibeta_inv_imp(T a, T b, T p, T q, const Policy& pol, T* py) std::swap(p, q); invert = !invert; } - if(pow(p, 1/a) < 0.5) + if (a < tools::min_value<T>()) + { + // Avoid spurious overflows for denorms: + if (p < 1) + { + x = 1; + y = 0; + } + else + { + x = 0; + y = 1; + } + } + else if(pow(p, 1/a) < 0.5) { #ifndef BOOST_NO_EXCEPTIONS try diff --git a/contrib/restricted/boost/math/include/boost/math/special_functions/next.hpp b/contrib/restricted/boost/math/include/boost/math/special_functions/next.hpp index 07c55e8d4a..1b56932423 100644 --- a/contrib/restricted/boost/math/include/boost/math/special_functions/next.hpp +++ b/contrib/restricted/boost/math/include/boost/math/special_functions/next.hpp @@ -82,8 +82,8 @@ inline T normalize_value(const T& val, const std::true_type&) } template <class T> -inline T get_smallest_value(std::true_type const&) -{ +inline T get_smallest_value(std::true_type const&) { + static_assert(std::numeric_limits<T>::is_specialized, "Type T must be specialized."); // // numeric_limits lies about denorms being present - particularly // when this can be turned on or off at runtime, as is the case @@ -106,11 +106,12 @@ inline T get_smallest_value(std::false_type const&) template <class T> inline T get_smallest_value() { -#if defined(BOOST_MSVC) && (BOOST_MSVC <= 1310) - return get_smallest_value<T>(std::integral_constant<bool, std::numeric_limits<T>::is_specialized && (std::numeric_limits<T>::has_denorm == 1)>()); -#else - return get_smallest_value<T>(std::integral_constant<bool, std::numeric_limits<T>::is_specialized && (std::numeric_limits<T>::has_denorm == std::denorm_present)>()); -#endif + return get_smallest_value<T>(std::integral_constant<bool, std::numeric_limits<T>::is_specialized>()); +} + +template <class T> +inline bool has_denorm_now() { + return get_smallest_value<T>() < tools::min_value<T>(); } // diff --git a/contrib/restricted/boost/math/include/boost/math/special_functions/round.hpp b/contrib/restricted/boost/math/include/boost/math/special_functions/round.hpp index b32c8e4347..e74acba85b 100644 --- a/contrib/restricted/boost/math/include/boost/math/special_functions/round.hpp +++ b/contrib/restricted/boost/math/include/boost/math/special_functions/round.hpp @@ -12,6 +12,7 @@ #endif #include <boost/math/tools/config.hpp> +#include <boost/math/ccmath/detail/config.hpp> #include <boost/math/policies/error_handling.hpp> #include <boost/math/special_functions/math_fwd.hpp> #include <boost/math/special_functions/fpclassify.hpp> @@ -19,11 +20,9 @@ #include <limits> #include <cmath> -#ifndef BOOST_NO_CXX17_IF_CONSTEXPR +#if !defined(BOOST_MATH_NO_CCMATH) && !defined(BOOST_MATH_NO_CONSTEXPR_DETECTION) #include <boost/math/ccmath/ldexp.hpp> -# if !defined(BOOST_MATH_NO_CONSTEXPR_DETECTION) # define BOOST_MATH_HAS_CONSTEXPR_LDEXP -# endif #endif namespace boost{ namespace math{ diff --git a/contrib/restricted/boost/math/include/boost/math/special_functions/trunc.hpp b/contrib/restricted/boost/math/include/boost/math/special_functions/trunc.hpp index 7ac33a92d3..a084de560b 100644 --- a/contrib/restricted/boost/math/include/boost/math/special_functions/trunc.hpp +++ b/contrib/restricted/boost/math/include/boost/math/special_functions/trunc.hpp @@ -14,14 +14,14 @@ #include <type_traits> #include <boost/math/special_functions/math_fwd.hpp> #include <boost/math/tools/config.hpp> +#include <boost/math/ccmath/detail/config.hpp> #include <boost/math/policies/error_handling.hpp> #include <boost/math/special_functions/fpclassify.hpp> +#include <boost/math/tools/is_constant_evaluated.hpp> -#ifndef BOOST_NO_CXX17_IF_CONSTEXPR +#if !defined(BOOST_MATH_NO_CCMATH) && !defined(BOOST_MATH_NO_CONSTEXPR_DETECTION) #include <boost/math/ccmath/ldexp.hpp> -# if !defined(BOOST_MATH_NO_CONSTEXPR_DETECTION) # define BOOST_MATH_HAS_CONSTEXPR_LDEXP -# endif #endif namespace boost{ namespace math{ namespace detail{ diff --git a/contrib/restricted/boost/math/include/boost/math/tools/big_constant.hpp b/contrib/restricted/boost/math/include/boost/math/tools/big_constant.hpp index bec93a1165..b26d8cf86d 100644 --- a/contrib/restricted/boost/math/include/boost/math/tools/big_constant.hpp +++ b/contrib/restricted/boost/math/include/boost/math/tools/big_constant.hpp @@ -60,7 +60,7 @@ inline T make_big_value(largest_float, const char* s, std::false_type const&, st } #else template <typename T> -inline T make_big_value(largest_float, const char* s, std::false_type const&, std::false_type const&) +inline T make_big_value(largest_float, const char*, std::false_type const&, std::false_type const&) { static_assert(sizeof(T) == 0, "Type is unsupported in standalone mode. Please disable and try again."); } diff --git a/contrib/restricted/boost/math/include/boost/math/tools/config.hpp b/contrib/restricted/boost/math/include/boost/math/tools/config.hpp index 4b4a8a66ba..584cbd5485 100644 --- a/contrib/restricted/boost/math/include/boost/math/tools/config.hpp +++ b/contrib/restricted/boost/math/include/boost/math/tools/config.hpp @@ -146,34 +146,17 @@ # define BOOST_MATH_EXEC_COMPATIBLE #endif -// Attributes from C++14 and newer -#ifdef __has_cpp_attribute - -// C++17 -#if (__cplusplus >= 201703L || _MSVC_LANG >= 201703L) -# if __has_cpp_attribute(maybe_unused) -# define BOOST_MATH_MAYBE_UNUSED [[maybe_unused]] -# endif -#endif - -#endif // Standalone config - -// If attributes are not defined make sure we don't have compiler errors -#ifndef BOOST_MATH_MAYBE_UNUSED -# define BOOST_MATH_MAYBE_UNUSED -#endif - // C++23 #if __cplusplus > 202002L || _MSVC_LANG > 202002L # if __GNUC__ >= 13 // libstdc++3 only defines to/from_chars for std::float128_t when one of these defines are set // otherwise we're right out of luck... # if defined(_GLIBCXX_LDOUBLE_IS_IEEE_BINARY128) || defined(_GLIBCXX_HAVE_FLOAT128_MATH) -# include <cstring> // std::strlen is used with from_chars -# include <charconv> -# error #include <stdfloat> -# define BOOST_MATH_USE_CHARCONV_FOR_CONVERSION -#endif +# include <cstring> // std::strlen is used with from_chars +# include <charconv> +# error #include <stdfloat> +# define BOOST_MATH_USE_CHARCONV_FOR_CONVERSION +# endif # endif #endif diff --git a/contrib/restricted/boost/math/include/boost/math/tools/promotion.hpp b/contrib/restricted/boost/math/include/boost/math/tools/promotion.hpp index 80b3d89a10..263f3a6205 100644 --- a/contrib/restricted/boost/math/include/boost/math/tools/promotion.hpp +++ b/contrib/restricted/boost/math/include/boost/math/tools/promotion.hpp @@ -26,8 +26,12 @@ #include <boost/math/tools/config.hpp> #include <type_traits> -#if __has_include(<stdfloat>) -# error #include <stdfloat> +#if defined __has_include +# if __cplusplus > 202002L || _MSVC_LANG > 202002L +# if __has_include (<stdfloat>) +# error #include <stdfloat> +# endif +# endif #endif namespace boost diff --git a/contrib/restricted/boost/math/include/boost/math/tools/roots.hpp b/contrib/restricted/boost/math/include/boost/math/tools/roots.hpp index 9fd7bee488..af46e2fcdd 100644 --- a/contrib/restricted/boost/math/include/boost/math/tools/roots.hpp +++ b/contrib/restricted/boost/math/include/boost/math/tools/roots.hpp @@ -275,14 +275,7 @@ T newton_raphson_iterate(F f, T guess, T min, T max, int digits, std::uintmax_t& if (fabs(delta * 2) > fabs(delta2)) { // Last two steps haven't converged. - T shift = (delta > 0) ? (result - min) / 2 : (result - max) / 2; - if ((result != 0) && (fabs(shift) > fabs(result))) - { - delta = sign(delta) * fabs(result) * 1.1f; // Protect against huge jumps! - //delta = sign(delta) * result; // Protect against huge jumps! Failed for negative result. https://github.com/boostorg/math/issues/216 - } - else - delta = shift; + delta = (delta > 0) ? (result - min) / 2 : (result - max) / 2; // reset delta1/2 so we don't take this branch next time round: delta1 = 3 * delta; delta2 = 3 * delta; @@ -596,7 +589,8 @@ namespace detail { #ifdef BOOST_MATH_INSTRUMENT std::cout << "Second order root iteration, delta = " << delta << ", residual = " << f0 << "\n"; #endif - T convergence = fabs(delta / delta2); + // We need to avoid delta/delta2 overflowing here: + T convergence = (fabs(delta2) > 1) || (fabs(tools::max_value<T>() * delta2) > fabs(delta)) ? fabs(delta / delta2) : tools::max_value<T>(); if ((convergence > 0.8) && (convergence < 2)) { // last two steps haven't converged. diff --git a/contrib/restricted/boost/math/ya.make b/contrib/restricted/boost/math/ya.make index a692f03ae5..660c3ca286 100644 --- a/contrib/restricted/boost/math/ya.make +++ b/contrib/restricted/boost/math/ya.make @@ -6,9 +6,9 @@ LICENSE(BSL-1.0) LICENSE_TEXTS(.yandex_meta/licenses.list.txt) -VERSION(1.83.0) +VERSION(1.84.0) -ORIGINAL_SOURCE(https://github.com/boostorg/math/archive/boost-1.83.0.tar.gz) +ORIGINAL_SOURCE(https://github.com/boostorg/math/archive/boost-1.84.0.tar.gz) PEERDIR( contrib/restricted/boost/assert |