boost/multiprecision
2
0
Fork 0
mirror of https://github.com/boostorg/multiprecision.git synced 2026-02-19 14:32:35 +00:00
Code Issues Projects Releases Wiki Activity

Fix variable precision mpfr/mpfi to give correct results when the precision changes.

Browse Source
This commit is contained in:
jzmaddock
2016-05-31 10:04:52 +01:00
parent 45008b5c55
commit 8e8f4fabeb
2 changed files with 321 additions and 9 deletions
Download Patch File Download Diff File Expand all files Collapse all files

138
include/boost/multiprecision/mpfi.hpp
View File

@@ -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

192
include/boost/multiprecision/mpfr.hpp
View File

@@ -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