1// (C) Copyright John Maddock 2005-2006.
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_FRACTION_INCLUDED
7#define BOOST_MATH_TOOLS_FRACTION_INCLUDED
8
9#ifdef _MSC_VER
10#pragma once
11#endif
12
13#include <boost/math/tools/precision.hpp>
14#include <boost/math/tools/complex.hpp>
15#include <type_traits>
16#include <cstdint>
17#include <cmath>
18
19namespace boost{ namespace math{ namespace tools{
20
21namespace detail
22{
23
24 template <typename T>
25 struct is_pair : public std::false_type{};
26
27 template <typename T, typename U>
28 struct is_pair<std::pair<T,U>> : public std::true_type{};
29
30 template <typename Gen>
31 struct fraction_traits_simple
32 {
33 using result_type = typename Gen::result_type;
34 using value_type = typename Gen::result_type;
35
36 static result_type a(const value_type&) BOOST_MATH_NOEXCEPT(value_type)
37 {
38 return 1;
39 }
40 static result_type b(const value_type& v) BOOST_MATH_NOEXCEPT(value_type)
41 {
42 return v;
43 }
44 };
45
46 template <typename Gen>
47 struct fraction_traits_pair
48 {
49 using value_type = typename Gen::result_type;
50 using result_type = typename value_type::first_type;
51
52 static result_type a(const value_type& v) BOOST_MATH_NOEXCEPT(value_type)
53 {
54 return v.first;
55 }
56 static result_type b(const value_type& v) BOOST_MATH_NOEXCEPT(value_type)
57 {
58 return v.second;
59 }
60 };
61
62 template <typename Gen>
63 struct fraction_traits
64 : public std::conditional<
65 is_pair<typename Gen::result_type>::value,
66 fraction_traits_pair<Gen>,
67 fraction_traits_simple<Gen>>::type
68 {
69 };
70
71 template <typename T, bool = is_complex_type<T>::value>
72 struct tiny_value
73 {
74 // For float, double, and long double, 1/min_value<T>() is finite.
75 // But for mpfr_float and cpp_bin_float, 1/min_value<T>() is inf.
76 // Multiply the min by 16 so that the reciprocal doesn't overflow.
77 static T get() {
78 return 16*tools::min_value<T>();
79 }
80 };
81 template <typename T>
82 struct tiny_value<T, true>
83 {
84 using value_type = typename T::value_type;
85 static T get() {
86 return 16*tools::min_value<value_type>();
87 }
88 };
89
90} // namespace detail
91
92//
93// continued_fraction_b
94// Evaluates:
95//
96// b0 + a1
97// ---------------
98// b1 + a2
99// ----------
100// b2 + a3
101// -----
102// b3 + ...
103//
104// Note that the first a0 returned by generator Gen is discarded.
105//
106template <typename Gen, typename U>
107inline typename detail::fraction_traits<Gen>::result_type continued_fraction_b(Gen& g, const U& factor, std::uintmax_t& max_terms)
108 noexcept(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
109{
110 BOOST_MATH_STD_USING // ADL of std names
111
112 using traits = detail::fraction_traits<Gen>;
113 using result_type = typename traits::result_type;
114 using value_type = typename traits::value_type;
115 using integer_type = typename integer_scalar_type<result_type>::type;
116 using scalar_type = typename scalar_type<result_type>::type;
117
118 integer_type const zero(0), one(1);
119
120 result_type tiny = detail::tiny_value<result_type>::get();
121 scalar_type terminator = abs(factor);
122
123 value_type v = g();
124
125 result_type f, C, D, delta;
126 f = traits::b(v);
127 if(f == zero)
128 f = tiny;
129 C = f;
130 D = 0;
131
132 std::uintmax_t counter(max_terms);
133 do{
134 v = g();
135 D = traits::b(v) + traits::a(v) * D;
136 if(D == result_type(0))
137 D = tiny;
138 C = traits::b(v) + traits::a(v) / C;
139 if(C == zero)
140 C = tiny;
141 D = one/D;
142 delta = C*D;
143 f = f * delta;
144 }while((abs(delta - one) > terminator) && --counter);
145
146 max_terms = max_terms - counter;
147
148 return f;
149}
150
151template <typename Gen, typename U>
152inline typename detail::fraction_traits<Gen>::result_type continued_fraction_b(Gen& g, const U& factor)
153 noexcept(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
154{
155 std::uintmax_t max_terms = (std::numeric_limits<std::uintmax_t>::max)();
156 return continued_fraction_b(g, factor, max_terms);
157}
158
159template <typename Gen>
160inline typename detail::fraction_traits<Gen>::result_type continued_fraction_b(Gen& g, int bits)
161 noexcept(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
162{
163 BOOST_MATH_STD_USING // ADL of std names
164
165 using traits = detail::fraction_traits<Gen>;
166 using result_type = typename traits::result_type;
167
168 result_type factor = ldexp(1.0f, 1 - bits); // 1 / pow(result_type(2), bits);
169 std::uintmax_t max_terms = (std::numeric_limits<std::uintmax_t>::max)();
170 return continued_fraction_b(g, factor, max_terms);
171}
172
173template <typename Gen>
174inline typename detail::fraction_traits<Gen>::result_type continued_fraction_b(Gen& g, int bits, std::uintmax_t& max_terms)
175 noexcept(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
176{
177 BOOST_MATH_STD_USING // ADL of std names
178
179 using traits = detail::fraction_traits<Gen>;
180 using result_type = typename traits::result_type;
181
182 result_type factor = ldexp(1.0f, 1 - bits); // 1 / pow(result_type(2), bits);
183 return continued_fraction_b(g, factor, max_terms);
184}
185
186//
187// continued_fraction_a
188// Evaluates:
189//
190// a1
191// ---------------
192// b1 + a2
193// ----------
194// b2 + a3
195// -----
196// b3 + ...
197//
198// Note that the first a1 and b1 returned by generator Gen are both used.
199//
200template <typename Gen, typename U>
201inline typename detail::fraction_traits<Gen>::result_type continued_fraction_a(Gen& g, const U& factor, std::uintmax_t& max_terms)
202 noexcept(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
203{
204 BOOST_MATH_STD_USING // ADL of std names
205
206 using traits = detail::fraction_traits<Gen>;
207 using result_type = typename traits::result_type;
208 using value_type = typename traits::value_type;
209 using integer_type = typename integer_scalar_type<result_type>::type;
210 using scalar_type = typename scalar_type<result_type>::type;
211
212 integer_type const zero(0), one(1);
213
214 result_type tiny = detail::tiny_value<result_type>::get();
215 scalar_type terminator = abs(factor);
216
217 value_type v = g();
218
219 result_type f, C, D, delta, a0;
220 f = traits::b(v);
221 a0 = traits::a(v);
222 if(f == zero)
223 f = tiny;
224 C = f;
225 D = 0;
226
227 std::uintmax_t counter(max_terms);
228
229 do{
230 v = g();
231 D = traits::b(v) + traits::a(v) * D;
232 if(D == zero)
233 D = tiny;
234 C = traits::b(v) + traits::a(v) / C;
235 if(C == zero)
236 C = tiny;
237 D = one/D;
238 delta = C*D;
239 f = f * delta;
240 }while((abs(delta - one) > terminator) && --counter);
241
242 max_terms = max_terms - counter;
243
244 return a0/f;
245}
246
247template <typename Gen, typename U>
248inline typename detail::fraction_traits<Gen>::result_type continued_fraction_a(Gen& g, const U& factor)
249 noexcept(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
250{
251 std::uintmax_t max_iter = (std::numeric_limits<std::uintmax_t>::max)();
252 return continued_fraction_a(g, factor, max_iter);
253}
254
255template <typename Gen>
256inline typename detail::fraction_traits<Gen>::result_type continued_fraction_a(Gen& g, int bits)
257 noexcept(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
258{
259 BOOST_MATH_STD_USING // ADL of std names
260
261 typedef detail::fraction_traits<Gen> traits;
262 typedef typename traits::result_type result_type;
263
264 result_type factor = ldexp(1.0f, 1-bits); // 1 / pow(result_type(2), bits);
265 std::uintmax_t max_iter = (std::numeric_limits<std::uintmax_t>::max)();
266
267 return continued_fraction_a(g, factor, max_iter);
268}
269
270template <typename Gen>
271inline typename detail::fraction_traits<Gen>::result_type continued_fraction_a(Gen& g, int bits, std::uintmax_t& max_terms)
272 noexcept(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
273{
274 BOOST_MATH_STD_USING // ADL of std names
275
276 using traits = detail::fraction_traits<Gen>;
277 using result_type = typename traits::result_type;
278
279 result_type factor = ldexp(1.0f, 1-bits); // 1 / pow(result_type(2), bits);
280 return continued_fraction_a(g, factor, max_terms);
281}
282
283} // namespace tools
284} // namespace math
285} // namespace boost
286
287#endif // BOOST_MATH_TOOLS_FRACTION_INCLUDED
288

source code of boost/libs/math/include/boost/math/tools/fraction.hpp