diff DEPENDENCIES/generic/include/boost/multiprecision/detail/generic_interconvert.hpp @ 101:c530137014c0

Update Boost headers (1.58.0)

--- a/DEPENDENCIES/generic/include/boost/multiprecision/detail/generic_interconvert.hpp	Fri Sep 04 12:01:02 2015 +0100
+++ b/DEPENDENCIES/generic/include/boost/multiprecision/detail/generic_interconvert.hpp	Mon Sep 07 11:12:49 2015 +0100
@@ -8,9 +8,25 @@
 
 #include <boost/multiprecision/detail/default_ops.hpp>
 
+#ifdef BOOST_MSVC
+#pragma warning(push)
+#pragma warning(disable:4127)
+#endif
+
 namespace boost{ namespace multiprecision{ namespace detail{
 
 template <class To, class From>
+inline To do_cast(const From & from)
+{
+   return static_cast<To>(from);
+}
+template <class To, class B, ::boost::multiprecision::expression_template_option et>
+inline To do_cast(const number<B, et>& from)
+{
+   return from.template convert_to<To>();
+}
+
+template <class To, class From>
 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*/)
 {
    using default_ops::eval_get_sign;
@@ -136,19 +152,19 @@
    //
    int c = eval_fpclassify(from);
 
-   if(c == FP_ZERO) 
+   if(c == (int)FP_ZERO)
    {
       to = ui_type(0);
       return;
    }
-   else if(c == FP_NAN)
+   else if(c == (int)FP_NAN)
    {
-      to = "nan";
+      to = static_cast<const char*>("nan");
       return;
    }
-   else if(c == FP_INFINITE)
+   else if(c == (int)FP_INFINITE)
    {
-      to = "inf";
+      to = static_cast<const char*>("inf");
       if(eval_get_sign(from) < 0)
          to.negate();
       return;
@@ -177,7 +193,7 @@
    typedef typename To::exponent_type to_exponent;
    if((e > (std::numeric_limits<to_exponent>::max)()) || (e < (std::numeric_limits<to_exponent>::min)()))
    {
-      to = "inf";
+      to = static_cast<const char*>("inf");
       if(eval_get_sign(from) < 0)
          to.negate();
       return;
@@ -210,19 +226,239 @@
    assign_components(to, n.backend(), d.backend());
 }
 
+template <class R, class LargeInteger>
+R safe_convert_to_float(const LargeInteger& i)
+{
+   using std::ldexp;
+   if(!i)
+      return R(0);
+   if(std::numeric_limits<R>::is_specialized && std::numeric_limits<R>::max_exponent)
+   {
+      LargeInteger val(i);
+      if(val.sign() < 0)
+         val = -val;
+      unsigned mb = msb(val);
+      if(mb >= std::numeric_limits<R>::max_exponent)
+      {
+         int scale_factor = (int)mb + 1 - std::numeric_limits<R>::max_exponent;
+         BOOST_ASSERT(scale_factor >= 1);
+         val >>= scale_factor;
+         R result = val.template convert_to<R>();
+         if(std::numeric_limits<R>::digits == 0 || std::numeric_limits<R>::digits >= std::numeric_limits<R>::max_exponent)
+         {
+            //
+            // Calculate and add on the remainder, only if there are more
+            // digits in the mantissa that the size of the exponent, in 
+            // other words if we are dropping digits in the conversion
+            // otherwise:
+            //
+            LargeInteger remainder(i);
+            remainder &= (LargeInteger(1) << scale_factor) - 1;
+            result += ldexp(safe_convert_to_float<R>(remainder), -scale_factor);
+         }
+         return i.sign() < 0 ? static_cast<R>(-result) : result;
+      }
+   }
+   return i.template convert_to<R>();
+}
+
+template <class To, class Integer>
+inline typename disable_if_c<is_number<To>::value || is_floating_point<To>::value>::type 
+   generic_convert_rational_to_float_imp(To& result, const Integer& n, const Integer& d, const mpl::true_&)
+{
+   //
+   // If we get here, then there's something about one type or the other
+   // that prevents an exactly rounded result from being calculated
+   // (or at least it's not clear how to implement such a thing).
+   //
+   using default_ops::eval_divide;
+   number<To> fn(safe_convert_to_float<number<To> >(n)), fd(safe_convert_to_float<number<To> >(d));
+   eval_divide(result, fn.backend(), fd.backend());
+}
+template <class To, class Integer>
+inline typename enable_if_c<is_number<To>::value || is_floating_point<To>::value>::type 
+   generic_convert_rational_to_float_imp(To& result, const Integer& n, const Integer& d, const mpl::true_&)
+{
+   //
+   // If we get here, then there's something about one type or the other
+   // that prevents an exactly rounded result from being calculated
+   // (or at least it's not clear how to implement such a thing).
+   //
+   To fd(safe_convert_to_float<To>(d));
+   result = safe_convert_to_float<To>(n);
+   result /= fd;
+}
+
+template <class To, class Integer>
+typename enable_if_c<is_number<To>::value || is_floating_point<To>::value>::type 
+   generic_convert_rational_to_float_imp(To& result, Integer& num, Integer& denom, const mpl::false_&)
+{
+   //
+   // If we get here, then the precision of type To is known, and the integer type is unbounded
+   // so we can use integer division plus manipulation of the remainder to get an exactly
+   // rounded result.
+   //
+   if(num == 0)
+   {
+      result = 0;
+      return;
+   }
+   bool s = false;
+   if(num < 0)
+   {
+      s = true;
+      num = -num;
+   }
+   int denom_bits = msb(denom);
+   int shift = std::numeric_limits<To>::digits + denom_bits - msb(num) + 1;
+   if(shift > 0)
+      num <<= shift;
+   else if(shift < 0)
+      denom <<= std::abs(shift);
+   Integer q, r;
+   divide_qr(num, denom, q, r);
+   int q_bits = msb(q);
+   if(q_bits == std::numeric_limits<To>::digits)
+   {
+      //
+      // Round up if 2 * r > denom:
+      //
+      r <<= 1;
+      int c = r.compare(denom);
+      if(c > 0)
+         ++q;
+      else if((c == 0) && (q & 1u))
+      {
+         ++q;
+      }
+   }
+   else
+   {
+      BOOST_ASSERT(q_bits == 1 + std::numeric_limits<To>::digits);
+      //
+      // We basically already have the rounding info:
+      //
+      if(q & 1u)
+      {
+         if(r || (q & 2u))
+            ++q;
+      }
+   }
+   using std::ldexp;
+   result = do_cast<To>(q);
+   result = ldexp(result, -shift);
+   if(s)
+      result = -result;
+}
+template <class To, class Integer>
+inline typename disable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
+   generic_convert_rational_to_float_imp(To& result, Integer& num, Integer& denom, const mpl::false_& tag)
+{
+   number<To> t;
+   generic_convert_rational_to_float_imp(t, num, denom, tag);
+   result = t.backend();
+}
+
 template <class To, class From>
-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*/)
+inline void generic_convert_rational_to_float(To& result, const From& f)
 {
-   typedef typename component_type<number<From> >::type   from_component_type;
-   using default_ops::eval_divide;
+   //
+   // Type From is always a Backend to number<>, or an
+   // instance of number<>, but we allow
+   // To to be either a Backend type, or a real number type,
+   // that way we can call this from generic conversions, and
+   // from specific conversions to built in types.
+   //
+   typedef typename mpl::if_c<is_number<From>::value, From, number<From> >::type actual_from_type;
+   typedef typename mpl::if_c<is_number<To>::value || is_floating_point<To>::value, To, number<To> >::type actual_to_type;
+   typedef typename component_type<actual_from_type>::type integer_type;
+   typedef mpl::bool_<!std::numeric_limits<integer_type>::is_specialized 
+                      || std::numeric_limits<integer_type>::is_bounded
+                      || !std::numeric_limits<actual_to_type>::is_specialized 
+                      || !std::numeric_limits<actual_to_type>::is_bounded
+                      || (std::numeric_limits<actual_to_type>::radix != 2)> dispatch_tag;
 
-   number<From> t(from);
-   from_component_type n(numerator(t)), d(denominator(t));
-   number<To> fn(n), fd(d);
-   eval_divide(to, fn.backend(), fd.backend());
+   integer_type n(numerator(static_cast<actual_from_type>(f))), d(denominator(static_cast<actual_from_type>(f)));
+   generic_convert_rational_to_float_imp(result, n, d, dispatch_tag());
+}
+
+template <class To, class From>
+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*/)
+{
+   generic_convert_rational_to_float(to, from);
+}
+
+template <class To, class From>
+void generic_interconvert_float2rational(To& to, const From& from, const mpl::int_<2>& /*radix*/)
+{
+   typedef typename mpl::front<typename To::unsigned_types>::type ui_type;
+   static const int shift = std::numeric_limits<long long>::digits;
+   typename From::exponent_type e;
+   typename component_type<number<To> >::type num, denom;
+   number<From> val(from);
+   val = frexp(val, &e);
+   while(val)
+   {
+      val = ldexp(val, shift);
+      e -= shift;
+      long long ll = boost::math::lltrunc(val);
+      val -= ll;
+      num <<= shift;
+      num += ll;
+   }
+   denom = ui_type(1u);
+   if(e < 0)
+      denom <<= -e;
+   else if(e > 0)
+      num <<= e;
+   assign_components(to, num.backend(), denom.backend());
+}
+
+template <class To, class From, int Radix>
+void generic_interconvert_float2rational(To& to, const From& from, const mpl::int_<Radix>& /*radix*/)
+{
+   //
+   // This is almost the same as the binary case above, but we have to use
+   // scalbn and ilogb rather than ldexp and frexp, we also only extract
+   // one Radix digit at a time which is terribly inefficient!
+   //
+   typedef typename mpl::front<typename To::unsigned_types>::type ui_type;
+   typename From::exponent_type e;
+   typename component_type<To>::type num, denom;
+   number<From> val(from);
+   e = ilogb(val);
+   val = scalbn(val, -e);
+   while(val)
+   {
+      long long ll = boost::math::lltrunc(val);
+      val -= ll;
+      val = scalbn(val, 1);
+      num *= Radix;
+      num += ll;
+      --e;
+   }
+   ++e;
+   denom = ui_type(Radix);
+   denom = pow(denom, abs(e));
+   if(e > 0)
+   {
+      num *= denom;
+      denom = 1;
+   }
+   assign_components(to, num, denom);
+}
+
+template <class To, class From>
+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*/)
+{
+   generic_interconvert_float2rational(to, from, mpl::int_<std::numeric_limits<number<From> >::radix>());
 }
 
 }}} // namespaces
 
+#ifdef BOOST_MSVC
+#pragma warning(pop)
+#endif
+
 #endif  // BOOST_MP_GENERIC_INTERCONVERT_HPP