From 8e8f4fabebac305fc360cb09c7b0bb1e439f4962 Mon Sep 17 00:00:00 2001
From: jzmaddock <john@johnmaddock.co.uk>
Date: Tue, 31 May 2016 10:04:52 +0100
Subject: [PATCH] Fix variable precision mpfr/mpfi to give correct results when
 the precision changes.

---
 include/boost/multiprecision/mpfi.hpp | 138 +++++++++++++++++-
 include/boost/multiprecision/mpfr.hpp | 192 +++++++++++++++++++++++++-
 2 files changed, 321 insertions(+), 9 deletions(-)

diff --git a/include/boost/multiprecision/mpfi.hpp b/include/boost/multiprecision/mpfi.hpp
index d78047e8..00a69298 100644
--- a/include/boost/multiprecision/mpfi.hpp
+++ b/include/boost/multiprecision/mpfi.hpp
@@ -12,11 +12,17 @@
 #include <boost/multiprecision/detail/big_lanczos.hpp>
 #include <boost/multiprecision/detail/digits.hpp>
 #include <boost/multiprecision/mpfr.hpp>
+#include <boost/multiprecision/logged_adaptor.hpp>
 #include <boost/math/constants/constants.hpp>
+#include <boost/functional/hash_fwd.hpp>
 #include <mpfi.h>
 #include <cmath>
 #include <algorithm>
 
+#ifndef BOOST_MULTIPRECISION_MPFI_DEFAULT_PRECISION
+#  define BOOST_MULTIPRECISION_MPFI_DEFAULT_PRECISION 20
+#endif
+
 namespace boost{
 namespace multiprecision{
 namespace backends{
@@ -315,7 +321,7 @@ protected:
    mpfi_t m_data;
    static unsigned& get_default_precision() BOOST_NOEXCEPT
    {
-      static unsigned val = 50;
+      static unsigned val = BOOST_MULTIPRECISION_MPFI_DEFAULT_PRECISION;
       return val;
    }
 };
@@ -1175,7 +1181,7 @@ inline int digits<boost::multiprecision::mpfi_float>()
    BOOST_NOEXCEPT
 #endif
 {
-   return boost::multiprecision::backends::detail::get_default_precision();
+   return  multiprecision::detail::digits10_2_2(boost::multiprecision::mpfi_float::default_precision());
 }
 template <>
 inline int digits<boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<0>, boost::multiprecision::et_off> >()
@@ -1183,9 +1189,109 @@ inline int digits<boost::multiprecision::number<boost::multiprecision::mpfi_floa
    BOOST_NOEXCEPT
 #endif
 {
-   return boost::multiprecision::backends::detail::get_default_precision();
+   return  multiprecision::detail::digits10_2_2(boost::multiprecision::mpfi_float::default_precision());
 }
 
+template <>
+inline boost::multiprecision::mpfi_float
+max_value<boost::multiprecision::mpfi_float>()
+{
+   boost::multiprecision::mpfi_float result(0.5);
+   mpfi_mul_2exp(result.backend().data(), result.backend().data(), mpfr_get_emax());
+   //BOOST_ASSERT(mpfi_number_p(result.backend().data()));
+   return result;
+}
+
+template <>
+inline boost::multiprecision::mpfi_float
+min_value<boost::multiprecision::mpfi_float>()
+{
+   boost::multiprecision::mpfi_float result(0.5);
+   mpfi_div_2exp(result.backend().data(), result.backend().data(), -mpfr_get_emin());
+   //BOOST_ASSERT(mpfi_number_p(result.backend().data()));
+   return result;
+}
+
+template <>
+inline boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<0>, boost::multiprecision::et_off>
+max_value<boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<0>, boost::multiprecision::et_off> >()
+{
+   boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<0>, boost::multiprecision::et_off> result(0.5);
+   mpfi_mul_2exp(result.backend().data(), result.backend().data(), mpfr_get_emax());
+   //BOOST_ASSERT(mpfi_number_p(result.backend().data()));
+   return result;
+}
+
+template <>
+inline boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<0>, boost::multiprecision::et_off>
+min_value<boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<0>, boost::multiprecision::et_off> >()
+{
+   boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<0>, boost::multiprecision::et_off> result(0.5);
+   mpfi_div_2exp(result.backend().data(), result.backend().data(), -mpfr_get_emin());
+   //BOOST_ASSERT(mpfi_number_p(result.backend().data()));
+   return result;
+}
+
+// mpfi gets used with logged_adaptor fairly often, so specialize for that use case as well:
+typedef boost::multiprecision::number<boost::multiprecision::backends::logged_adaptor<boost::multiprecision::mpfi_float::backend_type>, boost::multiprecision::et_on> logged_type1;
+typedef boost::multiprecision::number<boost::multiprecision::backends::logged_adaptor<boost::multiprecision::mpfi_float::backend_type>, boost::multiprecision::et_off> logged_type2;
+
+template <>
+inline int digits<logged_type1>()
+#ifdef BOOST_MATH_NOEXCEPT
+BOOST_NOEXCEPT
+#endif
+{
+   return  multiprecision::detail::digits10_2_2(logged_type1::default_precision());
+}
+template <>
+inline int digits<logged_type2 >()
+#ifdef BOOST_MATH_NOEXCEPT
+BOOST_NOEXCEPT
+#endif
+{
+   return  multiprecision::detail::digits10_2_2(logged_type1::default_precision());
+}
+
+template <>
+inline logged_type1
+max_value<logged_type1>()
+{
+   logged_type1 result(0.5);
+   mpfi_mul_2exp(result.backend().value().data(), result.backend().value().data(), mpfr_get_emax());
+   //BOOST_ASSERT(mpfi_number_p(result.backend().data()));
+   return result;
+}
+
+template <>
+inline logged_type1
+min_value<logged_type1>()
+{
+   logged_type1 result(0.5);
+   mpfi_div_2exp(result.backend().value().data(), result.backend().value().data(), -mpfr_get_emin());
+   //BOOST_ASSERT(mpfi_number_p(result.backend().data()));
+   return result;
+}
+
+template <>
+inline logged_type2
+max_value<logged_type2 >()
+{
+   logged_type2 result(0.5);
+   mpfi_mul_2exp(result.backend().value().data(), result.backend().value().data(), mpfr_get_emax());
+   //BOOST_ASSERT(mpfi_number_p(result.backend().data()));
+   return result;
+}
+
+template <>
+inline logged_type2
+min_value<logged_type2 >()
+{
+   logged_type2 result(0.5);
+   mpfi_div_2exp(result.backend().value().data(), result.backend().value().data(), -mpfr_get_emin());
+   //BOOST_ASSERT(mpfi_number_p(result.backend().data()));
+   return result;
+}
 } // namespace tools
 
 namespace constants{ namespace detail{
@@ -1239,13 +1345,19 @@ struct constant_pi<boost::multiprecision::number<boost::multiprecision::mpfi_flo
       }
       return result;
    }
+   static inline result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfi_const_pi(result.backend().data());
+      return result;
+   }
 };
 template<unsigned Digits10, boost::multiprecision::expression_template_option ExpressionTemplates>
 struct constant_ln_two<boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<Digits10>, ExpressionTemplates> >
 {
    typedef boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<Digits10>, ExpressionTemplates> result_type;
    template<int N>
-   static inline result_type const& get(const mpl::int_<N>&)
+   static inline result_type get(const mpl::int_<N>&)
    {
       mpfi_initializer<result_type>::force_instantiate();
       static result_type result;
@@ -1257,6 +1369,12 @@ struct constant_ln_two<boost::multiprecision::number<boost::multiprecision::mpfi
       }
       return result;
    }
+   static inline result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfi_const_log2(result.backend().data());
+      return result;
+   }
 };
 template<unsigned Digits10, boost::multiprecision::expression_template_option ExpressionTemplates>
 struct constant_euler<boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<Digits10>, ExpressionTemplates> >
@@ -1275,6 +1393,12 @@ struct constant_euler<boost::multiprecision::number<boost::multiprecision::mpfi_
       }
       return result;
    }
+   static inline result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfi_const_euler(result.backend().data());
+      return result;
+   }
 };
 template<unsigned Digits10, boost::multiprecision::expression_template_option ExpressionTemplates>
 struct constant_catalan<boost::multiprecision::number<boost::multiprecision::mpfi_float_backend<Digits10>, ExpressionTemplates> >
@@ -1293,6 +1417,12 @@ struct constant_catalan<boost::multiprecision::number<boost::multiprecision::mpf
       }
       return result;
    }
+   static inline result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfi_const_catalan(result.backend().data());
+      return result;
+   }
 };
 
 }} // namespaces
diff --git a/include/boost/multiprecision/mpfr.hpp b/include/boost/multiprecision/mpfr.hpp
index 82dd166c..5af90c09 100644
--- a/include/boost/multiprecision/mpfr.hpp
+++ b/include/boost/multiprecision/mpfr.hpp
@@ -7,6 +7,7 @@
 #define BOOST_MATH_BN_MPFR_HPP
 
 #include <boost/multiprecision/number.hpp>
+#include <boost/multiprecision/debug_adaptor.hpp>
 #include <boost/multiprecision/gmp.hpp>
 #include <boost/math/special_functions/fpclassify.hpp>
 #include <boost/cstdint.hpp>
@@ -16,6 +17,10 @@
 #include <cmath>
 #include <algorithm>
 
+#ifndef BOOST_MULTIPRECISION_MPFR_DEFAULT_PRECISION
+#  define BOOST_MULTIPRECISION_MPFR_DEFAULT_PRECISION 20
+#endif
+
 namespace boost{
 namespace multiprecision{
 
@@ -58,7 +63,6 @@ struct mpfr_cleanup
 template <bool b>
 typename mpfr_cleanup<b>::initializer const mpfr_cleanup<b>::init;
 
-inline long get_default_precision() { return 50; }
 
 template <unsigned digits10, mpfr_allocation_type AllocationType>
 struct mpfr_float_imp;
@@ -314,7 +318,7 @@ struct mpfr_float_imp<digits10, allocate_dynamic>
    }
    ~mpfr_float_imp() BOOST_NOEXCEPT
    {
-      if(m_data[0]._mpfr_d)
+       if(m_data[0]._mpfr_d)
          mpfr_clear(m_data);
       detail::mpfr_cleanup<true>::force_instantiate();
    }
@@ -360,7 +364,7 @@ protected:
    mpfr_t m_data;
    static unsigned& get_default_precision() BOOST_NOEXCEPT
    {
-      static unsigned val = 50;
+      static unsigned val = BOOST_MULTIPRECISION_MPFR_DEFAULT_PRECISION;
       return val;
    }
 };
@@ -1481,7 +1485,7 @@ inline int digits<boost::multiprecision::mpfr_float>()
    BOOST_NOEXCEPT
 #endif
 {
-   return boost::multiprecision::backends::detail::get_default_precision();
+   return multiprecision::detail::digits10_2_2(boost::multiprecision::mpfr_float::default_precision());
 }
 template <>
 inline int digits<boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<0>, boost::multiprecision::et_off> >()
@@ -1489,7 +1493,7 @@ inline int digits<boost::multiprecision::number<boost::multiprecision::mpfr_floa
    BOOST_NOEXCEPT
 #endif
 {
-   return boost::multiprecision::backends::detail::get_default_precision();
+   return multiprecision::detail::digits10_2_2(boost::multiprecision::mpfr_float::default_precision());
 }
 
 template <>
@@ -1532,6 +1536,51 @@ inline boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<0
    return result;
 }
 
+template <>
+inline int digits<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float::backend_type> >>()
+#ifdef BOOST_MATH_NOEXCEPT
+BOOST_NOEXCEPT
+#endif
+{
+   return multiprecision::detail::digits10_2_2(boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float::backend_type> >::default_precision());
+}
+template <>
+inline int digits<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<0> >, boost::multiprecision::et_off> >()
+#ifdef BOOST_MATH_NOEXCEPT
+BOOST_NOEXCEPT
+#endif
+{
+   return multiprecision::detail::digits10_2_2(boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float::backend_type> >::default_precision());
+}
+
+template <>
+inline boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float::backend_type> >
+max_value<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float::backend_type> >>()
+{
+   return max_value<boost::multiprecision::mpfr_float>().backend();
+}
+
+template <>
+inline boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float::backend_type> >
+min_value<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float::backend_type> >>()
+{
+   return min_value<boost::multiprecision::mpfr_float>().backend();
+}
+
+template <>
+inline boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<0> >, boost::multiprecision::et_off>
+max_value<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<0> >, boost::multiprecision::et_off> >()
+{
+   return max_value<boost::multiprecision::mpfr_float>().backend();
+}
+
+template <>
+inline boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<0> >, boost::multiprecision::et_off>
+min_value<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<0> >, boost::multiprecision::et_off> >()
+{
+   return min_value<boost::multiprecision::mpfr_float>().backend();
+}
+
 } // namespace tools
 
 namespace constants{ namespace detail{
@@ -1584,6 +1633,12 @@ struct constant_pi<boost::multiprecision::number<boost::multiprecision::mpfr_flo
       }
       return result;
    }
+   static inline const result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfr_const_pi(result.backend().data(), GMP_RNDN);
+      return result;
+   }
 };
 template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
 struct constant_ln_two<boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType>, ExpressionTemplates> >
@@ -1602,6 +1657,12 @@ struct constant_ln_two<boost::multiprecision::number<boost::multiprecision::mpfr
       }
       return result;
    }
+   static inline const result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfr_const_log2(result.backend().data(), GMP_RNDN);
+      return result;
+   }
 };
 template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
 struct constant_euler<boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType>, ExpressionTemplates> >
@@ -1620,6 +1681,12 @@ struct constant_euler<boost::multiprecision::number<boost::multiprecision::mpfr_
       }
       return result;
    }
+   static inline const result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfr_const_euler(result.backend().data(), GMP_RNDN);
+      return result;
+   }
 };
 template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
 struct constant_catalan<boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType>, ExpressionTemplates> >
@@ -1638,6 +1705,109 @@ struct constant_catalan<boost::multiprecision::number<boost::multiprecision::mpf
       }
       return result;
    }
+   static inline const result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfr_const_catalan(result.backend().data(), GMP_RNDN);
+      return result;
+   }
+};
+
+template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
+struct constant_pi<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> >
+{
+   typedef boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> result_type;
+   template<int N>
+   static inline const result_type& get(const mpl::int_<N>&)
+   {
+      detail::mpfr_constant_initializer<constant_pi<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> >, N>::force_instantiate();
+      static result_type result;
+      static bool has_init = false;
+      if(!has_init)
+      {
+         mpfr_const_pi(result.backend().value().data(), GMP_RNDN);
+         has_init = true;
+      }
+      return result;
+   }
+   static inline const result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfr_const_pi(result.backend().value().data(), GMP_RNDN);
+      return result;
+   }
+};
+template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
+struct constant_ln_two<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> >
+{
+   typedef boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> result_type;
+   template<int N>
+   static inline const result_type& get(const mpl::int_<N>&)
+   {
+      detail::mpfr_constant_initializer<constant_ln_two<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> >, N>::force_instantiate();
+      static result_type result;
+      static bool init = false;
+      if(!init)
+      {
+         mpfr_const_log2(result.backend().value().data(), GMP_RNDN);
+         init = true;
+      }
+      return result;
+   }
+   static inline const result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfr_const_log2(result.backend().value().data(), GMP_RNDN);
+      return result;
+   }
+};
+template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
+struct constant_euler<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> >
+{
+   typedef boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> result_type;
+   template<int N>
+   static inline const result_type& get(const mpl::int_<N>&)
+   {
+      detail::mpfr_constant_initializer<constant_euler<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> >, N>::force_instantiate();
+      static result_type result;
+      static bool init = false;
+      if(!init)
+      {
+         mpfr_const_euler(result.backend().value().data(), GMP_RNDN);
+         init = true;
+      }
+      return result;
+   }
+   static inline const result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfr_const_euler(result.backend().value().data(), GMP_RNDN);
+      return result;
+   }
+};
+template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
+struct constant_catalan<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> >
+{
+   typedef boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> result_type;
+   template<int N>
+   static inline const result_type& get(const mpl::int_<N>&)
+   {
+      detail::mpfr_constant_initializer<constant_catalan<boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> >, N>::force_instantiate();
+      static result_type result;
+      static bool init = false;
+      if(!init)
+      {
+         mpfr_const_catalan(result.backend().value().data(), GMP_RNDN);
+         init = true;
+      }
+      return result;
+   }
+   static inline const result_type get(const mpl::int_<0>&)
+   {
+      result_type result;
+      mpfr_const_catalan(result.backend().value().data(), GMP_RNDN);
+      return result;
+   }
 };
 
 }} // namespaces
@@ -1658,6 +1828,18 @@ inline boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<D
    return (boost::math::signbit)(a) != (boost::math::signbit)(b) ? boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType>, ExpressionTemplates>(-a) : a;
 }
 
+template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
+inline int signbit BOOST_PREVENT_MACRO_SUBSTITUTION(const boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates>& arg)
+{
+   return (arg.backend().value().data()[0]._mpfr_sign < 0) ? 1 : 0;
+}
+
+template<unsigned Digits10, boost::multiprecision::mpfr_allocation_type AllocateType, boost::multiprecision::expression_template_option ExpressionTemplates>
+inline boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates> copysign BOOST_PREVENT_MACRO_SUBSTITUTION(const boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType>, ExpressionTemplates>& a, const boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType>, ExpressionTemplates>& b)
+{
+   return (boost::math::signbit)(a) != (boost::math::signbit)(b) ? boost::multiprecision::number<boost::multiprecision::debug_adaptor<boost::multiprecision::mpfr_float_backend<Digits10, AllocateType> >, ExpressionTemplates>(-a) : a;
+}
+
 
 } // namespace math