ceph/src/boost/libs/math/example/barycentric_interpolation_example_2.cpp

   1
   2 // Copyright Nick Thompson, 2017
   3
   4 // Distributed under the Boost Software License, Version 1.0.
   5 // (See accompanying file LICENSE_1_0.txt or
   6 // copy at http://www.boost.org/LICENSE_1_0.txt).
   7
   8 #include <iostream>
   9 #include <limits>
  10 #include <map>
  11
  12 //[barycentric_rational_example2
  13
  14 /*`This further example shows how to use the iterator based constructor, and then uses the
  15 function object in our root finding algorithms to locate the points where the potential
  16 achieves a specific value.
  17 */
  18
  19 #include <boost/math/interpolators/barycentric_rational.hpp>
  20 #include <boost/range/adaptors.hpp>
  21 #include <boost/math/tools/roots.hpp>
  22
  23 int main()
  24 {
  25     // The lithium potential is given in Kohn's paper, Table I,
  26     // we could equally use an unordered_map, a list of tuples or pairs,
  27     // or a 2-dimentional array equally easily:
  28     std::map<double, double> r;
  29
  30     r[0.02] = 5.727;
  31     r[0.04] = 5.544;
  32     r[0.06] = 5.450;
  33     r[0.08] = 5.351;
  34     r[0.10] = 5.253;
  35     r[0.12] = 5.157;
  36     r[0.14] = 5.058;
  37     r[0.16] = 4.960;
  38     r[0.18] = 4.862;
  39     r[0.20] = 4.762;
  40     r[0.24] = 4.563;
  41     r[0.28] = 4.360;
  42     r[0.32] = 4.1584;
  43     r[0.36] = 3.9463;
  44     r[0.40] = 3.7360;
  45     r[0.44] = 3.5429;
  46     r[0.48] = 3.3797;
  47     r[0.52] = 3.2417;
  48     r[0.56] = 3.1209;
  49     r[0.60] = 3.0138;
  50     r[0.68] = 2.8342;
  51     r[0.76] = 2.6881;
  52     r[0.84] = 2.5662;
  53     r[0.92] = 2.4242;
  54     r[1.00] = 2.3766;
  55     r[1.08] = 2.3058;
  56     r[1.16] = 2.2458;
  57     r[1.24] = 2.2035;
  58     r[1.32] = 2.1661;
  59     r[1.40] = 2.1350;
  60     r[1.48] = 2.1090;
  61     r[1.64] = 2.0697;
  62     r[1.80] = 2.0466;
  63     r[1.96] = 2.0325;
  64     r[2.12] = 2.0288;
  65     r[2.28] = 2.0292;
  66     r[2.44] = 2.0228;
  67     r[2.60] = 2.0124;
  68     r[2.76] = 2.0065;
  69     r[2.92] = 2.0031;
  70     r[3.08] = 2.0015;
  71     r[3.24] = 2.0008;
  72     r[3.40] = 2.0004;
  73     r[3.56] = 2.0002;
  74     r[3.72] = 2.0001;
  75
  76     // Let's discover the absissa that will generate a potential of exactly 3.0,
  77     // start by creating 2 ranges for the x and y values:
  78     auto x_range = boost::adaptors::keys(r);
  79     auto y_range = boost::adaptors::values(r);
  80     boost::math::barycentric_rational<double> b(x_range.begin(), x_range.end(), y_range.begin());
  81     //
  82     // We'll use a lamda expression to provide the functor to our root finder, since we want
  83     // the abscissa value that yields 3, not zero.  We pass the functor b by value to the
  84     // lambda expression since barycentric_rational is trivial to copy.
  85     // Here we're using simple bisection to find the root:
  86     boost::uintmax_t iterations = (std::numeric_limits<boost::uintmax_t>::max)();
  87     double abscissa_3 = boost::math::tools::bisect([=](double x) { return b(x) - 3; }, 0.44, 1.24, boost::math::tools::eps_tolerance<double>(), iterations).first;
  88     std::cout << "Abscissa value that yields a potential of 3 = " << abscissa_3 << std::endl;
  89     std::cout << "Root was found in " << iterations << " iterations." << std::endl;
  90     //
  91     // However, we have a more efficient root finding algorithm than simple bisection:
  92     iterations = (std::numeric_limits<boost::uintmax_t>::max)();
  93     abscissa_3 = boost::math::tools::bracket_and_solve_root([=](double x) { return b(x) - 3; }, 0.6, 1.2, false, boost::math::tools::eps_tolerance<double>(), iterations).first;
  94     std::cout << "Abscissa value that yields a potential of 3 = " << abscissa_3 << std::endl;
  95     std::cout << "Root was found in " << iterations << " iterations." << std::endl;
  96 }
  97 //] [/barycentric_rational_example2]
  98
  99
 100 //[barycentric_rational_example2_out
 101 /*` Program output is:
 102 [pre
 103 Abscissa value that yields a potential of 3 = 0.604728
 104 Root was found in 54 iterations.
 105 Abscissa value that yields a potential of 3 = 0.604728
 106 Root was found in 10 iterations.
 107 ]
 108 */
 109 //]