annotate DEPENDENCIES/generic/include/boost/math/tools/rational.hpp @ 133:4acb5d8d80b6 tip

Don't fail environmental check if README.md exists (but .txt and no-suffix don't)
author Chris Cannam
date Tue, 30 Jul 2019 12:25:44 +0100
parents 2665513ce2d3
children
rev   line source
Chris@16 1 // (C) Copyright John Maddock 2006.
Chris@16 2 // Use, modification and distribution are subject to the
Chris@16 3 // Boost Software License, Version 1.0. (See accompanying file
Chris@16 4 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
Chris@16 5
Chris@16 6 #ifndef BOOST_MATH_TOOLS_RATIONAL_HPP
Chris@16 7 #define BOOST_MATH_TOOLS_RATIONAL_HPP
Chris@16 8
Chris@16 9 #ifdef _MSC_VER
Chris@16 10 #pragma once
Chris@16 11 #endif
Chris@16 12
Chris@16 13 #include <boost/array.hpp>
Chris@16 14 #include <boost/math/tools/config.hpp>
Chris@16 15 #include <boost/mpl/int.hpp>
Chris@16 16
Chris@16 17 #if BOOST_MATH_POLY_METHOD == 1
Chris@16 18 # define BOOST_HEADER() <BOOST_JOIN(boost/math/tools/detail/polynomial_horner1_, BOOST_MATH_MAX_POLY_ORDER).hpp>
Chris@16 19 # include BOOST_HEADER()
Chris@16 20 # undef BOOST_HEADER
Chris@16 21 #elif BOOST_MATH_POLY_METHOD == 2
Chris@16 22 # define BOOST_HEADER() <BOOST_JOIN(boost/math/tools/detail/polynomial_horner2_, BOOST_MATH_MAX_POLY_ORDER).hpp>
Chris@16 23 # include BOOST_HEADER()
Chris@16 24 # undef BOOST_HEADER
Chris@16 25 #elif BOOST_MATH_POLY_METHOD == 3
Chris@16 26 # define BOOST_HEADER() <BOOST_JOIN(boost/math/tools/detail/polynomial_horner3_, BOOST_MATH_MAX_POLY_ORDER).hpp>
Chris@16 27 # include BOOST_HEADER()
Chris@16 28 # undef BOOST_HEADER
Chris@16 29 #endif
Chris@16 30 #if BOOST_MATH_RATIONAL_METHOD == 1
Chris@16 31 # define BOOST_HEADER() <BOOST_JOIN(boost/math/tools/detail/rational_horner1_, BOOST_MATH_MAX_POLY_ORDER).hpp>
Chris@16 32 # include BOOST_HEADER()
Chris@16 33 # undef BOOST_HEADER
Chris@16 34 #elif BOOST_MATH_RATIONAL_METHOD == 2
Chris@16 35 # define BOOST_HEADER() <BOOST_JOIN(boost/math/tools/detail/rational_horner2_, BOOST_MATH_MAX_POLY_ORDER).hpp>
Chris@16 36 # include BOOST_HEADER()
Chris@16 37 # undef BOOST_HEADER
Chris@16 38 #elif BOOST_MATH_RATIONAL_METHOD == 3
Chris@16 39 # define BOOST_HEADER() <BOOST_JOIN(boost/math/tools/detail/rational_horner3_, BOOST_MATH_MAX_POLY_ORDER).hpp>
Chris@16 40 # include BOOST_HEADER()
Chris@16 41 # undef BOOST_HEADER
Chris@16 42 #endif
Chris@16 43
Chris@16 44 #if 0
Chris@16 45 //
Chris@16 46 // This just allows dependency trackers to find the headers
Chris@16 47 // used in the above PP-magic.
Chris@16 48 //
Chris@16 49 #include <boost/math/tools/detail/polynomial_horner1_2.hpp>
Chris@16 50 #include <boost/math/tools/detail/polynomial_horner1_3.hpp>
Chris@16 51 #include <boost/math/tools/detail/polynomial_horner1_4.hpp>
Chris@16 52 #include <boost/math/tools/detail/polynomial_horner1_5.hpp>
Chris@16 53 #include <boost/math/tools/detail/polynomial_horner1_6.hpp>
Chris@16 54 #include <boost/math/tools/detail/polynomial_horner1_7.hpp>
Chris@16 55 #include <boost/math/tools/detail/polynomial_horner1_8.hpp>
Chris@16 56 #include <boost/math/tools/detail/polynomial_horner1_9.hpp>
Chris@16 57 #include <boost/math/tools/detail/polynomial_horner1_10.hpp>
Chris@16 58 #include <boost/math/tools/detail/polynomial_horner1_11.hpp>
Chris@16 59 #include <boost/math/tools/detail/polynomial_horner1_12.hpp>
Chris@16 60 #include <boost/math/tools/detail/polynomial_horner1_13.hpp>
Chris@16 61 #include <boost/math/tools/detail/polynomial_horner1_14.hpp>
Chris@16 62 #include <boost/math/tools/detail/polynomial_horner1_15.hpp>
Chris@16 63 #include <boost/math/tools/detail/polynomial_horner1_16.hpp>
Chris@16 64 #include <boost/math/tools/detail/polynomial_horner1_17.hpp>
Chris@16 65 #include <boost/math/tools/detail/polynomial_horner1_18.hpp>
Chris@16 66 #include <boost/math/tools/detail/polynomial_horner1_19.hpp>
Chris@16 67 #include <boost/math/tools/detail/polynomial_horner1_20.hpp>
Chris@16 68 #include <boost/math/tools/detail/polynomial_horner2_2.hpp>
Chris@16 69 #include <boost/math/tools/detail/polynomial_horner2_3.hpp>
Chris@16 70 #include <boost/math/tools/detail/polynomial_horner2_4.hpp>
Chris@16 71 #include <boost/math/tools/detail/polynomial_horner2_5.hpp>
Chris@16 72 #include <boost/math/tools/detail/polynomial_horner2_6.hpp>
Chris@16 73 #include <boost/math/tools/detail/polynomial_horner2_7.hpp>
Chris@16 74 #include <boost/math/tools/detail/polynomial_horner2_8.hpp>
Chris@16 75 #include <boost/math/tools/detail/polynomial_horner2_9.hpp>
Chris@16 76 #include <boost/math/tools/detail/polynomial_horner2_10.hpp>
Chris@16 77 #include <boost/math/tools/detail/polynomial_horner2_11.hpp>
Chris@16 78 #include <boost/math/tools/detail/polynomial_horner2_12.hpp>
Chris@16 79 #include <boost/math/tools/detail/polynomial_horner2_13.hpp>
Chris@16 80 #include <boost/math/tools/detail/polynomial_horner2_14.hpp>
Chris@16 81 #include <boost/math/tools/detail/polynomial_horner2_15.hpp>
Chris@16 82 #include <boost/math/tools/detail/polynomial_horner2_16.hpp>
Chris@16 83 #include <boost/math/tools/detail/polynomial_horner2_17.hpp>
Chris@16 84 #include <boost/math/tools/detail/polynomial_horner2_18.hpp>
Chris@16 85 #include <boost/math/tools/detail/polynomial_horner2_19.hpp>
Chris@16 86 #include <boost/math/tools/detail/polynomial_horner2_20.hpp>
Chris@16 87 #include <boost/math/tools/detail/polynomial_horner3_2.hpp>
Chris@16 88 #include <boost/math/tools/detail/polynomial_horner3_3.hpp>
Chris@16 89 #include <boost/math/tools/detail/polynomial_horner3_4.hpp>
Chris@16 90 #include <boost/math/tools/detail/polynomial_horner3_5.hpp>
Chris@16 91 #include <boost/math/tools/detail/polynomial_horner3_6.hpp>
Chris@16 92 #include <boost/math/tools/detail/polynomial_horner3_7.hpp>
Chris@16 93 #include <boost/math/tools/detail/polynomial_horner3_8.hpp>
Chris@16 94 #include <boost/math/tools/detail/polynomial_horner3_9.hpp>
Chris@16 95 #include <boost/math/tools/detail/polynomial_horner3_10.hpp>
Chris@16 96 #include <boost/math/tools/detail/polynomial_horner3_11.hpp>
Chris@16 97 #include <boost/math/tools/detail/polynomial_horner3_12.hpp>
Chris@16 98 #include <boost/math/tools/detail/polynomial_horner3_13.hpp>
Chris@16 99 #include <boost/math/tools/detail/polynomial_horner3_14.hpp>
Chris@16 100 #include <boost/math/tools/detail/polynomial_horner3_15.hpp>
Chris@16 101 #include <boost/math/tools/detail/polynomial_horner3_16.hpp>
Chris@16 102 #include <boost/math/tools/detail/polynomial_horner3_17.hpp>
Chris@16 103 #include <boost/math/tools/detail/polynomial_horner3_18.hpp>
Chris@16 104 #include <boost/math/tools/detail/polynomial_horner3_19.hpp>
Chris@16 105 #include <boost/math/tools/detail/polynomial_horner3_20.hpp>
Chris@16 106 #include <boost/math/tools/detail/rational_horner1_2.hpp>
Chris@16 107 #include <boost/math/tools/detail/rational_horner1_3.hpp>
Chris@16 108 #include <boost/math/tools/detail/rational_horner1_4.hpp>
Chris@16 109 #include <boost/math/tools/detail/rational_horner1_5.hpp>
Chris@16 110 #include <boost/math/tools/detail/rational_horner1_6.hpp>
Chris@16 111 #include <boost/math/tools/detail/rational_horner1_7.hpp>
Chris@16 112 #include <boost/math/tools/detail/rational_horner1_8.hpp>
Chris@16 113 #include <boost/math/tools/detail/rational_horner1_9.hpp>
Chris@16 114 #include <boost/math/tools/detail/rational_horner1_10.hpp>
Chris@16 115 #include <boost/math/tools/detail/rational_horner1_11.hpp>
Chris@16 116 #include <boost/math/tools/detail/rational_horner1_12.hpp>
Chris@16 117 #include <boost/math/tools/detail/rational_horner1_13.hpp>
Chris@16 118 #include <boost/math/tools/detail/rational_horner1_14.hpp>
Chris@16 119 #include <boost/math/tools/detail/rational_horner1_15.hpp>
Chris@16 120 #include <boost/math/tools/detail/rational_horner1_16.hpp>
Chris@16 121 #include <boost/math/tools/detail/rational_horner1_17.hpp>
Chris@16 122 #include <boost/math/tools/detail/rational_horner1_18.hpp>
Chris@16 123 #include <boost/math/tools/detail/rational_horner1_19.hpp>
Chris@16 124 #include <boost/math/tools/detail/rational_horner1_20.hpp>
Chris@16 125 #include <boost/math/tools/detail/rational_horner2_2.hpp>
Chris@16 126 #include <boost/math/tools/detail/rational_horner2_3.hpp>
Chris@16 127 #include <boost/math/tools/detail/rational_horner2_4.hpp>
Chris@16 128 #include <boost/math/tools/detail/rational_horner2_5.hpp>
Chris@16 129 #include <boost/math/tools/detail/rational_horner2_6.hpp>
Chris@16 130 #include <boost/math/tools/detail/rational_horner2_7.hpp>
Chris@16 131 #include <boost/math/tools/detail/rational_horner2_8.hpp>
Chris@16 132 #include <boost/math/tools/detail/rational_horner2_9.hpp>
Chris@16 133 #include <boost/math/tools/detail/rational_horner2_10.hpp>
Chris@16 134 #include <boost/math/tools/detail/rational_horner2_11.hpp>
Chris@16 135 #include <boost/math/tools/detail/rational_horner2_12.hpp>
Chris@16 136 #include <boost/math/tools/detail/rational_horner2_13.hpp>
Chris@16 137 #include <boost/math/tools/detail/rational_horner2_14.hpp>
Chris@16 138 #include <boost/math/tools/detail/rational_horner2_15.hpp>
Chris@16 139 #include <boost/math/tools/detail/rational_horner2_16.hpp>
Chris@16 140 #include <boost/math/tools/detail/rational_horner2_17.hpp>
Chris@16 141 #include <boost/math/tools/detail/rational_horner2_18.hpp>
Chris@16 142 #include <boost/math/tools/detail/rational_horner2_19.hpp>
Chris@16 143 #include <boost/math/tools/detail/rational_horner2_20.hpp>
Chris@16 144 #include <boost/math/tools/detail/rational_horner3_2.hpp>
Chris@16 145 #include <boost/math/tools/detail/rational_horner3_3.hpp>
Chris@16 146 #include <boost/math/tools/detail/rational_horner3_4.hpp>
Chris@16 147 #include <boost/math/tools/detail/rational_horner3_5.hpp>
Chris@16 148 #include <boost/math/tools/detail/rational_horner3_6.hpp>
Chris@16 149 #include <boost/math/tools/detail/rational_horner3_7.hpp>
Chris@16 150 #include <boost/math/tools/detail/rational_horner3_8.hpp>
Chris@16 151 #include <boost/math/tools/detail/rational_horner3_9.hpp>
Chris@16 152 #include <boost/math/tools/detail/rational_horner3_10.hpp>
Chris@16 153 #include <boost/math/tools/detail/rational_horner3_11.hpp>
Chris@16 154 #include <boost/math/tools/detail/rational_horner3_12.hpp>
Chris@16 155 #include <boost/math/tools/detail/rational_horner3_13.hpp>
Chris@16 156 #include <boost/math/tools/detail/rational_horner3_14.hpp>
Chris@16 157 #include <boost/math/tools/detail/rational_horner3_15.hpp>
Chris@16 158 #include <boost/math/tools/detail/rational_horner3_16.hpp>
Chris@16 159 #include <boost/math/tools/detail/rational_horner3_17.hpp>
Chris@16 160 #include <boost/math/tools/detail/rational_horner3_18.hpp>
Chris@16 161 #include <boost/math/tools/detail/rational_horner3_19.hpp>
Chris@16 162 #include <boost/math/tools/detail/rational_horner3_20.hpp>
Chris@16 163 #endif
Chris@16 164
Chris@16 165 namespace boost{ namespace math{ namespace tools{
Chris@16 166
Chris@16 167 //
Chris@16 168 // Forward declaration to keep two phase lookup happy:
Chris@16 169 //
Chris@16 170 template <class T, class U>
Chris@16 171 U evaluate_polynomial(const T* poly, U const& z, std::size_t count);
Chris@16 172
Chris@16 173 namespace detail{
Chris@16 174
Chris@16 175 template <class T, class V, class Tag>
Chris@16 176 inline V evaluate_polynomial_c_imp(const T* a, const V& val, const Tag*)
Chris@16 177 {
Chris@16 178 return evaluate_polynomial(a, val, Tag::value);
Chris@16 179 }
Chris@16 180
Chris@16 181 } // namespace detail
Chris@16 182
Chris@16 183 //
Chris@16 184 // Polynomial evaluation with runtime size.
Chris@16 185 // This requires a for-loop which may be more expensive than
Chris@16 186 // the loop expanded versions above:
Chris@16 187 //
Chris@16 188 template <class T, class U>
Chris@16 189 inline U evaluate_polynomial(const T* poly, U const& z, std::size_t count)
Chris@16 190 {
Chris@16 191 BOOST_ASSERT(count > 0);
Chris@16 192 U sum = static_cast<U>(poly[count - 1]);
Chris@16 193 for(int i = static_cast<int>(count) - 2; i >= 0; --i)
Chris@16 194 {
Chris@16 195 sum *= z;
Chris@16 196 sum += static_cast<U>(poly[i]);
Chris@16 197 }
Chris@16 198 return sum;
Chris@16 199 }
Chris@16 200 //
Chris@16 201 // Compile time sized polynomials, just inline forwarders to the
Chris@16 202 // implementations above:
Chris@16 203 //
Chris@16 204 template <std::size_t N, class T, class V>
Chris@16 205 inline V evaluate_polynomial(const T(&a)[N], const V& val)
Chris@16 206 {
Chris@16 207 typedef mpl::int_<N> tag_type;
Chris@16 208 return detail::evaluate_polynomial_c_imp(static_cast<const T*>(a), val, static_cast<tag_type const*>(0));
Chris@16 209 }
Chris@16 210
Chris@16 211 template <std::size_t N, class T, class V>
Chris@16 212 inline V evaluate_polynomial(const boost::array<T,N>& a, const V& val)
Chris@16 213 {
Chris@16 214 typedef mpl::int_<N> tag_type;
Chris@16 215 return detail::evaluate_polynomial_c_imp(static_cast<const T*>(a.data()), val, static_cast<tag_type const*>(0));
Chris@16 216 }
Chris@16 217 //
Chris@16 218 // Even polynomials are trivial: just square the argument!
Chris@16 219 //
Chris@16 220 template <class T, class U>
Chris@16 221 inline U evaluate_even_polynomial(const T* poly, U z, std::size_t count)
Chris@16 222 {
Chris@16 223 return evaluate_polynomial(poly, U(z*z), count);
Chris@16 224 }
Chris@16 225
Chris@16 226 template <std::size_t N, class T, class V>
Chris@16 227 inline V evaluate_even_polynomial(const T(&a)[N], const V& z)
Chris@16 228 {
Chris@16 229 return evaluate_polynomial(a, V(z*z));
Chris@16 230 }
Chris@16 231
Chris@16 232 template <std::size_t N, class T, class V>
Chris@16 233 inline V evaluate_even_polynomial(const boost::array<T,N>& a, const V& z)
Chris@16 234 {
Chris@16 235 return evaluate_polynomial(a, V(z*z));
Chris@16 236 }
Chris@16 237 //
Chris@16 238 // Odd polynomials come next:
Chris@16 239 //
Chris@16 240 template <class T, class U>
Chris@16 241 inline U evaluate_odd_polynomial(const T* poly, U z, std::size_t count)
Chris@16 242 {
Chris@16 243 return poly[0] + z * evaluate_polynomial(poly+1, U(z*z), count-1);
Chris@16 244 }
Chris@16 245
Chris@16 246 template <std::size_t N, class T, class V>
Chris@16 247 inline V evaluate_odd_polynomial(const T(&a)[N], const V& z)
Chris@16 248 {
Chris@16 249 typedef mpl::int_<N-1> tag_type;
Chris@16 250 return a[0] + z * detail::evaluate_polynomial_c_imp(static_cast<const T*>(a) + 1, V(z*z), static_cast<tag_type const*>(0));
Chris@16 251 }
Chris@16 252
Chris@16 253 template <std::size_t N, class T, class V>
Chris@16 254 inline V evaluate_odd_polynomial(const boost::array<T,N>& a, const V& z)
Chris@16 255 {
Chris@16 256 typedef mpl::int_<N-1> tag_type;
Chris@16 257 return a[0] + z * detail::evaluate_polynomial_c_imp(static_cast<const T*>(a.data()) + 1, V(z*z), static_cast<tag_type const*>(0));
Chris@16 258 }
Chris@16 259
Chris@16 260 template <class T, class U, class V>
Chris@16 261 V evaluate_rational(const T* num, const U* denom, const V& z_, std::size_t count);
Chris@16 262
Chris@16 263 namespace detail{
Chris@16 264
Chris@16 265 template <class T, class U, class V, class Tag>
Chris@16 266 inline V evaluate_rational_c_imp(const T* num, const U* denom, const V& z, const Tag*)
Chris@16 267 {
Chris@16 268 return boost::math::tools::evaluate_rational(num, denom, z, Tag::value);
Chris@16 269 }
Chris@16 270
Chris@16 271 }
Chris@16 272 //
Chris@16 273 // Rational functions: numerator and denominator must be
Chris@16 274 // equal in size. These always have a for-loop and so may be less
Chris@16 275 // efficient than evaluating a pair of polynomials. However, there
Chris@16 276 // are some tricks we can use to prevent overflow that might otherwise
Chris@16 277 // occur in polynomial evaluation, if z is large. This is important
Chris@16 278 // in our Lanczos code for example.
Chris@16 279 //
Chris@16 280 template <class T, class U, class V>
Chris@16 281 V evaluate_rational(const T* num, const U* denom, const V& z_, std::size_t count)
Chris@16 282 {
Chris@16 283 V z(z_);
Chris@16 284 V s1, s2;
Chris@16 285 if(z <= 1)
Chris@16 286 {
Chris@16 287 s1 = static_cast<V>(num[count-1]);
Chris@16 288 s2 = static_cast<V>(denom[count-1]);
Chris@16 289 for(int i = (int)count - 2; i >= 0; --i)
Chris@16 290 {
Chris@16 291 s1 *= z;
Chris@16 292 s2 *= z;
Chris@16 293 s1 += num[i];
Chris@16 294 s2 += denom[i];
Chris@16 295 }
Chris@16 296 }
Chris@16 297 else
Chris@16 298 {
Chris@16 299 z = 1 / z;
Chris@16 300 s1 = static_cast<V>(num[0]);
Chris@16 301 s2 = static_cast<V>(denom[0]);
Chris@16 302 for(unsigned i = 1; i < count; ++i)
Chris@16 303 {
Chris@16 304 s1 *= z;
Chris@16 305 s2 *= z;
Chris@16 306 s1 += num[i];
Chris@16 307 s2 += denom[i];
Chris@16 308 }
Chris@16 309 }
Chris@16 310 return s1 / s2;
Chris@16 311 }
Chris@16 312
Chris@16 313 template <std::size_t N, class T, class U, class V>
Chris@16 314 inline V evaluate_rational(const T(&a)[N], const U(&b)[N], const V& z)
Chris@16 315 {
Chris@16 316 return detail::evaluate_rational_c_imp(a, b, z, static_cast<const mpl::int_<N>*>(0));
Chris@16 317 }
Chris@16 318
Chris@16 319 template <std::size_t N, class T, class U, class V>
Chris@16 320 inline V evaluate_rational(const boost::array<T,N>& a, const boost::array<U,N>& b, const V& z)
Chris@16 321 {
Chris@16 322 return detail::evaluate_rational_c_imp(a.data(), b.data(), z, static_cast<mpl::int_<N>*>(0));
Chris@16 323 }
Chris@16 324
Chris@16 325 } // namespace tools
Chris@16 326 } // namespace math
Chris@16 327 } // namespace boost
Chris@16 328
Chris@16 329 #endif // BOOST_MATH_TOOLS_RATIONAL_HPP
Chris@16 330
Chris@16 331
Chris@16 332
Chris@16 333