| 1 | /////////////////////////////////////////////////////////////// |
| 2 | // Copyright 2012-2020 John Maddock. |
| 3 | // Copyright 2020 Madhur Chauhan. |
| 4 | // Copyright 2021 Matt Borland. |
| 5 | // Distributed under the Boost Software License, Version 1.0. |
| 6 | // (See accompanying file LICENSE_1_0.txt or copy at |
| 7 | // https://www.boost.org/LICENSE_1_0.txt) |
| 8 | // |
| 9 | // Comparison operators for cpp_int_backend: |
| 10 | // |
| 11 | #ifndef BOOST_MP_CPP_INT_MISC_HPP |
| 12 | #define BOOST_MP_CPP_INT_MISC_HPP |
| 13 | |
| 14 | #include <boost/multiprecision/detail/standalone_config.hpp> |
| 15 | #include <boost/multiprecision/detail/number_base.hpp> |
| 16 | #include <boost/multiprecision/cpp_int/cpp_int_config.hpp> |
| 17 | #include <boost/multiprecision/detail/float128_functions.hpp> |
| 18 | #include <boost/multiprecision/detail/assert.hpp> |
| 19 | #include <boost/multiprecision/detail/constexpr.hpp> |
| 20 | #include <boost/multiprecision/detail/bitscan.hpp> // lsb etc |
| 21 | #include <boost/multiprecision/detail/hash.hpp> |
| 22 | #include <boost/multiprecision/detail/no_exceptions_support.hpp> |
| 23 | #include <numeric> // std::gcd |
| 24 | #include <type_traits> |
| 25 | #include <stdexcept> |
| 26 | #include <cmath> |
| 27 | |
| 28 | #ifndef BOOST_MP_STANDALONE |
| 29 | #include <boost/integer/common_factor_rt.hpp> |
| 30 | #endif |
| 31 | |
| 32 | #ifdef BOOST_MP_MATH_AVAILABLE |
| 33 | #include <boost/math/special_functions/next.hpp> |
| 34 | #endif |
| 35 | |
| 36 | #ifdef BOOST_MSVC |
| 37 | #pragma warning(push) |
| 38 | #pragma warning(disable : 4702) |
| 39 | #pragma warning(disable : 4127) // conditional expression is constant |
| 40 | #pragma warning(disable : 4146) // unary minus operator applied to unsigned type, result still unsigned |
| 41 | #endif |
| 42 | |
| 43 | // Forward decleration of gcd and lcm functions |
| 44 | namespace boost { namespace multiprecision { namespace detail { |
| 45 | |
| 46 | template <typename T> |
| 47 | inline BOOST_CXX14_CONSTEXPR T constexpr_gcd(T a, T b) noexcept; |
| 48 | |
| 49 | template <typename T> |
| 50 | inline BOOST_CXX14_CONSTEXPR T constexpr_lcm(T a, T b) noexcept; |
| 51 | |
| 52 | }}} // namespace boost::multiprecision::detail |
| 53 | |
| 54 | namespace boost { namespace multiprecision { namespace backends { |
| 55 | |
| 56 | template <class T, bool has_limits = std::numeric_limits<T>::is_specialized> |
| 57 | struct numeric_limits_workaround : public std::numeric_limits<T> |
| 58 | { |
| 59 | }; |
| 60 | template <class R> |
| 61 | struct numeric_limits_workaround<R, false> |
| 62 | { |
| 63 | static constexpr unsigned digits = ~static_cast<R>(0) < 0 ? sizeof(R) * CHAR_BIT - 1 : sizeof(R) * CHAR_BIT; |
| 64 | static constexpr R (min)(){ return (static_cast<R>(-1) < 0) ? static_cast<R>(1) << digits : 0; } |
| 65 | static constexpr R (max)() { return (static_cast<R>(-1) < 0) ? ~(static_cast<R>(1) << digits) : ~static_cast<R>(0); } |
| 66 | }; |
| 67 | |
| 68 | template <class R, class CppInt> |
| 69 | BOOST_MP_CXX14_CONSTEXPR void check_in_range(const CppInt& val, const std::integral_constant<int, checked>&) |
| 70 | { |
| 71 | using cast_type = typename boost::multiprecision::detail::canonical<R, CppInt>::type; |
| 72 | |
| 73 | if (val.sign()) |
| 74 | { |
| 75 | BOOST_IF_CONSTEXPR (boost::multiprecision::detail::is_signed<R>::value == false) |
| 76 | BOOST_MP_THROW_EXCEPTION(std::range_error("Attempt to assign a negative value to an unsigned type." )); |
| 77 | if (val.compare(static_cast<cast_type>((numeric_limits_workaround<R>::min)())) < 0) |
| 78 | BOOST_MP_THROW_EXCEPTION(std::overflow_error("Could not convert to the target type - -value is out of range." )); |
| 79 | } |
| 80 | else |
| 81 | { |
| 82 | if (val.compare(static_cast<cast_type>((numeric_limits_workaround<R>::max)())) > 0) |
| 83 | BOOST_MP_THROW_EXCEPTION(std::overflow_error("Could not convert to the target type - -value is out of range." )); |
| 84 | } |
| 85 | } |
| 86 | template <class R, class CppInt> |
| 87 | inline BOOST_MP_CXX14_CONSTEXPR void check_in_range(const CppInt& /*val*/, const std::integral_constant<int, unchecked>&) noexcept {} |
| 88 | |
| 89 | inline BOOST_MP_CXX14_CONSTEXPR void check_is_negative(const std::integral_constant<bool, true>&) noexcept {} |
| 90 | inline void check_is_negative(const std::integral_constant<bool, false>&) |
| 91 | { |
| 92 | BOOST_MP_THROW_EXCEPTION(std::range_error("Attempt to assign a negative value to an unsigned type." )); |
| 93 | } |
| 94 | |
| 95 | template <class Integer> |
| 96 | inline BOOST_MP_CXX14_CONSTEXPR Integer negate_integer(Integer i, const std::integral_constant<bool, true>&) noexcept |
| 97 | { |
| 98 | return -i; |
| 99 | } |
| 100 | template <class Integer> |
| 101 | inline BOOST_MP_CXX14_CONSTEXPR Integer negate_integer(Integer i, const std::integral_constant<bool, false>&) noexcept |
| 102 | { |
| 103 | return ~(i - 1); |
| 104 | } |
| 105 | |
| 106 | template <class R, std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 107 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_integral<R>::value && !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, void>::type |
| 108 | eval_convert_to(R* result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& backend) |
| 109 | { |
| 110 | using checked_type = std::integral_constant<int, Checked1>; |
| 111 | check_in_range<R>(backend, checked_type()); |
| 112 | |
| 113 | BOOST_IF_CONSTEXPR(numeric_limits_workaround<R>::digits < cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits) |
| 114 | { |
| 115 | if ((backend.sign() && boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value) && (1 + static_cast<boost::multiprecision::limb_type>((std::numeric_limits<R>::max)()) <= backend.limbs()[0])) |
| 116 | { |
| 117 | *result = (numeric_limits_workaround<R>::min)(); |
| 118 | return; |
| 119 | } |
| 120 | else if (boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value && !backend.sign() && static_cast<boost::multiprecision::limb_type>((std::numeric_limits<R>::max)()) <= backend.limbs()[0]) |
| 121 | { |
| 122 | *result = (numeric_limits_workaround<R>::max)(); |
| 123 | return; |
| 124 | } |
| 125 | else |
| 126 | *result = static_cast<R>(backend.limbs()[0]); |
| 127 | } |
| 128 | else |
| 129 | *result = static_cast<R>(backend.limbs()[0]); |
| 130 | |
| 131 | BOOST_IF_CONSTEXPR(numeric_limits_workaround<R>::digits > cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits) |
| 132 | { |
| 133 | std::size_t shift = cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 134 | std::size_t i = 1u; |
| 135 | |
| 136 | while ((i < backend.size()) && (shift < static_cast<unsigned>(numeric_limits_workaround<R>::digits - cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits))) |
| 137 | { |
| 138 | *result += static_cast<R>(backend.limbs()[i]) << shift; |
| 139 | shift += cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 140 | ++i; |
| 141 | } |
| 142 | // |
| 143 | // We have one more limb to extract, but may not need all the bits, so treat this as a special case: |
| 144 | // |
| 145 | if (i < backend.size()) |
| 146 | { |
| 147 | const limb_type mask = ((numeric_limits_workaround<R>::digits - shift) == cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits) ? ~static_cast<limb_type>(0) : static_cast<limb_type>(static_cast<limb_type>(1u) << (numeric_limits_workaround<R>::digits - shift)) - 1u; |
| 148 | const limb_type limb_at_index_masked = static_cast<limb_type>(backend.limbs()[i] & mask); |
| 149 | |
| 150 | *result = static_cast<R>(*result + static_cast<R>(static_cast<R>(limb_at_index_masked) << shift)); |
| 151 | |
| 152 | if ((backend.limbs()[i] & static_cast<limb_type>(~mask)) || (i + 1 < backend.size())) |
| 153 | { |
| 154 | // Overflow: |
| 155 | if (backend.sign()) |
| 156 | { |
| 157 | check_is_negative(boost::multiprecision::detail::is_signed<R>()); |
| 158 | *result = (numeric_limits_workaround<R>::min)(); |
| 159 | } |
| 160 | else if (boost::multiprecision::detail::is_signed<R>::value) |
| 161 | *result = (numeric_limits_workaround<R>::max)(); |
| 162 | return; |
| 163 | } |
| 164 | } |
| 165 | } |
| 166 | else if (backend.size() > 1) |
| 167 | { |
| 168 | // Overflow: |
| 169 | if (backend.sign()) |
| 170 | { |
| 171 | check_is_negative(boost::multiprecision::detail::is_signed<R>()); |
| 172 | *result = (numeric_limits_workaround<R>::min)(); |
| 173 | } |
| 174 | else if (boost::multiprecision::detail::is_signed<R>::value) |
| 175 | *result = (numeric_limits_workaround<R>::max)(); |
| 176 | return; |
| 177 | } |
| 178 | if (backend.sign()) |
| 179 | { |
| 180 | check_is_negative(std::integral_constant<bool, boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value>()); |
| 181 | *result = negate_integer(*result, std::integral_constant<bool, boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value>()); |
| 182 | } |
| 183 | } |
| 184 | |
| 185 | template <class R, std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 186 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<std::is_floating_point<R>::value && !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, void>::type |
| 187 | eval_convert_to(R* result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& backend) noexcept(boost::multiprecision::detail::is_arithmetic<R>::value && std::numeric_limits<R>::has_infinity) |
| 188 | { |
| 189 | BOOST_MP_FLOAT128_USING using std::ldexp; |
| 190 | if (eval_is_zero(backend)) |
| 191 | { |
| 192 | *result = 0.0f; |
| 193 | return; |
| 194 | } |
| 195 | |
| 196 | #ifdef BOOST_HAS_FLOAT128 |
| 197 | std::ptrdiff_t bits_to_keep = static_cast<std::ptrdiff_t>(std::is_same<R, float128_type>::value ? 113 : std::numeric_limits<R>::digits); |
| 198 | #else |
| 199 | std::ptrdiff_t bits_to_keep = static_cast<std::ptrdiff_t>(std::numeric_limits<R>::digits); |
| 200 | #endif |
| 201 | std::ptrdiff_t bits = static_cast<std::ptrdiff_t>(eval_msb_imp(backend) + 1); |
| 202 | |
| 203 | if (bits > bits_to_keep) |
| 204 | { |
| 205 | // Extract the bits we need, and then manually round the result: |
| 206 | *result = 0.0f; |
| 207 | typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::const_limb_pointer p = backend.limbs(); |
| 208 | limb_type mask = ~static_cast<limb_type>(0u); |
| 209 | std::size_t index = backend.size() - 1; |
| 210 | std::size_t shift = cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits * index; |
| 211 | while (bits_to_keep > 0) |
| 212 | { |
| 213 | if (bits_to_keep < (std::ptrdiff_t)cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits) |
| 214 | { |
| 215 | if(index != backend.size() - 1) |
| 216 | { |
| 217 | const std::ptrdiff_t left_shift_amount = static_cast<std::ptrdiff_t>(static_cast<std::ptrdiff_t>(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits) - bits_to_keep); |
| 218 | |
| 219 | mask <<= left_shift_amount; |
| 220 | } |
| 221 | else |
| 222 | { |
| 223 | std::ptrdiff_t bits_in_first_limb = static_cast<std::ptrdiff_t>(bits % static_cast<std::ptrdiff_t>(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits)); |
| 224 | if (bits_in_first_limb == 0) |
| 225 | bits_in_first_limb = cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 226 | if (bits_in_first_limb > bits_to_keep) |
| 227 | mask <<= bits_in_first_limb - bits_to_keep; |
| 228 | } |
| 229 | } |
| 230 | *result += ldexp(static_cast<R>(p[index] & mask), static_cast<int>(shift)); |
| 231 | shift -= cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 232 | |
| 233 | const bool bits_has_non_zero_remainder = (bits % static_cast<std::ptrdiff_t>(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits) != 0); |
| 234 | |
| 235 | bits_to_keep -= ((index == backend.size() - 1) && bits_has_non_zero_remainder) |
| 236 | ? bits % static_cast<std::ptrdiff_t>(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits) |
| 237 | : static_cast<std::ptrdiff_t>(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits); |
| 238 | --index; |
| 239 | } |
| 240 | // Perform rounding: |
| 241 | bits -= 1 + std::numeric_limits<R>::digits; |
| 242 | if (eval_bit_test(backend, static_cast<unsigned>(bits))) |
| 243 | { |
| 244 | if ((eval_lsb_imp(backend) < static_cast<std::size_t>(bits)) || eval_bit_test(backend, static_cast<std::size_t>(bits + 1))) |
| 245 | { |
| 246 | #ifdef BOOST_MP_MATH_AVAILABLE |
| 247 | BOOST_IF_CONSTEXPR(std::numeric_limits<R>::has_infinity) |
| 248 | { |
| 249 | // Must NOT throw: |
| 250 | *result = boost::math::float_next(*result, boost::math::policies::make_policy(boost::math::policies::overflow_error<boost::math::policies::ignore_error>())); |
| 251 | } |
| 252 | else |
| 253 | { |
| 254 | *result = boost::math::float_next(*result); |
| 255 | } |
| 256 | #else |
| 257 | using std::nextafter; BOOST_MP_FLOAT128_USING |
| 258 | *result = nextafter(*result, *result * 2); |
| 259 | #endif |
| 260 | } |
| 261 | } |
| 262 | } |
| 263 | else |
| 264 | { |
| 265 | typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::const_limb_pointer p = backend.limbs(); |
| 266 | std::size_t shift = cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 267 | *result = static_cast<R>(*p); |
| 268 | for (std::size_t i = 1; i < backend.size(); ++i) |
| 269 | { |
| 270 | *result += static_cast<R>(ldexp(static_cast<long double>(p[i]), static_cast<int>(shift))); |
| 271 | shift += cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 272 | } |
| 273 | } |
| 274 | if (backend.sign()) |
| 275 | *result = -*result; |
| 276 | } |
| 277 | |
| 278 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 279 | BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, bool>::type |
| 280 | eval_is_zero(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) noexcept |
| 281 | { |
| 282 | return (val.size() == 1) && (val.limbs()[0] == 0); |
| 283 | } |
| 284 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 285 | BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, int>::type |
| 286 | eval_get_sign(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) noexcept |
| 287 | { |
| 288 | return eval_is_zero(val) ? 0 : val.sign() ? -1 : 1; |
| 289 | } |
| 290 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 291 | BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 292 | eval_abs(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value)) |
| 293 | { |
| 294 | result = val; |
| 295 | result.sign(false); |
| 296 | } |
| 297 | |
| 298 | // |
| 299 | // Get the location of the least-significant-bit: |
| 300 | // |
| 301 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 302 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, std::size_t>::type |
| 303 | eval_lsb_imp(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a) |
| 304 | { |
| 305 | // |
| 306 | // Find the index of the least significant limb that is non-zero: |
| 307 | // |
| 308 | std::size_t index = 0; |
| 309 | while (!a.limbs()[index] && (index < a.size())) |
| 310 | ++index; |
| 311 | // |
| 312 | // Find the index of the least significant bit within that limb: |
| 313 | // |
| 314 | std::size_t result = boost::multiprecision::detail::find_lsb(a.limbs()[index]); |
| 315 | |
| 316 | return result + index * cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 317 | } |
| 318 | |
| 319 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 320 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, std::size_t>::type |
| 321 | eval_lsb(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a) |
| 322 | { |
| 323 | using default_ops::eval_get_sign; |
| 324 | if (eval_get_sign(a) == 0) |
| 325 | { |
| 326 | BOOST_MP_THROW_EXCEPTION(std::domain_error("No bits were set in the operand." )); |
| 327 | } |
| 328 | if (a.sign()) |
| 329 | { |
| 330 | BOOST_MP_THROW_EXCEPTION(std::domain_error("Testing individual bits in negative values is not supported - results are undefined." )); |
| 331 | } |
| 332 | return eval_lsb_imp(a); |
| 333 | } |
| 334 | |
| 335 | // |
| 336 | // Get the location of the most-significant-bit: |
| 337 | // |
| 338 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 339 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, std::size_t>::type |
| 340 | eval_msb_imp(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a) |
| 341 | { |
| 342 | // |
| 343 | // Find the index of the most significant bit that is non-zero: |
| 344 | // |
| 345 | return (a.size() - 1) * cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits + boost::multiprecision::detail::find_msb(a.limbs()[a.size() - 1]); |
| 346 | } |
| 347 | |
| 348 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 349 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, std::size_t>::type |
| 350 | eval_msb(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a) |
| 351 | { |
| 352 | using default_ops::eval_get_sign; |
| 353 | if (eval_get_sign(a) == 0) |
| 354 | { |
| 355 | BOOST_MP_THROW_EXCEPTION(std::domain_error("No bits were set in the operand." )); |
| 356 | } |
| 357 | if (a.sign()) |
| 358 | { |
| 359 | BOOST_MP_THROW_EXCEPTION(std::domain_error("Testing individual bits in negative values is not supported - results are undefined." )); |
| 360 | } |
| 361 | return eval_msb_imp(a); |
| 362 | } |
| 363 | |
| 364 | #ifdef BOOST_GCC |
| 365 | // |
| 366 | // We really shouldn't need to be disabling this warning, but it really does appear to be |
| 367 | // spurious. The warning appears only when in release mode, and asserts are on. |
| 368 | // |
| 369 | #pragma GCC diagnostic push |
| 370 | #pragma GCC diagnostic ignored "-Warray-bounds" |
| 371 | #endif |
| 372 | |
| 373 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 374 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, bool>::type |
| 375 | eval_bit_test(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val, std::size_t index) noexcept |
| 376 | { |
| 377 | std::size_t offset = index / cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 378 | std::size_t shift = index % cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 379 | limb_type mask = shift ? limb_type(1u) << shift : limb_type(1u); |
| 380 | if (offset >= val.size()) |
| 381 | return false; |
| 382 | return val.limbs()[offset] & mask ? true : false; |
| 383 | } |
| 384 | |
| 385 | #ifdef BOOST_GCC |
| 386 | #pragma GCC diagnostic pop |
| 387 | #endif |
| 388 | |
| 389 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 390 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 391 | eval_bit_set(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val, std::size_t index) |
| 392 | { |
| 393 | std::size_t offset = index / cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 394 | std::size_t shift = index % cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 395 | limb_type mask = shift ? limb_type(1u) << shift : limb_type(1u); |
| 396 | if (offset >= val.size()) |
| 397 | { |
| 398 | std::size_t os = val.size(); |
| 399 | val.resize(offset + 1, offset + 1); |
| 400 | if (offset >= val.size()) |
| 401 | return; // fixed precision overflow |
| 402 | for (std::size_t i = os; i <= offset; ++i) |
| 403 | val.limbs()[i] = 0; |
| 404 | } |
| 405 | val.limbs()[offset] |= mask; |
| 406 | } |
| 407 | |
| 408 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 409 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 410 | eval_bit_unset(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val, std::size_t index) noexcept |
| 411 | { |
| 412 | std::size_t offset = index / cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 413 | std::size_t shift = index % cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 414 | limb_type mask = shift ? limb_type(1u) << shift : limb_type(1u); |
| 415 | if (offset >= val.size()) |
| 416 | return; |
| 417 | val.limbs()[offset] &= ~mask; |
| 418 | val.normalize(); |
| 419 | } |
| 420 | |
| 421 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 422 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 423 | eval_bit_flip(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val, std::size_t index) |
| 424 | { |
| 425 | std::size_t offset = index / cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 426 | std::size_t shift = index % cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits; |
| 427 | limb_type mask = shift ? limb_type(1u) << shift : limb_type(1u); |
| 428 | if (offset >= val.size()) |
| 429 | { |
| 430 | std::size_t os = val.size(); |
| 431 | val.resize(offset + 1, offset + 1); |
| 432 | if (offset >= val.size()) |
| 433 | return; // fixed precision overflow |
| 434 | for (std::size_t i = os; i <= offset; ++i) |
| 435 | val.limbs()[i] = 0; |
| 436 | } |
| 437 | val.limbs()[offset] ^= mask; |
| 438 | val.normalize(); |
| 439 | } |
| 440 | |
| 441 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 442 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 443 | eval_qr( |
| 444 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& x, |
| 445 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& y, |
| 446 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& q, |
| 447 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& r) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value)) |
| 448 | { |
| 449 | divide_unsigned_helper(&q, x, y, r); |
| 450 | q.sign(x.sign() != y.sign()); |
| 451 | r.sign(x.sign()); |
| 452 | } |
| 453 | |
| 454 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 455 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 456 | eval_qr( |
| 457 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& x, |
| 458 | limb_type y, |
| 459 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& q, |
| 460 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& r) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value)) |
| 461 | { |
| 462 | divide_unsigned_helper(&q, x, y, r); |
| 463 | q.sign(x.sign()); |
| 464 | r.sign(x.sign()); |
| 465 | } |
| 466 | |
| 467 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class U> |
| 468 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_integral<U>::value>::type eval_qr( |
| 469 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& x, |
| 470 | U y, |
| 471 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& q, |
| 472 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& r) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value)) |
| 473 | { |
| 474 | using default_ops::eval_qr; |
| 475 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t; |
| 476 | t = y; |
| 477 | eval_qr(x, t, q, r); |
| 478 | } |
| 479 | |
| 480 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class Integer> |
| 481 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_unsigned<Integer>::value && !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, Integer>::type |
| 482 | eval_integer_modulus(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, Integer mod) |
| 483 | { |
| 484 | BOOST_IF_CONSTEXPR (sizeof(Integer) <= sizeof(limb_type)) |
| 485 | { |
| 486 | if (mod <= (std::numeric_limits<limb_type>::max)()) |
| 487 | { |
| 488 | const std::ptrdiff_t n = a.size(); |
| 489 | const double_limb_type two_n_mod = static_cast<limb_type>(1u) + (~static_cast<limb_type>(0u) - mod) % mod; |
| 490 | limb_type res = a.limbs()[n - 1] % mod; |
| 491 | |
| 492 | for (std::ptrdiff_t i = n - 2; i >= 0; --i) |
| 493 | res = static_cast<limb_type>((res * two_n_mod + a.limbs()[i]) % mod); |
| 494 | return res; |
| 495 | } |
| 496 | else |
| 497 | return default_ops::eval_integer_modulus(a, mod); |
| 498 | } |
| 499 | else |
| 500 | { |
| 501 | return default_ops::eval_integer_modulus(a, mod); |
| 502 | } |
| 503 | } |
| 504 | |
| 505 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class Integer> |
| 506 | BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_signed<Integer>::value && boost::multiprecision::detail::is_integral<Integer>::value && !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, Integer>::type |
| 507 | eval_integer_modulus(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& x, Integer val) |
| 508 | { |
| 509 | return eval_integer_modulus(x, boost::multiprecision::detail::unsigned_abs(val)); |
| 510 | } |
| 511 | |
| 512 | BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR limb_type eval_gcd(limb_type u, limb_type v) |
| 513 | { |
| 514 | // boundary cases |
| 515 | if (!u || !v) |
| 516 | return u | v; |
| 517 | #if (defined(__cpp_lib_gcd_lcm) && (__cpp_lib_gcd_lcm >= 201606L)) |
| 518 | return std::gcd(m: u, n: v); |
| 519 | #else |
| 520 | std::size_t shift = boost::multiprecision::detail::find_lsb(u | v); |
| 521 | u >>= boost::multiprecision::detail::find_lsb(u); |
| 522 | do |
| 523 | { |
| 524 | v >>= boost::multiprecision::detail::find_lsb(v); |
| 525 | if (u > v) |
| 526 | std_constexpr::swap(u, v); |
| 527 | v -= u; |
| 528 | } while (v); |
| 529 | return u << shift; |
| 530 | #endif |
| 531 | } |
| 532 | |
| 533 | inline BOOST_MP_CXX14_CONSTEXPR double_limb_type eval_gcd(double_limb_type u, double_limb_type v) |
| 534 | { |
| 535 | #if (defined(__cpp_lib_gcd_lcm) && (__cpp_lib_gcd_lcm >= 201606L)) && (!defined(BOOST_HAS_INT128) || !defined(__STRICT_ANSI__)) |
| 536 | return std::gcd(m: u, n: v); |
| 537 | #else |
| 538 | if (u == 0) |
| 539 | return v; |
| 540 | |
| 541 | std::size_t shift = boost::multiprecision::detail::find_lsb(u | v); |
| 542 | u >>= boost::multiprecision::detail::find_lsb(u); |
| 543 | do |
| 544 | { |
| 545 | v >>= boost::multiprecision::detail::find_lsb(v); |
| 546 | if (u > v) |
| 547 | std_constexpr::swap(u, v); |
| 548 | v -= u; |
| 549 | } while (v); |
| 550 | return u << shift; |
| 551 | #endif |
| 552 | } |
| 553 | |
| 554 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 555 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 556 | eval_gcd( |
| 557 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, |
| 558 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, |
| 559 | limb_type b) |
| 560 | { |
| 561 | int s = eval_get_sign(a); |
| 562 | if (!b || !s) |
| 563 | { |
| 564 | result = a; |
| 565 | *result.limbs() |= b; |
| 566 | } |
| 567 | else |
| 568 | { |
| 569 | eval_modulus(result, a, b); |
| 570 | limb_type& res = *result.limbs(); |
| 571 | res = eval_gcd(u: res, v: b); |
| 572 | } |
| 573 | result.sign(false); |
| 574 | } |
| 575 | |
| 576 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 577 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 578 | eval_gcd( |
| 579 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, |
| 580 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, |
| 581 | double_limb_type b) |
| 582 | { |
| 583 | int s = eval_get_sign(a); |
| 584 | if (!b || !s) |
| 585 | { |
| 586 | if (!s) |
| 587 | result = b; |
| 588 | else |
| 589 | result = a; |
| 590 | return; |
| 591 | } |
| 592 | double_limb_type res = 0; |
| 593 | if(a.sign() == 0) |
| 594 | res = eval_integer_modulus(a, b); |
| 595 | else |
| 596 | { |
| 597 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t(a); |
| 598 | t.negate(); |
| 599 | res = eval_integer_modulus(t, b); |
| 600 | } |
| 601 | res = eval_gcd(u: res, v: b); |
| 602 | result = res; |
| 603 | result.sign(false); |
| 604 | } |
| 605 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 606 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 607 | eval_gcd( |
| 608 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, |
| 609 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, |
| 610 | signed_double_limb_type v) |
| 611 | { |
| 612 | eval_gcd(result, a, static_cast<double_limb_type>(v < 0 ? -v : v)); |
| 613 | } |
| 614 | // |
| 615 | // These 2 overloads take care of gcd against an (unsigned) short etc: |
| 616 | // |
| 617 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class Integer> |
| 618 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_unsigned<Integer>::value && (sizeof(Integer) <= sizeof(limb_type)) && !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 619 | eval_gcd( |
| 620 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, |
| 621 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, |
| 622 | const Integer& v) |
| 623 | { |
| 624 | eval_gcd(result, a, static_cast<limb_type>(v)); |
| 625 | } |
| 626 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class Integer> |
| 627 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_signed<Integer>::value && boost::multiprecision::detail::is_integral<Integer>::value && (sizeof(Integer) <= sizeof(limb_type)) && !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 628 | eval_gcd( |
| 629 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, |
| 630 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, |
| 631 | const Integer& v) |
| 632 | { |
| 633 | eval_gcd(result, a, static_cast<limb_type>(v < 0 ? -v : v)); |
| 634 | } |
| 635 | // |
| 636 | // What follows is Lehmer's GCD algorithm: |
| 637 | // Essentially this uses the leading digit(s) of U and V |
| 638 | // only to run a "simulated" Euclid algorithm. It stops |
| 639 | // when the calculated quotient differs from what would have been |
| 640 | // the true quotient. At that point the cosequences are used to |
| 641 | // calculate the new U and V. A nice lucid description appears |
| 642 | // in "An Analysis of Lehmer's Euclidean GCD Algorithm", |
| 643 | // by Jonathan Sorenson. https://www.researchgate.net/publication/2424634_An_Analysis_of_Lehmer%27s_Euclidean_GCD_Algorithm |
| 644 | // DOI: 10.1145/220346.220378. |
| 645 | // |
| 646 | // There are two versions of this algorithm here, and both are "double digit" |
| 647 | // variations: which is to say if there are k bits per limb, then they extract |
| 648 | // 2k bits into a double_limb_type and then run the algorithm on that. The first |
| 649 | // version is a straightforward version of the algorithm, and is designed for |
| 650 | // situations where double_limb_type is a native integer (for example where |
| 651 | // limb_type is a 32-bit integer on a 64-bit machine). For 32-bit limbs it |
| 652 | // reduces the size of U by about 30 bits per call. The second is a more complex |
| 653 | // version for situations where double_limb_type is a synthetic type: for example |
| 654 | // __int128. For 64 bit limbs it reduces the size of U by about 62 bits per call. |
| 655 | // |
| 656 | // The complexity of the algorithm given by Sorenson is roughly O(ln^2(N)) for |
| 657 | // two N bit numbers. |
| 658 | // |
| 659 | // The original double-digit version of the algorithm is described in: |
| 660 | // |
| 661 | // "A Double Digit Lehmer-Euclid Algorithm for Finding the GCD of Long Integers", |
| 662 | // Tudor Jebelean, J Symbolic Computation, 1995 (19), 145. |
| 663 | // |
| 664 | #ifndef BOOST_HAS_INT128 |
| 665 | // |
| 666 | // When double_limb_type is a native integer type then we should just use it and not worry about the consequences. |
| 667 | // This can eliminate approximately a full limb with each call. |
| 668 | // |
| 669 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class Storage> |
| 670 | void eval_gcd_lehmer(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& U, cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& V, std::size_t lu, Storage& storage) |
| 671 | { |
| 672 | // |
| 673 | // Extract the leading 2 * bits_per_limb bits from U and V: |
| 674 | // |
| 675 | std::size_t h = lu % bits_per_limb; |
| 676 | double_limb_type u = (static_cast<double_limb_type>((U.limbs()[U.size() - 1])) << bits_per_limb) | U.limbs()[U.size() - 2]; |
| 677 | double_limb_type v = (static_cast<double_limb_type>((V.size() < U.size() ? 0 : V.limbs()[V.size() - 1])) << bits_per_limb) | V.limbs()[U.size() - 2]; |
| 678 | if (h) |
| 679 | { |
| 680 | u <<= bits_per_limb - h; |
| 681 | u |= U.limbs()[U.size() - 3] >> h; |
| 682 | v <<= bits_per_limb - h; |
| 683 | v |= V.limbs()[U.size() - 3] >> h; |
| 684 | } |
| 685 | // |
| 686 | // Co-sequences x an y: we need only the last 3 values of these, |
| 687 | // the first 2 values are known correct, the third gets checked |
| 688 | // in each loop operation, and we terminate when they go wrong. |
| 689 | // |
| 690 | // x[i+0] is positive for even i. |
| 691 | // y[i+0] is positive for odd i. |
| 692 | // |
| 693 | // However we track only absolute values here: |
| 694 | // |
| 695 | double_limb_type x[3] = {1, 0}; |
| 696 | double_limb_type y[3] = {0, 1}; |
| 697 | std::size_t i = 0; |
| 698 | |
| 699 | #ifdef BOOST_MP_GCD_DEBUG |
| 700 | cpp_int UU, VV; |
| 701 | UU = U; |
| 702 | VV = V; |
| 703 | #endif |
| 704 | |
| 705 | while (true) |
| 706 | { |
| 707 | double_limb_type q = u / v; |
| 708 | x[2] = x[0] + q * x[1]; |
| 709 | y[2] = y[0] + q * y[1]; |
| 710 | double_limb_type tu = u; |
| 711 | u = v; |
| 712 | v = tu - q * v; |
| 713 | ++i; |
| 714 | // |
| 715 | // We must make sure that y[2] occupies a single limb otherwise |
| 716 | // the multiprecision multiplications below would be much more expensive. |
| 717 | // This can sometimes lose us one iteration, but is worth it for improved |
| 718 | // calculation efficiency. |
| 719 | // |
| 720 | if (y[2] >> bits_per_limb) |
| 721 | break; |
| 722 | // |
| 723 | // These are Jebelean's exact termination conditions: |
| 724 | // |
| 725 | if ((i & 1u) == 0) |
| 726 | { |
| 727 | BOOST_MP_ASSERT(u > v); |
| 728 | if ((v < x[2]) || ((u - v) < (y[2] + y[1]))) |
| 729 | break; |
| 730 | } |
| 731 | else |
| 732 | { |
| 733 | BOOST_MP_ASSERT(u > v); |
| 734 | if ((v < y[2]) || ((u - v) < (x[2] + x[1]))) |
| 735 | break; |
| 736 | } |
| 737 | #ifdef BOOST_MP_GCD_DEBUG |
| 738 | BOOST_MP_ASSERT(q == UU / VV); |
| 739 | UU %= VV; |
| 740 | UU.swap(VV); |
| 741 | #endif |
| 742 | x[0] = x[1]; |
| 743 | x[1] = x[2]; |
| 744 | y[0] = y[1]; |
| 745 | y[1] = y[2]; |
| 746 | } |
| 747 | if (i == 1) |
| 748 | { |
| 749 | // No change to U and V we've stalled! |
| 750 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t; |
| 751 | eval_modulus(t, U, V); |
| 752 | U.swap(V); |
| 753 | V.swap(t); |
| 754 | return; |
| 755 | } |
| 756 | // |
| 757 | // Update U and V. |
| 758 | // We have: |
| 759 | // |
| 760 | // U = x[0]U + y[0]V and |
| 761 | // V = x[1]U + y[1]V. |
| 762 | // |
| 763 | // But since we track only absolute values of x and y |
| 764 | // we have to take account of the implied signs and perform |
| 765 | // the appropriate subtraction depending on the whether i is |
| 766 | // even or odd: |
| 767 | // |
| 768 | std::size_t ts = U.size() + 1; |
| 769 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t1(storage, ts), t2(storage, ts), t3(storage, ts); |
| 770 | eval_multiply(t1, U, static_cast<limb_type>(x[0])); |
| 771 | eval_multiply(t2, V, static_cast<limb_type>(y[0])); |
| 772 | eval_multiply(t3, U, static_cast<limb_type>(x[1])); |
| 773 | if ((i & 1u) == 0) |
| 774 | { |
| 775 | if (x[0] == 0) |
| 776 | U = t2; |
| 777 | else |
| 778 | { |
| 779 | BOOST_MP_ASSERT(t2.compare(t1) >= 0); |
| 780 | eval_subtract(U, t2, t1); |
| 781 | BOOST_MP_ASSERT(U.sign() == false); |
| 782 | } |
| 783 | } |
| 784 | else |
| 785 | { |
| 786 | BOOST_MP_ASSERT(t1.compare(t2) >= 0); |
| 787 | eval_subtract(U, t1, t2); |
| 788 | BOOST_MP_ASSERT(U.sign() == false); |
| 789 | } |
| 790 | eval_multiply(t2, V, static_cast<limb_type>(y[1])); |
| 791 | if (i & 1u) |
| 792 | { |
| 793 | if (x[1] == 0) |
| 794 | V = t2; |
| 795 | else |
| 796 | { |
| 797 | BOOST_MP_ASSERT(t2.compare(t3) >= 0); |
| 798 | eval_subtract(V, t2, t3); |
| 799 | BOOST_MP_ASSERT(V.sign() == false); |
| 800 | } |
| 801 | } |
| 802 | else |
| 803 | { |
| 804 | BOOST_MP_ASSERT(t3.compare(t2) >= 0); |
| 805 | eval_subtract(V, t3, t2); |
| 806 | BOOST_MP_ASSERT(V.sign() == false); |
| 807 | } |
| 808 | BOOST_MP_ASSERT(U.compare(V) >= 0); |
| 809 | BOOST_MP_ASSERT(lu > eval_msb(U)); |
| 810 | #ifdef BOOST_MP_GCD_DEBUG |
| 811 | |
| 812 | BOOST_MP_ASSERT(UU == U); |
| 813 | BOOST_MP_ASSERT(VV == V); |
| 814 | |
| 815 | extern std::size_t total_lehmer_gcd_calls; |
| 816 | extern std::size_t total_lehmer_gcd_bits_saved; |
| 817 | extern std::size_t total_lehmer_gcd_cycles; |
| 818 | |
| 819 | ++total_lehmer_gcd_calls; |
| 820 | total_lehmer_gcd_bits_saved += lu - eval_msb(U); |
| 821 | total_lehmer_gcd_cycles += i; |
| 822 | #endif |
| 823 | if (lu < 2048) |
| 824 | { |
| 825 | // |
| 826 | // Since we have stripped all common powers of 2 from U and V at the start |
| 827 | // if either are even at this point, we can remove stray powers of 2 now. |
| 828 | // Note that it is not possible for *both* U and V to be even at this point. |
| 829 | // |
| 830 | // This has an adverse effect on performance for high bit counts, but has |
| 831 | // a significant positive effect for smaller counts. |
| 832 | // |
| 833 | if ((U.limbs()[0] & 1u) == 0) |
| 834 | { |
| 835 | eval_right_shift(U, eval_lsb(U)); |
| 836 | if (U.compare(V) < 0) |
| 837 | U.swap(V); |
| 838 | } |
| 839 | else if ((V.limbs()[0] & 1u) == 0) |
| 840 | { |
| 841 | eval_right_shift(V, eval_lsb(V)); |
| 842 | } |
| 843 | } |
| 844 | storage.deallocate(ts * 3); |
| 845 | } |
| 846 | |
| 847 | #else |
| 848 | // |
| 849 | // This branch is taken when double_limb_type is a synthetic type with no native hardware support. |
| 850 | // For example __int128. The assumption is that add/subtract/multiply of double_limb_type are efficient, |
| 851 | // but that division is very slow. |
| 852 | // |
| 853 | // We begin with a specialized routine for division. |
| 854 | // We know that most of the time this is called the result will be 1. |
| 855 | // For small limb counts, this almost doubles the performance of Lehmer's routine! |
| 856 | // |
| 857 | BOOST_FORCEINLINE void divide_subtract(double_limb_type& q, double_limb_type& u, const double_limb_type& v) |
| 858 | { |
| 859 | BOOST_MP_ASSERT(q == 1); // precondition on entry. |
| 860 | u -= v; |
| 861 | while (u >= v) |
| 862 | { |
| 863 | u -= v; |
| 864 | if (++q > 30) |
| 865 | { |
| 866 | double_limb_type t = u / v; |
| 867 | u -= t * v; |
| 868 | q += t; |
| 869 | } |
| 870 | } |
| 871 | } |
| 872 | |
| 873 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class Storage> |
| 874 | void eval_gcd_lehmer(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& U, cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& V, std::size_t lu, Storage& storage) |
| 875 | { |
| 876 | // |
| 877 | // Extract the leading 2*bits_per_limb bits from U and V: |
| 878 | // |
| 879 | std::size_t h = lu % bits_per_limb; |
| 880 | double_limb_type u, v; |
| 881 | if (h) |
| 882 | { |
| 883 | u = (static_cast<double_limb_type>((U.limbs()[U.size() - 1])) << bits_per_limb) | U.limbs()[U.size() - 2]; |
| 884 | v = (static_cast<double_limb_type>((V.size() < U.size() ? 0 : V.limbs()[V.size() - 1])) << bits_per_limb) | V.limbs()[U.size() - 2]; |
| 885 | u <<= bits_per_limb - h; |
| 886 | u |= U.limbs()[U.size() - 3] >> h; |
| 887 | v <<= bits_per_limb - h; |
| 888 | v |= V.limbs()[U.size() - 3] >> h; |
| 889 | } |
| 890 | else |
| 891 | { |
| 892 | u = (static_cast<double_limb_type>(U.limbs()[U.size() - 1]) << bits_per_limb) | U.limbs()[U.size() - 2]; |
| 893 | v = (static_cast<double_limb_type>(V.limbs()[U.size() - 1]) << bits_per_limb) | V.limbs()[U.size() - 2]; |
| 894 | } |
| 895 | // |
| 896 | // Cosequences are stored as limb_types, we take care not to overflow these: |
| 897 | // |
| 898 | // x[i+0] is positive for even i. |
| 899 | // y[i+0] is positive for odd i. |
| 900 | // |
| 901 | // However we track only absolute values here: |
| 902 | // |
| 903 | limb_type x[3] = { 1, 0 }; |
| 904 | limb_type y[3] = { 0, 1 }; |
| 905 | std::size_t i = 0; |
| 906 | |
| 907 | #ifdef BOOST_MP_GCD_DEBUG |
| 908 | cpp_int UU, VV; |
| 909 | UU = U; |
| 910 | VV = V; |
| 911 | #endif |
| 912 | // |
| 913 | // We begine by running a single digit version of Lehmer's algorithm, we still have |
| 914 | // to track u and v at double precision, but this adds only a tiny performance penalty. |
| 915 | // What we gain is fast division, and fast termination testing. |
| 916 | // When you see static_cast<limb_type>(u >> bits_per_limb) here, this is really just |
| 917 | // a direct access to the upper bits_per_limb of the double limb type. For __int128 |
| 918 | // this is simple a load of the upper 64 bits and the "shift" is optimised away. |
| 919 | // |
| 920 | double_limb_type old_u, old_v; |
| 921 | while (true) |
| 922 | { |
| 923 | limb_type q = static_cast<limb_type>(u >> bits_per_limb) / static_cast<limb_type>(v >> bits_per_limb); |
| 924 | x[2] = x[0] + q * x[1]; |
| 925 | y[2] = y[0] + q * y[1]; |
| 926 | double_limb_type tu = u; |
| 927 | old_u = u; |
| 928 | old_v = v; |
| 929 | u = v; |
| 930 | double_limb_type t = q * v; |
| 931 | if (tu < t) |
| 932 | { |
| 933 | ++i; |
| 934 | break; |
| 935 | } |
| 936 | v = tu - t; |
| 937 | ++i; |
| 938 | BOOST_MP_ASSERT((u <= v) || (t / q == old_v)); |
| 939 | if (u <= v) |
| 940 | { |
| 941 | // We've gone terribly wrong, probably numeric overflow: |
| 942 | break; |
| 943 | } |
| 944 | if ((i & 1u) == 0) |
| 945 | { |
| 946 | if ((static_cast<limb_type>(v >> bits_per_limb) < x[2]) || ((static_cast<limb_type>(u >> bits_per_limb) - static_cast<limb_type>(v >> bits_per_limb)) < (y[2] + y[1]))) |
| 947 | break; |
| 948 | } |
| 949 | else |
| 950 | { |
| 951 | if ((static_cast<limb_type>(v >> bits_per_limb) < y[2]) || ((static_cast<limb_type>(u >> bits_per_limb) - static_cast<limb_type>(v >> bits_per_limb)) < (x[2] + x[1]))) |
| 952 | break; |
| 953 | } |
| 954 | #ifdef BOOST_MP_GCD_DEBUG |
| 955 | BOOST_MP_ASSERT(q == UU / VV); |
| 956 | UU %= VV; |
| 957 | UU.swap(VV); |
| 958 | #endif |
| 959 | x[0] = x[1]; |
| 960 | x[1] = x[2]; |
| 961 | y[0] = y[1]; |
| 962 | y[1] = y[2]; |
| 963 | } |
| 964 | // |
| 965 | // We get here when the single digit algorithm has gone wrong, back up i, u and v: |
| 966 | // |
| 967 | --i; |
| 968 | u = old_u; |
| 969 | v = old_v; |
| 970 | // |
| 971 | // Now run the full double-digit algorithm: |
| 972 | // |
| 973 | while (true) |
| 974 | { |
| 975 | double_limb_type q = 1u; |
| 976 | double_limb_type tt = u; |
| 977 | divide_subtract(q, u, v); |
| 978 | std::swap(a&: u, b&: v); |
| 979 | tt = y[0] + q * static_cast<double_limb_type>(y[1]); |
| 980 | // |
| 981 | // If calculation of y[2] would overflow a single limb, then we *must* terminate. |
| 982 | // Note that x[2] < y[2] so there is no need to check that as well: |
| 983 | // |
| 984 | if (tt >> bits_per_limb) |
| 985 | { |
| 986 | ++i; |
| 987 | break; |
| 988 | } |
| 989 | x[2] = static_cast<limb_type>(x[0] + static_cast<double_limb_type>(q * x[1])); |
| 990 | y[2] = static_cast<limb_type>(tt); |
| 991 | ++i; |
| 992 | if ((i & 1u) == 0) |
| 993 | { |
| 994 | BOOST_MP_ASSERT(u > v); |
| 995 | if ((v < x[2]) || ((u - v) < (static_cast<double_limb_type>(y[2]) + y[1]))) |
| 996 | break; |
| 997 | } |
| 998 | else |
| 999 | { |
| 1000 | BOOST_MP_ASSERT(u > v); |
| 1001 | if ((v < y[2]) || ((u - v) < (static_cast<double_limb_type>(x[2]) + x[1]))) |
| 1002 | break; |
| 1003 | } |
| 1004 | #ifdef BOOST_MP_GCD_DEBUG |
| 1005 | BOOST_MP_ASSERT(q == UU / VV); |
| 1006 | UU %= VV; |
| 1007 | UU.swap(VV); |
| 1008 | #endif |
| 1009 | x[0] = x[1]; |
| 1010 | x[1] = x[2]; |
| 1011 | y[0] = y[1]; |
| 1012 | y[1] = y[2]; |
| 1013 | } |
| 1014 | if (i == 1) |
| 1015 | { |
| 1016 | // No change to U and V we've stalled! |
| 1017 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t; |
| 1018 | eval_modulus(t, U, V); |
| 1019 | U.swap(V); |
| 1020 | V.swap(t); |
| 1021 | return; |
| 1022 | } |
| 1023 | // |
| 1024 | // Update U and V. |
| 1025 | // We have: |
| 1026 | // |
| 1027 | // U = x[0]U + y[0]V and |
| 1028 | // V = x[1]U + y[1]V. |
| 1029 | // |
| 1030 | // But since we track only absolute values of x and y |
| 1031 | // we have to take account of the implied signs and perform |
| 1032 | // the appropriate subtraction depending on the whether i is |
| 1033 | // even or odd: |
| 1034 | // |
| 1035 | std::size_t ts = U.size() + 1; |
| 1036 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t1(storage, ts), t2(storage, ts), t3(storage, ts); |
| 1037 | eval_multiply(t1, U, x[0]); |
| 1038 | eval_multiply(t2, V, y[0]); |
| 1039 | eval_multiply(t3, U, x[1]); |
| 1040 | if ((i & 1u) == 0) |
| 1041 | { |
| 1042 | if (x[0] == 0) |
| 1043 | U = t2; |
| 1044 | else |
| 1045 | { |
| 1046 | BOOST_MP_ASSERT(t2.compare(t1) >= 0); |
| 1047 | eval_subtract(U, t2, t1); |
| 1048 | BOOST_MP_ASSERT(U.sign() == false); |
| 1049 | } |
| 1050 | } |
| 1051 | else |
| 1052 | { |
| 1053 | BOOST_MP_ASSERT(t1.compare(t2) >= 0); |
| 1054 | eval_subtract(U, t1, t2); |
| 1055 | BOOST_MP_ASSERT(U.sign() == false); |
| 1056 | } |
| 1057 | eval_multiply(t2, V, y[1]); |
| 1058 | if (i & 1u) |
| 1059 | { |
| 1060 | if (x[1] == 0) |
| 1061 | V = t2; |
| 1062 | else |
| 1063 | { |
| 1064 | BOOST_MP_ASSERT(t2.compare(t3) >= 0); |
| 1065 | eval_subtract(V, t2, t3); |
| 1066 | BOOST_MP_ASSERT(V.sign() == false); |
| 1067 | } |
| 1068 | } |
| 1069 | else |
| 1070 | { |
| 1071 | BOOST_MP_ASSERT(t3.compare(t2) >= 0); |
| 1072 | eval_subtract(V, t3, t2); |
| 1073 | BOOST_MP_ASSERT(V.sign() == false); |
| 1074 | } |
| 1075 | BOOST_MP_ASSERT(U.compare(V) >= 0); |
| 1076 | BOOST_MP_ASSERT(lu > eval_msb(U)); |
| 1077 | #ifdef BOOST_MP_GCD_DEBUG |
| 1078 | |
| 1079 | BOOST_MP_ASSERT(UU == U); |
| 1080 | BOOST_MP_ASSERT(VV == V); |
| 1081 | |
| 1082 | extern std::size_t total_lehmer_gcd_calls; |
| 1083 | extern std::size_t total_lehmer_gcd_bits_saved; |
| 1084 | extern std::size_t total_lehmer_gcd_cycles; |
| 1085 | |
| 1086 | ++total_lehmer_gcd_calls; |
| 1087 | total_lehmer_gcd_bits_saved += lu - eval_msb(U); |
| 1088 | total_lehmer_gcd_cycles += i; |
| 1089 | #endif |
| 1090 | if (lu < 2048) |
| 1091 | { |
| 1092 | // |
| 1093 | // Since we have stripped all common powers of 2 from U and V at the start |
| 1094 | // if either are even at this point, we can remove stray powers of 2 now. |
| 1095 | // Note that it is not possible for *both* U and V to be even at this point. |
| 1096 | // |
| 1097 | // This has an adverse effect on performance for high bit counts, but has |
| 1098 | // a significant positive effect for smaller counts. |
| 1099 | // |
| 1100 | if ((U.limbs()[0] & 1u) == 0) |
| 1101 | { |
| 1102 | eval_right_shift(U, eval_lsb(U)); |
| 1103 | if (U.compare(V) < 0) |
| 1104 | U.swap(V); |
| 1105 | } |
| 1106 | else if ((V.limbs()[0] & 1u) == 0) |
| 1107 | { |
| 1108 | eval_right_shift(V, eval_lsb(V)); |
| 1109 | } |
| 1110 | } |
| 1111 | storage.deallocate(ts * 3); |
| 1112 | } |
| 1113 | |
| 1114 | #endif |
| 1115 | |
| 1116 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1117 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 1118 | eval_gcd( |
| 1119 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, |
| 1120 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, |
| 1121 | const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& b) |
| 1122 | { |
| 1123 | using default_ops::eval_get_sign; |
| 1124 | using default_ops::eval_is_zero; |
| 1125 | using default_ops::eval_lsb; |
| 1126 | |
| 1127 | if (a.size() == 1) |
| 1128 | { |
| 1129 | eval_gcd(result, b, *a.limbs()); |
| 1130 | return; |
| 1131 | } |
| 1132 | if (b.size() == 1) |
| 1133 | { |
| 1134 | eval_gcd(result, a, *b.limbs()); |
| 1135 | return; |
| 1136 | } |
| 1137 | std::size_t temp_size = (std::max)(a.size(), b.size()) + 1; |
| 1138 | typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::scoped_shared_storage storage(a, temp_size * 6); |
| 1139 | |
| 1140 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> U(storage, temp_size); |
| 1141 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> V(storage, temp_size); |
| 1142 | cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t(storage, temp_size); |
| 1143 | U = a; |
| 1144 | V = b; |
| 1145 | |
| 1146 | int s = eval_get_sign(U); |
| 1147 | |
| 1148 | /* GCD(0,x) := x */ |
| 1149 | if (s < 0) |
| 1150 | { |
| 1151 | U.negate(); |
| 1152 | } |
| 1153 | else if (s == 0) |
| 1154 | { |
| 1155 | result = V; |
| 1156 | return; |
| 1157 | } |
| 1158 | s = eval_get_sign(V); |
| 1159 | if (s < 0) |
| 1160 | { |
| 1161 | V.negate(); |
| 1162 | } |
| 1163 | else if (s == 0) |
| 1164 | { |
| 1165 | result = U; |
| 1166 | return; |
| 1167 | } |
| 1168 | // |
| 1169 | // Remove common factors of 2: |
| 1170 | // |
| 1171 | std::size_t us = eval_lsb(U); |
| 1172 | std::size_t vs = eval_lsb(V); |
| 1173 | std::size_t shift = (std::min)(a: us, b: vs); |
| 1174 | if (us) |
| 1175 | eval_right_shift(U, us); |
| 1176 | if (vs) |
| 1177 | eval_right_shift(V, vs); |
| 1178 | |
| 1179 | if (U.compare(V) < 0) |
| 1180 | U.swap(V); |
| 1181 | |
| 1182 | while (!eval_is_zero(V)) |
| 1183 | { |
| 1184 | if (U.size() <= 2) |
| 1185 | { |
| 1186 | // |
| 1187 | // Special case: if V has no more than 2 limbs |
| 1188 | // then we can reduce U and V to a pair of integers and perform |
| 1189 | // direct integer gcd: |
| 1190 | // |
| 1191 | if (U.size() == 1) |
| 1192 | U = eval_gcd(*V.limbs(), *U.limbs()); |
| 1193 | else |
| 1194 | { |
| 1195 | double_limb_type i = U.limbs()[0] | (static_cast<double_limb_type>(U.limbs()[1]) << sizeof(limb_type) * CHAR_BIT); |
| 1196 | double_limb_type j = (V.size() == 1) ? *V.limbs() : V.limbs()[0] | (static_cast<double_limb_type>(V.limbs()[1]) << sizeof(limb_type) * CHAR_BIT); |
| 1197 | U = eval_gcd(u: i, v: j); |
| 1198 | } |
| 1199 | break; |
| 1200 | } |
| 1201 | std::size_t lu = eval_msb(U) + 1; |
| 1202 | std::size_t lv = eval_msb(V) + 1; |
| 1203 | #ifndef BOOST_MP_NO_CONSTEXPR_DETECTION |
| 1204 | if (!BOOST_MP_IS_CONST_EVALUATED(lu) && (lu - lv <= bits_per_limb / 2)) |
| 1205 | #else |
| 1206 | if (lu - lv <= bits_per_limb / 2) |
| 1207 | #endif |
| 1208 | { |
| 1209 | eval_gcd_lehmer(U, V, lu, storage); |
| 1210 | } |
| 1211 | else |
| 1212 | { |
| 1213 | eval_modulus(t, U, V); |
| 1214 | U.swap(V); |
| 1215 | V.swap(t); |
| 1216 | } |
| 1217 | } |
| 1218 | result = U; |
| 1219 | if (shift) |
| 1220 | eval_left_shift(result, shift); |
| 1221 | } |
| 1222 | // |
| 1223 | // Now again for trivial backends: |
| 1224 | // |
| 1225 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1226 | BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type |
| 1227 | eval_gcd(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& b) noexcept |
| 1228 | { |
| 1229 | *result.limbs() = boost::multiprecision::detail::constexpr_gcd(*a.limbs(), *b.limbs()); |
| 1230 | result.sign(false); |
| 1231 | } |
| 1232 | // This one is only enabled for unchecked cpp_int's, for checked int's we need the checking in the default version: |
| 1233 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1234 | BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && (Checked1 == unchecked)>::type |
| 1235 | eval_lcm(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& b) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value)) |
| 1236 | { |
| 1237 | *result.limbs() = boost::multiprecision::detail::constexpr_lcm(*a.limbs(), *b.limbs()); |
| 1238 | result.normalize(); // result may overflow the specified number of bits |
| 1239 | result.sign(false); |
| 1240 | } |
| 1241 | |
| 1242 | inline void conversion_overflow(const std::integral_constant<int, checked>&) |
| 1243 | { |
| 1244 | BOOST_MP_THROW_EXCEPTION(std::overflow_error("Overflow in conversion to narrower type" )); |
| 1245 | } |
| 1246 | inline BOOST_MP_CXX14_CONSTEXPR void conversion_overflow(const std::integral_constant<int, unchecked>&) {} |
| 1247 | |
| 1248 | #if defined(__clang__) && defined(__MINGW32__) |
| 1249 | // |
| 1250 | // clang-11 on Mingw segfaults on conversion of __int128 -> float. |
| 1251 | // See: https://bugs.llvm.org/show_bug.cgi?id=48941 |
| 1252 | // These workarounds pass everything through an intermediate uint64_t. |
| 1253 | // |
| 1254 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1255 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if< |
| 1256 | is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_signed_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && std::is_same<typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::local_limb_type, double_limb_type>::value>::type |
| 1257 | eval_convert_to(float* result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) |
| 1258 | { |
| 1259 | float f = static_cast<std::uint64_t>((*val.limbs()) >> 64); |
| 1260 | *result = std::ldexp(f, 64); |
| 1261 | *result += static_cast<std::uint64_t>((*val.limbs())); |
| 1262 | if(val.sign()) |
| 1263 | *result = -*result; |
| 1264 | } |
| 1265 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1266 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if< |
| 1267 | is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_signed_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && std::is_same<typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::local_limb_type, double_limb_type>::value>::type |
| 1268 | eval_convert_to(double* result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) |
| 1269 | { |
| 1270 | float f = static_cast<std::uint64_t>((*val.limbs()) >> 64); |
| 1271 | *result = std::ldexp(f, 64); |
| 1272 | *result += static_cast<std::uint64_t>((*val.limbs())); |
| 1273 | if(val.sign()) |
| 1274 | *result = -*result; |
| 1275 | } |
| 1276 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1277 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if< |
| 1278 | is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_signed_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && std::is_same<typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::local_limb_type, double_limb_type>::value>::type |
| 1279 | eval_convert_to(long double* result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) |
| 1280 | { |
| 1281 | float f = static_cast<std::uint64_t>((*val.limbs()) >> 64); |
| 1282 | *result = std::ldexp(f, 64); |
| 1283 | *result += static_cast<std::uint64_t>((*val.limbs())); |
| 1284 | if(val.sign()) |
| 1285 | *result = -*result; |
| 1286 | } |
| 1287 | #endif |
| 1288 | |
| 1289 | template <class R, std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1290 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if< |
| 1291 | is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_signed_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && std::is_convertible<typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::local_limb_type, R>::value>::type |
| 1292 | eval_convert_to(R* result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) |
| 1293 | { |
| 1294 | BOOST_IF_CONSTEXPR(std::numeric_limits<R>::is_specialized) |
| 1295 | { |
| 1296 | using common_type = typename std::common_type<R, typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::local_limb_type>::type; |
| 1297 | |
| 1298 | if (static_cast<common_type>(*val.limbs()) > static_cast<common_type>((std::numeric_limits<R>::max)())) |
| 1299 | { |
| 1300 | if (val.isneg()) |
| 1301 | { |
| 1302 | check_is_negative(std::integral_constant < bool, (boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value) || (number_category<R>::value == number_kind_floating_point) > ()); |
| 1303 | if (static_cast<common_type>(*val.limbs()) > -static_cast<common_type>((std::numeric_limits<R>::min)())) |
| 1304 | conversion_overflow(typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::checked_type()); |
| 1305 | *result = (std::numeric_limits<R>::min)(); |
| 1306 | } |
| 1307 | else |
| 1308 | { |
| 1309 | conversion_overflow(typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::checked_type()); |
| 1310 | *result = boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value ? (std::numeric_limits<R>::max)() : static_cast<R>(*val.limbs()); |
| 1311 | } |
| 1312 | } |
| 1313 | else |
| 1314 | { |
| 1315 | *result = static_cast<R>(*val.limbs()); |
| 1316 | if (val.isneg()) |
| 1317 | { |
| 1318 | check_is_negative(std::integral_constant < bool, (boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value) || (number_category<R>::value == number_kind_floating_point) > ()); |
| 1319 | *result = negate_integer(*result, std::integral_constant < bool, is_signed_number<R>::value || (number_category<R>::value == number_kind_floating_point) > ()); |
| 1320 | } |
| 1321 | } |
| 1322 | } |
| 1323 | else |
| 1324 | { |
| 1325 | *result = static_cast<R>(*val.limbs()); |
| 1326 | if (val.isneg()) |
| 1327 | { |
| 1328 | check_is_negative(std::integral_constant<bool, (boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value) || (number_category<R>::value == number_kind_floating_point) > ()); |
| 1329 | *result = negate_integer(*result, std::integral_constant<bool, is_signed_number<R>::value || (number_category<R>::value == number_kind_floating_point) > ()); |
| 1330 | } |
| 1331 | } |
| 1332 | } |
| 1333 | |
| 1334 | template <class R, std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1335 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if< |
| 1336 | is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_unsigned_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && std::is_convertible<typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::local_limb_type, R>::value>::type |
| 1337 | eval_convert_to(R* result, const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) |
| 1338 | { |
| 1339 | BOOST_IF_CONSTEXPR(std::numeric_limits<R>::is_specialized) |
| 1340 | { |
| 1341 | using common_type = typename std::common_type<R, typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::local_limb_type>::type; |
| 1342 | |
| 1343 | if(static_cast<common_type>(*val.limbs()) > static_cast<common_type>((std::numeric_limits<R>::max)())) |
| 1344 | { |
| 1345 | conversion_overflow(typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::checked_type()); |
| 1346 | *result = boost::multiprecision::detail::is_signed<R>::value && boost::multiprecision::detail::is_integral<R>::value ? (std::numeric_limits<R>::max)() : static_cast<R>(*val.limbs()); |
| 1347 | } |
| 1348 | else |
| 1349 | *result = static_cast<R>(*val.limbs()); |
| 1350 | } |
| 1351 | else |
| 1352 | *result = static_cast<R>(*val.limbs()); |
| 1353 | } |
| 1354 | |
| 1355 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1356 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, std::size_t>::type |
| 1357 | eval_lsb(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a) |
| 1358 | { |
| 1359 | using default_ops::eval_get_sign; |
| 1360 | if (eval_get_sign(a) == 0) |
| 1361 | { |
| 1362 | BOOST_MP_THROW_EXCEPTION(std::domain_error("No bits were set in the operand." )); |
| 1363 | } |
| 1364 | if (a.sign()) |
| 1365 | { |
| 1366 | BOOST_MP_THROW_EXCEPTION(std::domain_error("Testing individual bits in negative values is not supported - results are undefined." )); |
| 1367 | } |
| 1368 | // |
| 1369 | // Find the index of the least significant bit within that limb: |
| 1370 | // |
| 1371 | return boost::multiprecision::detail::find_lsb(*a.limbs()); |
| 1372 | } |
| 1373 | |
| 1374 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1375 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, std::size_t>::type |
| 1376 | eval_msb_imp(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a) |
| 1377 | { |
| 1378 | // |
| 1379 | // Find the index of the least significant bit within that limb: |
| 1380 | // |
| 1381 | return boost::multiprecision::detail::find_msb(*a.limbs()); |
| 1382 | } |
| 1383 | |
| 1384 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1385 | inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value, std::size_t>::type |
| 1386 | eval_msb(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a) |
| 1387 | { |
| 1388 | using default_ops::eval_get_sign; |
| 1389 | if (eval_get_sign(a) == 0) |
| 1390 | { |
| 1391 | BOOST_MP_THROW_EXCEPTION(std::domain_error("No bits were set in the operand." )); |
| 1392 | } |
| 1393 | if (a.sign()) |
| 1394 | { |
| 1395 | BOOST_MP_THROW_EXCEPTION(std::domain_error("Testing individual bits in negative values is not supported - results are undefined." )); |
| 1396 | } |
| 1397 | return eval_msb_imp(a); |
| 1398 | } |
| 1399 | |
| 1400 | template <std::size_t MinBits1, std::size_t MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1> |
| 1401 | inline BOOST_MP_CXX14_CONSTEXPR std::size_t hash_value(const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& val) noexcept |
| 1402 | { |
| 1403 | std::size_t result = 0; |
| 1404 | for (std::size_t i = 0; i < val.size(); ++i) |
| 1405 | { |
| 1406 | boost::multiprecision::detail::hash_combine(result, val.limbs()[i]); |
| 1407 | } |
| 1408 | boost::multiprecision::detail::hash_combine(result, val.sign()); |
| 1409 | return result; |
| 1410 | } |
| 1411 | |
| 1412 | #ifdef BOOST_MSVC |
| 1413 | #pragma warning(pop) |
| 1414 | #endif |
| 1415 | |
| 1416 | } // Namespace backends |
| 1417 | |
| 1418 | namespace detail { |
| 1419 | |
| 1420 | #ifndef BOOST_MP_STANDALONE |
| 1421 | template <typename T> |
| 1422 | inline BOOST_CXX14_CONSTEXPR T constexpr_gcd(T a, T b) noexcept |
| 1423 | { |
| 1424 | return boost::integer::gcd(a, b); |
| 1425 | } |
| 1426 | |
| 1427 | template <typename T> |
| 1428 | inline BOOST_CXX14_CONSTEXPR T constexpr_lcm(T a, T b) noexcept |
| 1429 | { |
| 1430 | return boost::integer::lcm(a, b); |
| 1431 | } |
| 1432 | |
| 1433 | #else |
| 1434 | |
| 1435 | template <typename T> |
| 1436 | inline BOOST_CXX14_CONSTEXPR T constexpr_gcd(T a, T b) noexcept |
| 1437 | { |
| 1438 | return boost::multiprecision::backends::eval_gcd(a, b); |
| 1439 | } |
| 1440 | |
| 1441 | template <typename T> |
| 1442 | inline BOOST_CXX14_CONSTEXPR T constexpr_lcm(T a, T b) noexcept |
| 1443 | { |
| 1444 | const T ab_gcd = boost::multiprecision::detail::constexpr_gcd(a, b); |
| 1445 | return (a * b) / ab_gcd; |
| 1446 | } |
| 1447 | |
| 1448 | #endif // BOOST_MP_STANDALONE |
| 1449 | |
| 1450 | } |
| 1451 | |
| 1452 | }} // Namespace boost::multiprecision |
| 1453 | |
| 1454 | #endif |
| 1455 | |