1 // (C) Copyright John Maddock 2014.
2 // Use, modification and distribution are subject to the
3 // Boost Software License, Version 1.0. (See accompanying file
4 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
5
6 #ifndef BOOST_MATH_TOOLS_MP_T
7 #define BOOST_MATH_TOOLS_MP_T
8
9 #ifndef BOOST_MATH_PRECISION
10 #define BOOST_MATH_PRECISION 1000
11 #endif
12
13 #if defined(BOOST_MATH_USE_MPFR)
14
15 #include <boost/multiprecision/mpfr.hpp>
16
17 typedef boost::multiprecision::number<boost::multiprecision::mpfr_float_backend<BOOST_MATH_PRECISION *301L / 1000L> > mp_t;
18
19 #else
20
21 #include <boost/multiprecision/cpp_bin_float.hpp>
22
23 typedef boost::multiprecision::number<boost::multiprecision::cpp_bin_float<BOOST_MATH_PRECISION, boost::multiprecision::digit_base_2> > mp_t;
24
25 #endif
26
ConvPrec(mp_t arg,int digits)27 inline mp_t ConvPrec(mp_t arg, int digits)
28 {
29 int e1, e2;
30 mp_t t = frexp(arg, &e1);
31 t = frexp(floor(ldexp(t, digits + 1)), &e2);
32 return ldexp(t, e1);
33 }
34
relative_error(mp_t a,mp_t b)35 inline mp_t relative_error(mp_t a, mp_t b)
36 {
37 mp_t min_val = boost::math::tools::min_value<mp_t>();
38 mp_t max_val = boost::math::tools::max_value<mp_t>();
39
40 if((a != 0) && (b != 0))
41 {
42 // mp_tODO: use isfinite:
43 if(fabs(b) >= max_val)
44 {
45 if(fabs(a) >= max_val)
46 return 0; // one infinity is as good as another!
47 }
48 // If the result is denormalised, treat all denorms as equivalent:
49 if((a < min_val) && (a > 0))
50 a = min_val;
51 else if((a > -min_val) && (a < 0))
52 a = -min_val;
53 if((b < min_val) && (b > 0))
54 b = min_val;
55 else if((b > -min_val) && (b < 0))
56 b = -min_val;
57 return (std::max)(fabs((a - b) / a), fabs((a - b) / b));
58 }
59
60 // Handle special case where one or both are zero:
61 if(min_val == 0)
62 return fabs(a - b);
63 if(fabs(a) < min_val)
64 a = min_val;
65 if(fabs(b) < min_val)
66 b = min_val;
67 return (std::max)(fabs((a - b) / a), fabs((a - b) / b));
68 }
69
70
71 #endif // BOOST_MATH_TOOLS_MP_T
72