123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287 |
- #ifndef BOOST_MATH_TOOLS_FRACTION_INCLUDED
- #define BOOST_MATH_TOOLS_FRACTION_INCLUDED
- #ifdef _MSC_VER
- #pragma once
- #endif
- #include <boost/math/tools/precision.hpp>
- #include <boost/math/tools/complex.hpp>
- #include <type_traits>
- #include <cstdint>
- #include <cmath>
- namespace boost{ namespace math{ namespace tools{
- namespace detail
- {
- template <typename T>
- struct is_pair : public std::false_type{};
- template <typename T, typename U>
- struct is_pair<std::pair<T,U>> : public std::true_type{};
- template <typename Gen>
- struct fraction_traits_simple
- {
- using result_type = typename Gen::result_type;
- using value_type = typename Gen::result_type;
- static result_type a(const value_type&) BOOST_MATH_NOEXCEPT(value_type)
- {
- return 1;
- }
- static result_type b(const value_type& v) BOOST_MATH_NOEXCEPT(value_type)
- {
- return v;
- }
- };
- template <typename Gen>
- struct fraction_traits_pair
- {
- using value_type = typename Gen::result_type;
- using result_type = typename value_type::first_type;
- static result_type a(const value_type& v) BOOST_MATH_NOEXCEPT(value_type)
- {
- return v.first;
- }
- static result_type b(const value_type& v) BOOST_MATH_NOEXCEPT(value_type)
- {
- return v.second;
- }
- };
- template <typename Gen>
- struct fraction_traits
- : public std::conditional<
- is_pair<typename Gen::result_type>::value,
- fraction_traits_pair<Gen>,
- fraction_traits_simple<Gen>>::type
- {
- };
- template <typename T, bool = is_complex_type<T>::value>
- struct tiny_value
- {
-
-
-
- static T get() {
- return 16*tools::min_value<T>();
- }
- };
- template <typename T>
- struct tiny_value<T, true>
- {
- using value_type = typename T::value_type;
- static T get() {
- return 16*tools::min_value<value_type>();
- }
- };
- }
- template <typename Gen, typename U>
- inline typename detail::fraction_traits<Gen>::result_type continued_fraction_b(Gen& g, const U& factor, std::uintmax_t& max_terms)
- BOOST_NOEXCEPT_IF(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
- {
- BOOST_MATH_STD_USING
- using traits = detail::fraction_traits<Gen>;
- using result_type = typename traits::result_type;
- using value_type = typename traits::value_type;
- using integer_type = typename integer_scalar_type<result_type>::type;
- using scalar_type = typename scalar_type<result_type>::type;
- integer_type const zero(0), one(1);
- result_type tiny = detail::tiny_value<result_type>::get();
- scalar_type terminator = abs(factor);
- value_type v = g();
- result_type f, C, D, delta;
- f = traits::b(v);
- if(f == zero)
- f = tiny;
- C = f;
- D = 0;
- std::uintmax_t counter(max_terms);
- do{
- v = g();
- D = traits::b(v) + traits::a(v) * D;
- if(D == result_type(0))
- D = tiny;
- C = traits::b(v) + traits::a(v) / C;
- if(C == zero)
- C = tiny;
- D = one/D;
- delta = C*D;
- f = f * delta;
- }while((abs(delta - one) > terminator) && --counter);
- max_terms = max_terms - counter;
- return f;
- }
- template <typename Gen, typename U>
- inline typename detail::fraction_traits<Gen>::result_type continued_fraction_b(Gen& g, const U& factor)
- BOOST_NOEXCEPT_IF(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
- {
- std::uintmax_t max_terms = (std::numeric_limits<std::uintmax_t>::max)();
- return continued_fraction_b(g, factor, max_terms);
- }
- template <typename Gen>
- inline typename detail::fraction_traits<Gen>::result_type continued_fraction_b(Gen& g, int bits)
- BOOST_NOEXCEPT_IF(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
- {
- BOOST_MATH_STD_USING
- using traits = detail::fraction_traits<Gen>;
- using result_type = typename traits::result_type;
- result_type factor = ldexp(1.0f, 1 - bits);
- std::uintmax_t max_terms = (std::numeric_limits<std::uintmax_t>::max)();
- return continued_fraction_b(g, factor, max_terms);
- }
- template <typename Gen>
- inline typename detail::fraction_traits<Gen>::result_type continued_fraction_b(Gen& g, int bits, std::uintmax_t& max_terms)
- BOOST_NOEXCEPT_IF(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
- {
- BOOST_MATH_STD_USING
- using traits = detail::fraction_traits<Gen>;
- using result_type = typename traits::result_type;
- result_type factor = ldexp(1.0f, 1 - bits);
- return continued_fraction_b(g, factor, max_terms);
- }
- template <typename Gen, typename U>
- inline typename detail::fraction_traits<Gen>::result_type continued_fraction_a(Gen& g, const U& factor, std::uintmax_t& max_terms)
- BOOST_NOEXCEPT_IF(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
- {
- BOOST_MATH_STD_USING
- using traits = detail::fraction_traits<Gen>;
- using result_type = typename traits::result_type;
- using value_type = typename traits::value_type;
- using integer_type = typename integer_scalar_type<result_type>::type;
- using scalar_type = typename scalar_type<result_type>::type;
- integer_type const zero(0), one(1);
- result_type tiny = detail::tiny_value<result_type>::get();
- scalar_type terminator = abs(factor);
- value_type v = g();
- result_type f, C, D, delta, a0;
- f = traits::b(v);
- a0 = traits::a(v);
- if(f == zero)
- f = tiny;
- C = f;
- D = 0;
- std::uintmax_t counter(max_terms);
- do{
- v = g();
- D = traits::b(v) + traits::a(v) * D;
- if(D == zero)
- D = tiny;
- C = traits::b(v) + traits::a(v) / C;
- if(C == zero)
- C = tiny;
- D = one/D;
- delta = C*D;
- f = f * delta;
- }while((abs(delta - one) > terminator) && --counter);
- max_terms = max_terms - counter;
- return a0/f;
- }
- template <typename Gen, typename U>
- inline typename detail::fraction_traits<Gen>::result_type continued_fraction_a(Gen& g, const U& factor)
- BOOST_NOEXCEPT_IF(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
- {
- std::uintmax_t max_iter = (std::numeric_limits<std::uintmax_t>::max)();
- return continued_fraction_a(g, factor, max_iter);
- }
- template <typename Gen>
- inline typename detail::fraction_traits<Gen>::result_type continued_fraction_a(Gen& g, int bits)
- BOOST_NOEXCEPT_IF(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
- {
- BOOST_MATH_STD_USING
- typedef detail::fraction_traits<Gen> traits;
- typedef typename traits::result_type result_type;
- result_type factor = ldexp(1.0f, 1-bits);
- std::uintmax_t max_iter = (std::numeric_limits<std::uintmax_t>::max)();
- return continued_fraction_a(g, factor, max_iter);
- }
- template <typename Gen>
- inline typename detail::fraction_traits<Gen>::result_type continued_fraction_a(Gen& g, int bits, std::uintmax_t& max_terms)
- BOOST_NOEXCEPT_IF(BOOST_MATH_IS_FLOAT(typename detail::fraction_traits<Gen>::result_type) && noexcept(std::declval<Gen>()()))
- {
- BOOST_MATH_STD_USING
- using traits = detail::fraction_traits<Gen>;
- using result_type = typename traits::result_type;
- result_type factor = ldexp(1.0f, 1-bits);
- return continued_fraction_a(g, factor, max_terms);
- }
- }
- }
- }
- #endif
|