1 // (C) Copyright John Maddock 2018.
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)
6 #include <boost/math/tools/fraction.hpp>
9 #include <boost/multiprecision/cpp_complex.hpp>
13 struct golden_ratio_fraction
15 typedef T result_type
;
17 result_type
operator()()
35 typedef std::pair
<T
, T
> result_type
;
37 std::pair
<T
, T
> operator()()
40 return std::make_pair(a
, b
);
48 tan_fraction
<T
> fract(a
);
49 return a
/ continued_fraction_b(fract
, std::numeric_limits
<T
>::epsilon());
54 struct expint_fraction
56 typedef std::pair
<T
, T
> result_type
;
57 expint_fraction(unsigned n_
, T z_
) : b(z_
+ T(n_
)), i(-1), n(n_
) {}
58 std::pair
<T
, T
> operator()()
60 std::pair
<T
, T
> result
= std::make_pair(-static_cast<T
>((i
+ 1) * (n
+ i
)), b
);
73 inline std::complex<T
> expint_as_fraction(unsigned n
, std::complex<T
> const& z
)
75 boost::uintmax_t max_iter
= 1000;
76 expint_fraction
<std::complex<T
> > f(n
, z
);
77 std::complex<T
> result
= boost::math::tools::continued_fraction_b(
79 std::complex<T
>(std::numeric_limits
<T
>::epsilon()),
81 result
= exp(-z
) / result
;
85 //[cf_upper_gamma_fraction
87 struct upper_incomplete_gamma_fract
90 typedef typename
T::value_type scalar_type
;
94 typedef std::pair
<T
, T
> result_type
;
96 upper_incomplete_gamma_fract(T a1
, T z1
)
97 : z(z1
- a1
+ scalar_type(1)), a(a1
), k(0)
101 result_type
operator()()
105 return result_type(scalar_type(k
) * (a
- scalar_type(k
)), z
);
111 inline std::complex<T
> gamma_Q_as_fraction(const std::complex<T
>& a
, const std::complex<T
>& z
)
113 upper_incomplete_gamma_fract
<std::complex<T
> > f(a
, z
);
114 std::complex<T
> eps(std::numeric_limits
<T
>::epsilon());
115 return pow(z
, a
) / (exp(z
) *(z
- a
+ T(1) + boost::math::tools::continued_fraction_a(f
, eps
)));
118 inline boost::multiprecision::cpp_complex_50
gamma_Q_as_fraction(const boost::multiprecision::cpp_complex_50
& a
, const boost::multiprecision::cpp_complex_50
& z
)
120 upper_incomplete_gamma_fract
<boost::multiprecision::cpp_complex_50
> f(a
, z
);
121 boost::multiprecision::cpp_complex_50
eps(std::numeric_limits
<boost::multiprecision::cpp_complex_50::value_type
>::epsilon());
122 return pow(z
, a
) / (exp(z
) * (z
- a
+ 1 + boost::math::tools::continued_fraction_a(f
, eps
)));
128 using namespace boost::math::tools
;
131 golden_ratio_fraction
<double> func
;
132 double gr
= continued_fraction_a(
134 std::numeric_limits
<double>::epsilon());
135 std::cout
<< "The golden ratio is: " << gr
<< std::endl
;
138 std::cout
<< tan(0.5) << std::endl
;
140 std::complex<double> arg(3, 2);
141 std::cout
<< expint_as_fraction(5, arg
) << std::endl
;
143 std::complex<double> a(3, 3), z(3, 2);
144 std::cout
<< gamma_Q_as_fraction(a
, z
) << std::endl
;
146 boost::multiprecision::cpp_complex_50
am(3, 3), zm(3, 2);
147 std::cout
<< gamma_Q_as_fraction(am
, zm
) << std::endl
;