math/special_functions/bessel_iterators.hpp

0001
0002 ///////////////////////////////////////////////////////////////////////////////
0003 //  Copyright 2018 John Maddock
0004 //  Distributed under the Boost
0005 //  Software License, Version 1.0. (See accompanying file
0006 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
0007
0008 #ifndef BOOST_MATH_BESSEL_ITERATORS_HPP
0009 #define BOOST_MATH_BESSEL_ITERATORS_HPP
0010
0011 #include <boost/math/tools/recurrence.hpp>
0012
0013 namespace boost {
0014    namespace math {
0015       namespace detail {
0016
0017          template <class T>
0018          struct bessel_jy_recurrence
0019          {
0020             bessel_jy_recurrence(T v, T z) : v(v), z(z) {}
0021             boost::math::tuple<T, T, T> operator()(int k)
0022             {
0023                return boost::math::tuple<T, T, T>(1, -2 * (v + k) / z, 1);
0024             }
0025
0026             T v, z;
0027          };
0028          template <class T>
0029          struct bessel_ik_recurrence
0030          {
0031             bessel_ik_recurrence(T v, T z) : v(v), z(z) {}
0032             boost::math::tuple<T, T, T> operator()(int k)
0033             {
0034                return boost::math::tuple<T, T, T>(1, -2 * (v + k) / z, -1);
0035             }
0036
0037             T v, z;
0038          };
0039       } // namespace detail
0040
0041       template <class T, class Policy = boost::math::policies::policy<> >
0042       struct bessel_j_backwards_iterator
0043       {
0044          typedef std::ptrdiff_t difference_type;
0045          typedef T value_type;
0046          typedef T* pointer;
0047          typedef T& reference;
0048          typedef std::input_iterator_tag iterator_category;
0049
0050          bessel_j_backwards_iterator(const T& v, const T& x)
0051             : it(detail::bessel_jy_recurrence<T>(v, x), boost::math::cyl_bessel_j(v, x, Policy()))
0052          {
0053             if(v < 0)
0054                boost::math::policies::raise_domain_error("bessel_j_backwards_iterator<%1%>", "Order must be > 0 stable backwards recurrence but got %1%", v, Policy());
0055          }
0056
0057          bessel_j_backwards_iterator(const T& v, const T& x, const T& J_v)
0058             : it(detail::bessel_jy_recurrence<T>(v, x), J_v)
0059          {
0060             if(v < 0)
0061                boost::math::policies::raise_domain_error("bessel_j_backwards_iterator<%1%>", "Order must be > 0 stable backwards recurrence but got %1%", v, Policy());
0062          }
0063          bessel_j_backwards_iterator(const T& v, const T& x, const T& J_v_plus_1, const T& J_v)
0064             : it(detail::bessel_jy_recurrence<T>(v, x), J_v_plus_1, J_v)
0065          {
0066             if (v < -1)
0067                boost::math::policies::raise_domain_error("bessel_j_backwards_iterator<%1%>", "Order must be > 0 stable backwards recurrence but got %1%", v, Policy());
0068          }
0069
0070          bessel_j_backwards_iterator& operator++()
0071          {
0072             ++it;
0073             return *this;
0074          }
0075
0076          bessel_j_backwards_iterator operator++(int)
0077          {
0078             bessel_j_backwards_iterator t(*this);
0079             ++(*this);
0080             return t;
0081          }
0082
0083          T operator*() { return *it; }
0084
0085       private:
0086          boost::math::tools::backward_recurrence_iterator< detail::bessel_jy_recurrence<T> > it;
0087       };
0088
0089       template <class T, class Policy = boost::math::policies::policy<> >
0090       struct bessel_i_backwards_iterator
0091       {
0092          typedef std::ptrdiff_t difference_type;
0093          typedef T value_type;
0094          typedef T* pointer;
0095          typedef T& reference;
0096          typedef std::input_iterator_tag iterator_category;
0097
0098          bessel_i_backwards_iterator(const T& v, const T& x)
0099             : it(detail::bessel_ik_recurrence<T>(v, x), boost::math::cyl_bessel_i(v, x, Policy()))
0100          {
0101             if(v < -1)
0102                boost::math::policies::raise_domain_error("bessel_i_backwards_iterator<%1%>", "Order must be > 0 stable backwards recurrence but got %1%", v, Policy());
0103          }
0104          bessel_i_backwards_iterator(const T& v, const T& x, const T& I_v)
0105             : it(detail::bessel_ik_recurrence<T>(v, x), I_v)
0106          {
0107             if(v < -1)
0108                boost::math::policies::raise_domain_error("bessel_i_backwards_iterator<%1%>", "Order must be > 0 stable backwards recurrence but got %1%", v, Policy());
0109          }
0110          bessel_i_backwards_iterator(const T& v, const T& x, const T& I_v_plus_1, const T& I_v)
0111             : it(detail::bessel_ik_recurrence<T>(v, x), I_v_plus_1, I_v)
0112          {
0113             if(v < -1)
0114                boost::math::policies::raise_domain_error("bessel_i_backwards_iterator<%1%>", "Order must be > 0 stable backwards recurrence but got %1%", v, Policy());
0115          }
0116
0117          bessel_i_backwards_iterator& operator++()
0118          {
0119             ++it;
0120             return *this;
0121          }
0122
0123          bessel_i_backwards_iterator operator++(int)
0124          {
0125             bessel_i_backwards_iterator t(*this);
0126             ++(*this);
0127             return t;
0128          }
0129
0130          T operator*() { return *it; }
0131
0132       private:
0133          boost::math::tools::backward_recurrence_iterator< detail::bessel_ik_recurrence<T> > it;
0134       };
0135
0136       template <class T, class Policy = boost::math::policies::policy<> >
0137       struct bessel_i_forwards_iterator
0138       {
0139          typedef std::ptrdiff_t difference_type;
0140          typedef T value_type;
0141          typedef T* pointer;
0142          typedef T& reference;
0143          typedef std::input_iterator_tag iterator_category;
0144
0145          bessel_i_forwards_iterator(const T& v, const T& x)
0146             : it(detail::bessel_ik_recurrence<T>(v, x), boost::math::cyl_bessel_i(v, x, Policy()))
0147          {
0148             if(v > 1)
0149                boost::math::policies::raise_domain_error("bessel_i_forwards_iterator<%1%>", "Order must be < 0 stable forwards recurrence but got %1%", v, Policy());
0150          }
0151          bessel_i_forwards_iterator(const T& v, const T& x, const T& I_v)
0152             : it(detail::bessel_ik_recurrence<T>(v, x), I_v)
0153          {
0154             if (v > 1)
0155                boost::math::policies::raise_domain_error("bessel_i_forwards_iterator<%1%>", "Order must be < 0 stable forwards recurrence but got %1%", v, Policy());
0156          }
0157          bessel_i_forwards_iterator(const T& v, const T& x, const T& I_v_plus_1, const T& I_v)
0158             : it(detail::bessel_ik_recurrence<T>(v, x), I_v_plus_1, I_v)
0159          {
0160             if (v > 1)
0161                boost::math::policies::raise_domain_error("bessel_i_forwards_iterator<%1%>", "Order must be < 0 stable forwards recurrence but got %1%", v, Policy());
0162          }
0163
0164          bessel_i_forwards_iterator& operator++()
0165          {
0166             ++it;
0167             return *this;
0168          }
0169
0170          bessel_i_forwards_iterator operator++(int)
0171          {
0172             bessel_i_forwards_iterator t(*this);
0173             ++(*this);
0174             return t;
0175          }
0176
0177          T operator*() { return *it; }
0178
0179       private:
0180          boost::math::tools::forward_recurrence_iterator< detail::bessel_ik_recurrence<T> > it;
0181       };
0182
0183    }
0184 } // namespaces
0185
0186 #endif // BOOST_MATH_BESSEL_ITERATORS_HPP