annotate DEPENDENCIES/generic/include/boost/multiprecision/detail/generic_interconvert.hpp @ 125:34e428693f5d vext

Vext -> Repoint
author Chris Cannam
date Thu, 14 Jun 2018 11:15:39 +0100
parents c530137014c0
children
rev   line source
Chris@16 1 ///////////////////////////////////////////////////////////////////////////////
Chris@16 2 // Copyright 2011 John Maddock. Distributed under the Boost
Chris@16 3 // Software License, Version 1.0. (See accompanying file
Chris@16 4 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
Chris@16 5
Chris@16 6 #ifndef BOOST_MP_GENERIC_INTERCONVERT_HPP
Chris@16 7 #define BOOST_MP_GENERIC_INTERCONVERT_HPP
Chris@16 8
Chris@16 9 #include <boost/multiprecision/detail/default_ops.hpp>
Chris@16 10
Chris@101 11 #ifdef BOOST_MSVC
Chris@101 12 #pragma warning(push)
Chris@101 13 #pragma warning(disable:4127)
Chris@101 14 #endif
Chris@101 15
Chris@16 16 namespace boost{ namespace multiprecision{ namespace detail{
Chris@16 17
Chris@16 18 template <class To, class From>
Chris@101 19 inline To do_cast(const From & from)
Chris@101 20 {
Chris@101 21 return static_cast<To>(from);
Chris@101 22 }
Chris@101 23 template <class To, class B, ::boost::multiprecision::expression_template_option et>
Chris@101 24 inline To do_cast(const number<B, et>& from)
Chris@101 25 {
Chris@101 26 return from.template convert_to<To>();
Chris@101 27 }
Chris@101 28
Chris@101 29 template <class To, class From>
Chris@16 30 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_floating_point>& /*to_type*/, const mpl::int_<number_kind_integer>& /*from_type*/)
Chris@16 31 {
Chris@16 32 using default_ops::eval_get_sign;
Chris@16 33 using default_ops::eval_bitwise_and;
Chris@16 34 using default_ops::eval_convert_to;
Chris@16 35 using default_ops::eval_right_shift;
Chris@16 36 using default_ops::eval_ldexp;
Chris@16 37 using default_ops::eval_add;
Chris@16 38 // smallest unsigned type handled natively by "From" is likely to be it's limb_type:
Chris@16 39 typedef typename canonical<unsigned char, From>::type limb_type;
Chris@16 40 // get the corresponding type that we can assign to "To":
Chris@16 41 typedef typename canonical<limb_type, To>::type to_type;
Chris@16 42 From t(from);
Chris@16 43 bool is_neg = eval_get_sign(t) < 0;
Chris@16 44 if(is_neg)
Chris@16 45 t.negate();
Chris@16 46 // Pick off the first limb:
Chris@16 47 limb_type limb;
Chris@16 48 limb_type mask = ~static_cast<limb_type>(0);
Chris@16 49 From fl;
Chris@16 50 eval_bitwise_and(fl, t, mask);
Chris@16 51 eval_convert_to(&limb, fl);
Chris@16 52 to = static_cast<to_type>(limb);
Chris@16 53 eval_right_shift(t, std::numeric_limits<limb_type>::digits);
Chris@16 54 //
Chris@16 55 // Then keep picking off more limbs until "t" is zero:
Chris@16 56 //
Chris@16 57 To l;
Chris@16 58 unsigned shift = std::numeric_limits<limb_type>::digits;
Chris@16 59 while(!eval_is_zero(t))
Chris@16 60 {
Chris@16 61 eval_bitwise_and(fl, t, mask);
Chris@16 62 eval_convert_to(&limb, fl);
Chris@16 63 l = static_cast<to_type>(limb);
Chris@16 64 eval_right_shift(t, std::numeric_limits<limb_type>::digits);
Chris@16 65 eval_ldexp(l, l, shift);
Chris@16 66 eval_add(to, l);
Chris@16 67 shift += std::numeric_limits<limb_type>::digits;
Chris@16 68 }
Chris@16 69 //
Chris@16 70 // Finish off by setting the sign:
Chris@16 71 //
Chris@16 72 if(is_neg)
Chris@16 73 to.negate();
Chris@16 74 }
Chris@16 75
Chris@16 76 template <class To, class From>
Chris@16 77 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_integer>& /*to_type*/, const mpl::int_<number_kind_integer>& /*from_type*/)
Chris@16 78 {
Chris@16 79 using default_ops::eval_get_sign;
Chris@16 80 using default_ops::eval_bitwise_and;
Chris@16 81 using default_ops::eval_convert_to;
Chris@16 82 using default_ops::eval_right_shift;
Chris@16 83 using default_ops::eval_left_shift;
Chris@16 84 using default_ops::eval_bitwise_or;
Chris@16 85 using default_ops::eval_is_zero;
Chris@16 86 // smallest unsigned type handled natively by "From" is likely to be it's limb_type:
Chris@16 87 typedef typename canonical<unsigned char, From>::type limb_type;
Chris@16 88 // get the corresponding type that we can assign to "To":
Chris@16 89 typedef typename canonical<limb_type, To>::type to_type;
Chris@16 90 From t(from);
Chris@16 91 bool is_neg = eval_get_sign(t) < 0;
Chris@16 92 if(is_neg)
Chris@16 93 t.negate();
Chris@16 94 // Pick off the first limb:
Chris@16 95 limb_type limb;
Chris@16 96 limb_type mask = static_cast<limb_type>(~static_cast<limb_type>(0));
Chris@16 97 From fl;
Chris@16 98 eval_bitwise_and(fl, t, mask);
Chris@16 99 eval_convert_to(&limb, fl);
Chris@16 100 to = static_cast<to_type>(limb);
Chris@16 101 eval_right_shift(t, std::numeric_limits<limb_type>::digits);
Chris@16 102 //
Chris@16 103 // Then keep picking off more limbs until "t" is zero:
Chris@16 104 //
Chris@16 105 To l;
Chris@16 106 unsigned shift = std::numeric_limits<limb_type>::digits;
Chris@16 107 while(!eval_is_zero(t))
Chris@16 108 {
Chris@16 109 eval_bitwise_and(fl, t, mask);
Chris@16 110 eval_convert_to(&limb, fl);
Chris@16 111 l = static_cast<to_type>(limb);
Chris@16 112 eval_right_shift(t, std::numeric_limits<limb_type>::digits);
Chris@16 113 eval_left_shift(l, shift);
Chris@16 114 eval_bitwise_or(to, l);
Chris@16 115 shift += std::numeric_limits<limb_type>::digits;
Chris@16 116 }
Chris@16 117 //
Chris@16 118 // Finish off by setting the sign:
Chris@16 119 //
Chris@16 120 if(is_neg)
Chris@16 121 to.negate();
Chris@16 122 }
Chris@16 123
Chris@16 124 template <class To, class From>
Chris@16 125 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_floating_point>& /*to_type*/, const mpl::int_<number_kind_floating_point>& /*from_type*/)
Chris@16 126 {
Chris@16 127 #ifdef BOOST_MSVC
Chris@16 128 #pragma warning(push)
Chris@16 129 #pragma warning(disable:4127)
Chris@16 130 #endif
Chris@16 131 //
Chris@16 132 // The code here only works when the radix of "From" is 2, we could try shifting by other
Chris@16 133 // radixes but it would complicate things.... use a string conversion when the radix is other
Chris@16 134 // than 2:
Chris@16 135 //
Chris@16 136 if(std::numeric_limits<number<From> >::radix != 2)
Chris@16 137 {
Chris@16 138 to = from.str(0, std::ios_base::fmtflags()).c_str();
Chris@16 139 return;
Chris@16 140 }
Chris@16 141
Chris@16 142
Chris@16 143 typedef typename canonical<unsigned char, To>::type ui_type;
Chris@16 144
Chris@16 145 using default_ops::eval_fpclassify;
Chris@16 146 using default_ops::eval_add;
Chris@16 147 using default_ops::eval_subtract;
Chris@16 148 using default_ops::eval_convert_to;
Chris@16 149
Chris@16 150 //
Chris@16 151 // First classify the input, then handle the special cases:
Chris@16 152 //
Chris@16 153 int c = eval_fpclassify(from);
Chris@16 154
Chris@101 155 if(c == (int)FP_ZERO)
Chris@16 156 {
Chris@16 157 to = ui_type(0);
Chris@16 158 return;
Chris@16 159 }
Chris@101 160 else if(c == (int)FP_NAN)
Chris@16 161 {
Chris@101 162 to = static_cast<const char*>("nan");
Chris@16 163 return;
Chris@16 164 }
Chris@101 165 else if(c == (int)FP_INFINITE)
Chris@16 166 {
Chris@101 167 to = static_cast<const char*>("inf");
Chris@16 168 if(eval_get_sign(from) < 0)
Chris@16 169 to.negate();
Chris@16 170 return;
Chris@16 171 }
Chris@16 172
Chris@16 173 typename From::exponent_type e;
Chris@16 174 From f, term;
Chris@16 175 to = ui_type(0);
Chris@16 176
Chris@16 177 eval_frexp(f, from, &e);
Chris@16 178
Chris@16 179 static const int shift = std::numeric_limits<boost::intmax_t>::digits - 1;
Chris@16 180
Chris@16 181 while(!eval_is_zero(f))
Chris@16 182 {
Chris@16 183 // extract int sized bits from f:
Chris@16 184 eval_ldexp(f, f, shift);
Chris@16 185 eval_floor(term, f);
Chris@16 186 e -= shift;
Chris@16 187 eval_ldexp(to, to, shift);
Chris@16 188 typename boost::multiprecision::detail::canonical<boost::intmax_t, To>::type ll;
Chris@16 189 eval_convert_to(&ll, term);
Chris@16 190 eval_add(to, ll);
Chris@16 191 eval_subtract(f, term);
Chris@16 192 }
Chris@16 193 typedef typename To::exponent_type to_exponent;
Chris@16 194 if((e > (std::numeric_limits<to_exponent>::max)()) || (e < (std::numeric_limits<to_exponent>::min)()))
Chris@16 195 {
Chris@101 196 to = static_cast<const char*>("inf");
Chris@16 197 if(eval_get_sign(from) < 0)
Chris@16 198 to.negate();
Chris@16 199 return;
Chris@16 200 }
Chris@16 201 eval_ldexp(to, to, static_cast<to_exponent>(e));
Chris@16 202 #ifdef BOOST_MSVC
Chris@16 203 #pragma warning(pop)
Chris@16 204 #endif
Chris@16 205 }
Chris@16 206
Chris@16 207 template <class To, class From>
Chris@16 208 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_rational>& /*to_type*/, const mpl::int_<number_kind_rational>& /*from_type*/)
Chris@16 209 {
Chris@16 210 typedef typename component_type<number<To> >::type to_component_type;
Chris@16 211
Chris@16 212 number<From> t(from);
Chris@16 213 to_component_type n(numerator(t)), d(denominator(t));
Chris@16 214 using default_ops::assign_components;
Chris@16 215 assign_components(to, n.backend(), d.backend());
Chris@16 216 }
Chris@16 217
Chris@16 218 template <class To, class From>
Chris@16 219 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_rational>& /*to_type*/, const mpl::int_<number_kind_integer>& /*from_type*/)
Chris@16 220 {
Chris@16 221 typedef typename component_type<number<To> >::type to_component_type;
Chris@16 222
Chris@16 223 number<From> t(from);
Chris@16 224 to_component_type n(t), d(1);
Chris@16 225 using default_ops::assign_components;
Chris@16 226 assign_components(to, n.backend(), d.backend());
Chris@16 227 }
Chris@16 228
Chris@101 229 template <class R, class LargeInteger>
Chris@101 230 R safe_convert_to_float(const LargeInteger& i)
Chris@101 231 {
Chris@101 232 using std::ldexp;
Chris@101 233 if(!i)
Chris@101 234 return R(0);
Chris@101 235 if(std::numeric_limits<R>::is_specialized && std::numeric_limits<R>::max_exponent)
Chris@101 236 {
Chris@101 237 LargeInteger val(i);
Chris@101 238 if(val.sign() < 0)
Chris@101 239 val = -val;
Chris@101 240 unsigned mb = msb(val);
Chris@101 241 if(mb >= std::numeric_limits<R>::max_exponent)
Chris@101 242 {
Chris@101 243 int scale_factor = (int)mb + 1 - std::numeric_limits<R>::max_exponent;
Chris@101 244 BOOST_ASSERT(scale_factor >= 1);
Chris@101 245 val >>= scale_factor;
Chris@101 246 R result = val.template convert_to<R>();
Chris@101 247 if(std::numeric_limits<R>::digits == 0 || std::numeric_limits<R>::digits >= std::numeric_limits<R>::max_exponent)
Chris@101 248 {
Chris@101 249 //
Chris@101 250 // Calculate and add on the remainder, only if there are more
Chris@101 251 // digits in the mantissa that the size of the exponent, in
Chris@101 252 // other words if we are dropping digits in the conversion
Chris@101 253 // otherwise:
Chris@101 254 //
Chris@101 255 LargeInteger remainder(i);
Chris@101 256 remainder &= (LargeInteger(1) << scale_factor) - 1;
Chris@101 257 result += ldexp(safe_convert_to_float<R>(remainder), -scale_factor);
Chris@101 258 }
Chris@101 259 return i.sign() < 0 ? static_cast<R>(-result) : result;
Chris@101 260 }
Chris@101 261 }
Chris@101 262 return i.template convert_to<R>();
Chris@101 263 }
Chris@101 264
Chris@101 265 template <class To, class Integer>
Chris@101 266 inline typename disable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
Chris@101 267 generic_convert_rational_to_float_imp(To& result, const Integer& n, const Integer& d, const mpl::true_&)
Chris@101 268 {
Chris@101 269 //
Chris@101 270 // If we get here, then there's something about one type or the other
Chris@101 271 // that prevents an exactly rounded result from being calculated
Chris@101 272 // (or at least it's not clear how to implement such a thing).
Chris@101 273 //
Chris@101 274 using default_ops::eval_divide;
Chris@101 275 number<To> fn(safe_convert_to_float<number<To> >(n)), fd(safe_convert_to_float<number<To> >(d));
Chris@101 276 eval_divide(result, fn.backend(), fd.backend());
Chris@101 277 }
Chris@101 278 template <class To, class Integer>
Chris@101 279 inline typename enable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
Chris@101 280 generic_convert_rational_to_float_imp(To& result, const Integer& n, const Integer& d, const mpl::true_&)
Chris@101 281 {
Chris@101 282 //
Chris@101 283 // If we get here, then there's something about one type or the other
Chris@101 284 // that prevents an exactly rounded result from being calculated
Chris@101 285 // (or at least it's not clear how to implement such a thing).
Chris@101 286 //
Chris@101 287 To fd(safe_convert_to_float<To>(d));
Chris@101 288 result = safe_convert_to_float<To>(n);
Chris@101 289 result /= fd;
Chris@101 290 }
Chris@101 291
Chris@101 292 template <class To, class Integer>
Chris@101 293 typename enable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
Chris@101 294 generic_convert_rational_to_float_imp(To& result, Integer& num, Integer& denom, const mpl::false_&)
Chris@101 295 {
Chris@101 296 //
Chris@101 297 // If we get here, then the precision of type To is known, and the integer type is unbounded
Chris@101 298 // so we can use integer division plus manipulation of the remainder to get an exactly
Chris@101 299 // rounded result.
Chris@101 300 //
Chris@101 301 if(num == 0)
Chris@101 302 {
Chris@101 303 result = 0;
Chris@101 304 return;
Chris@101 305 }
Chris@101 306 bool s = false;
Chris@101 307 if(num < 0)
Chris@101 308 {
Chris@101 309 s = true;
Chris@101 310 num = -num;
Chris@101 311 }
Chris@101 312 int denom_bits = msb(denom);
Chris@101 313 int shift = std::numeric_limits<To>::digits + denom_bits - msb(num) + 1;
Chris@101 314 if(shift > 0)
Chris@101 315 num <<= shift;
Chris@101 316 else if(shift < 0)
Chris@101 317 denom <<= std::abs(shift);
Chris@101 318 Integer q, r;
Chris@101 319 divide_qr(num, denom, q, r);
Chris@101 320 int q_bits = msb(q);
Chris@101 321 if(q_bits == std::numeric_limits<To>::digits)
Chris@101 322 {
Chris@101 323 //
Chris@101 324 // Round up if 2 * r > denom:
Chris@101 325 //
Chris@101 326 r <<= 1;
Chris@101 327 int c = r.compare(denom);
Chris@101 328 if(c > 0)
Chris@101 329 ++q;
Chris@101 330 else if((c == 0) && (q & 1u))
Chris@101 331 {
Chris@101 332 ++q;
Chris@101 333 }
Chris@101 334 }
Chris@101 335 else
Chris@101 336 {
Chris@101 337 BOOST_ASSERT(q_bits == 1 + std::numeric_limits<To>::digits);
Chris@101 338 //
Chris@101 339 // We basically already have the rounding info:
Chris@101 340 //
Chris@101 341 if(q & 1u)
Chris@101 342 {
Chris@101 343 if(r || (q & 2u))
Chris@101 344 ++q;
Chris@101 345 }
Chris@101 346 }
Chris@101 347 using std::ldexp;
Chris@101 348 result = do_cast<To>(q);
Chris@101 349 result = ldexp(result, -shift);
Chris@101 350 if(s)
Chris@101 351 result = -result;
Chris@101 352 }
Chris@101 353 template <class To, class Integer>
Chris@101 354 inline typename disable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
Chris@101 355 generic_convert_rational_to_float_imp(To& result, Integer& num, Integer& denom, const mpl::false_& tag)
Chris@101 356 {
Chris@101 357 number<To> t;
Chris@101 358 generic_convert_rational_to_float_imp(t, num, denom, tag);
Chris@101 359 result = t.backend();
Chris@101 360 }
Chris@101 361
Chris@16 362 template <class To, class From>
Chris@101 363 inline void generic_convert_rational_to_float(To& result, const From& f)
Chris@16 364 {
Chris@101 365 //
Chris@101 366 // Type From is always a Backend to number<>, or an
Chris@101 367 // instance of number<>, but we allow
Chris@101 368 // To to be either a Backend type, or a real number type,
Chris@101 369 // that way we can call this from generic conversions, and
Chris@101 370 // from specific conversions to built in types.
Chris@101 371 //
Chris@101 372 typedef typename mpl::if_c<is_number<From>::value, From, number<From> >::type actual_from_type;
Chris@101 373 typedef typename mpl::if_c<is_number<To>::value || is_floating_point<To>::value, To, number<To> >::type actual_to_type;
Chris@101 374 typedef typename component_type<actual_from_type>::type integer_type;
Chris@101 375 typedef mpl::bool_<!std::numeric_limits<integer_type>::is_specialized
Chris@101 376 || std::numeric_limits<integer_type>::is_bounded
Chris@101 377 || !std::numeric_limits<actual_to_type>::is_specialized
Chris@101 378 || !std::numeric_limits<actual_to_type>::is_bounded
Chris@101 379 || (std::numeric_limits<actual_to_type>::radix != 2)> dispatch_tag;
Chris@16 380
Chris@101 381 integer_type n(numerator(static_cast<actual_from_type>(f))), d(denominator(static_cast<actual_from_type>(f)));
Chris@101 382 generic_convert_rational_to_float_imp(result, n, d, dispatch_tag());
Chris@101 383 }
Chris@101 384
Chris@101 385 template <class To, class From>
Chris@101 386 inline void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_floating_point>& /*to_type*/, const mpl::int_<number_kind_rational>& /*from_type*/)
Chris@101 387 {
Chris@101 388 generic_convert_rational_to_float(to, from);
Chris@101 389 }
Chris@101 390
Chris@101 391 template <class To, class From>
Chris@101 392 void generic_interconvert_float2rational(To& to, const From& from, const mpl::int_<2>& /*radix*/)
Chris@101 393 {
Chris@101 394 typedef typename mpl::front<typename To::unsigned_types>::type ui_type;
Chris@101 395 static const int shift = std::numeric_limits<long long>::digits;
Chris@101 396 typename From::exponent_type e;
Chris@101 397 typename component_type<number<To> >::type num, denom;
Chris@101 398 number<From> val(from);
Chris@101 399 val = frexp(val, &e);
Chris@101 400 while(val)
Chris@101 401 {
Chris@101 402 val = ldexp(val, shift);
Chris@101 403 e -= shift;
Chris@101 404 long long ll = boost::math::lltrunc(val);
Chris@101 405 val -= ll;
Chris@101 406 num <<= shift;
Chris@101 407 num += ll;
Chris@101 408 }
Chris@101 409 denom = ui_type(1u);
Chris@101 410 if(e < 0)
Chris@101 411 denom <<= -e;
Chris@101 412 else if(e > 0)
Chris@101 413 num <<= e;
Chris@101 414 assign_components(to, num.backend(), denom.backend());
Chris@101 415 }
Chris@101 416
Chris@101 417 template <class To, class From, int Radix>
Chris@101 418 void generic_interconvert_float2rational(To& to, const From& from, const mpl::int_<Radix>& /*radix*/)
Chris@101 419 {
Chris@101 420 //
Chris@101 421 // This is almost the same as the binary case above, but we have to use
Chris@101 422 // scalbn and ilogb rather than ldexp and frexp, we also only extract
Chris@101 423 // one Radix digit at a time which is terribly inefficient!
Chris@101 424 //
Chris@101 425 typedef typename mpl::front<typename To::unsigned_types>::type ui_type;
Chris@101 426 typename From::exponent_type e;
Chris@101 427 typename component_type<To>::type num, denom;
Chris@101 428 number<From> val(from);
Chris@101 429 e = ilogb(val);
Chris@101 430 val = scalbn(val, -e);
Chris@101 431 while(val)
Chris@101 432 {
Chris@101 433 long long ll = boost::math::lltrunc(val);
Chris@101 434 val -= ll;
Chris@101 435 val = scalbn(val, 1);
Chris@101 436 num *= Radix;
Chris@101 437 num += ll;
Chris@101 438 --e;
Chris@101 439 }
Chris@101 440 ++e;
Chris@101 441 denom = ui_type(Radix);
Chris@101 442 denom = pow(denom, abs(e));
Chris@101 443 if(e > 0)
Chris@101 444 {
Chris@101 445 num *= denom;
Chris@101 446 denom = 1;
Chris@101 447 }
Chris@101 448 assign_components(to, num, denom);
Chris@101 449 }
Chris@101 450
Chris@101 451 template <class To, class From>
Chris@101 452 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_rational>& /*to_type*/, const mpl::int_<number_kind_floating_point>& /*from_type*/)
Chris@101 453 {
Chris@101 454 generic_interconvert_float2rational(to, from, mpl::int_<std::numeric_limits<number<From> >::radix>());
Chris@16 455 }
Chris@16 456
Chris@16 457 }}} // namespaces
Chris@16 458
Chris@101 459 #ifdef BOOST_MSVC
Chris@101 460 #pragma warning(pop)
Chris@101 461 #endif
Chris@101 462
Chris@16 463 #endif // BOOST_MP_GENERIC_INTERCONVERT_HPP
Chris@16 464