libs/math/doc/html/math_toolkit/barycentric.html

   1 <html>
   2 <head>
   3 <meta http-equiv="Content-Type" content="text/html; charset=US-ASCII">
   4 <title>Barycentric Rational Interpolation</title>
   5 <link rel="stylesheet" href="../math.css" type="text/css">
   6 <meta name="generator" content="DocBook XSL Stylesheets V1.79.1">
   7 <link rel="home" href="../index.html" title="Math Toolkit 2.11.0">
   8 <link rel="up" href="../interpolation.html" title="Chapter&#160;12.&#160;Interpolation">
   9 <link rel="prev" href="whittaker_shannon.html" title="Whittaker-Shannon interpolation">
  10 <link rel="next" href="vector_barycentric.html" title="Vector-valued Barycentric Rational Interpolation">
  11 </head>
  12 <body bgcolor="white" text="black" link="#0000FF" vlink="#840084" alink="#0000FF">
  13 <table cellpadding="2" width="100%"><tr>
  14 <td valign="top"><img alt="Boost C++ Libraries" width="277" height="86" src="../../../../../boost.png"></td>
  15 <td align="center"><a href="../../../../../index.html">Home</a></td>
  16 <td align="center"><a href="../../../../../libs/libraries.htm">Libraries</a></td>
  17 <td align="center"><a href="http://www.boost.org/users/people.html">People</a></td>
  18 <td align="center"><a href="http://www.boost.org/users/faq.html">FAQ</a></td>
  19 <td align="center"><a href="../../../../../more/index.htm">More</a></td>
  20 </tr></table>
  21 <hr>
  22 <div class="spirit-nav">
  23 <a accesskey="p" href="whittaker_shannon.html"><img src="../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../interpolation.html"><img src="../../../../../doc/src/images/up.png" alt="Up"></a><a accesskey="h" href="../index.html"><img src="../../../../../doc/src/images/home.png" alt="Home"></a><a accesskey="n" href="vector_barycentric.html"><img src="../../../../../doc/src/images/next.png" alt="Next"></a>
  24 </div>
  25 <div class="section">
  26 <div class="titlepage"><div><div><h2 class="title" style="clear: both">
  27 <a name="math_toolkit.barycentric"></a><a class="link" href="barycentric.html" title="Barycentric Rational Interpolation">Barycentric Rational Interpolation</a>
  28 </h2></div></div></div>
  29 <h4>
  30 <a name="math_toolkit.barycentric.h0"></a>
  31       <span class="phrase"><a name="math_toolkit.barycentric.synopsis"></a></span><a class="link" href="barycentric.html#math_toolkit.barycentric.synopsis">Synopsis</a>
  32     </h4>
  33 <pre class="programlisting"><span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">math</span><span class="special">/</span><span class="identifier">interpolators</span><span class="special">/</span><span class="identifier">barycentric_rational</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
  34
  35 <span class="keyword">namespace</span> <span class="identifier">boost</span><span class="special">{</span> <span class="keyword">namespace</span> <span class="identifier">math</span><span class="special">{</span>
  36     <span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Real</span><span class="special">&gt;</span>
  37     <span class="keyword">class</span> <span class="identifier">barycentric_rational</span>
  38     <span class="special">{</span>
  39     <span class="keyword">public</span><span class="special">:</span>
  40         <span class="keyword">template</span> <span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">InputIterator1</span><span class="special">,</span> <span class="keyword">class</span> <span class="identifier">InputIterator2</span><span class="special">&gt;</span>
  41         <span class="identifier">barycentric_rational</span><span class="special">(</span><span class="identifier">InputIterator1</span> <span class="identifier">start_x</span><span class="special">,</span> <span class="identifier">InputIterator1</span> <span class="identifier">end_x</span><span class="special">,</span> <span class="identifier">InputIterator2</span> <span class="identifier">start_y</span><span class="special">,</span> <span class="identifier">size_t</span> <span class="identifier">approximation_order</span> <span class="special">=</span> <span class="number">3</span><span class="special">);</span>
  42
  43         <span class="identifier">barycentric_rational</span><span class="special">(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="identifier">Real</span><span class="special">&gt;&amp;&amp;</span> <span class="identifier">x</span><span class="special">,</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="identifier">Real</span><span class="special">&gt;&amp;&amp;</span> <span class="identifier">y</span><span class="special">,</span> <span class="identifier">size_t</span> <span class="identifier">order</span> <span class="special">=</span> <span class="number">3</span><span class="special">);</span>
  44
  45         <span class="identifier">barycentric_rational</span><span class="special">(</span><span class="keyword">const</span> <span class="identifier">Real</span><span class="special">*</span> <span class="keyword">const</span> <span class="identifier">x</span><span class="special">,</span> <span class="keyword">const</span> <span class="identifier">Real</span><span class="special">*</span> <span class="keyword">const</span> <span class="identifier">y</span><span class="special">,</span> <span class="identifier">size_t</span> <span class="identifier">n</span><span class="special">,</span> <span class="identifier">size_t</span> <span class="identifier">approximation_order</span> <span class="special">=</span> <span class="number">3</span><span class="special">);</span>
  46
  47         <span class="identifier">Real</span> <span class="keyword">operator</span><span class="special">()(</span><span class="identifier">Real</span> <span class="identifier">x</span><span class="special">)</span> <span class="keyword">const</span><span class="special">;</span>
  48
  49         <span class="identifier">Real</span> <span class="identifier">prime</span><span class="special">(</span><span class="identifier">Real</span> <span class="identifier">x</span><span class="special">)</span> <span class="keyword">const</span><span class="special">;</span>
  50
  51         <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="identifier">Real</span><span class="special">&gt;&amp;&amp;</span> <span class="identifier">return_x</span><span class="special">();</span>
  52
  53         <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="identifier">Real</span><span class="special">&gt;&amp;&amp;</span> <span class="identifier">return_y</span><span class="special">();</span>
  54     <span class="special">};</span>
  55
  56 <span class="special">}}</span>
  57 </pre>
  58 <h4>
  59 <a name="math_toolkit.barycentric.h1"></a>
  60       <span class="phrase"><a name="math_toolkit.barycentric.description"></a></span><a class="link" href="barycentric.html#math_toolkit.barycentric.description">Description</a>
  61     </h4>
  62 <p>
  63       Barycentric rational interpolation is a high-accuracy interpolation method
  64       for non-uniformly spaced samples. It requires &#119926;(<span class="emphasis"><em>N</em></span>) time
  65       for construction, and &#119926;(<span class="emphasis"><em>N</em></span>) time for each evaluation. Linear
  66       time evaluation is not optimal; for instance the cubic B-spline can be evaluated
  67       in constant time. However, using the cubic B-spline requires uniformly-spaced
  68       samples, which are not always available.
  69     </p>
  70 <p>
  71       Use of the class requires a vector of independent variables <code class="computeroutput"><span class="identifier">x</span><span class="special">[</span><span class="number">0</span><span class="special">]</span></code>,
  72       <code class="computeroutput"><span class="identifier">x</span><span class="special">[</span><span class="number">1</span><span class="special">]</span></code>, .... <code class="computeroutput"><span class="identifier">x</span><span class="special">[</span><span class="identifier">n</span><span class="special">-</span><span class="number">1</span><span class="special">]</span></code>
  73       where <code class="computeroutput"><span class="identifier">x</span><span class="special">[</span><span class="identifier">i</span><span class="special">+</span><span class="number">1</span><span class="special">]</span> <span class="special">&gt;</span> <span class="identifier">x</span><span class="special">[</span><span class="identifier">i</span><span class="special">]</span></code>,
  74       and a vector of dependent variables <code class="computeroutput"><span class="identifier">y</span><span class="special">[</span><span class="number">0</span><span class="special">]</span></code>,
  75       <code class="computeroutput"><span class="identifier">y</span><span class="special">[</span><span class="number">1</span><span class="special">]</span></code>, ... , <code class="computeroutput"><span class="identifier">y</span><span class="special">[</span><span class="identifier">n</span><span class="special">-</span><span class="number">1</span><span class="special">]</span></code>.
  76       The call is trivial:
  77     </p>
  78 <pre class="programlisting"><span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">x</span><span class="special">(</span><span class="number">500</span><span class="special">);</span>
  79 <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">y</span><span class="special">(</span><span class="number">500</span><span class="special">);</span>
  80 <span class="comment">// populate x, y, then:</span>
  81 <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">barycentric_rational</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">interpolant</span><span class="special">(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">move</span><span class="special">(</span><span class="identifier">x</span><span class="special">),</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">move</span><span class="special">(</span><span class="identifier">y</span><span class="special">));</span>
  82 </pre>
  83 <p>
  84       This implicitly calls the constructor with approximation order 3, and hence
  85       the accuracy is &#119926;(<span class="emphasis"><em>h</em></span><sup>4</sup>). In general, if you require an approximation
  86       order <span class="emphasis"><em>d</em></span>, then the error is &#119926;(<span class="emphasis"><em>h</em></span><sup><span class="emphasis"><em>d</em></span>+1</sup>).
  87       A call to the constructor with an explicit approximation order is demonstrated
  88       below
  89     </p>
  90 <pre class="programlisting"><span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">barycentric_rational</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">interpolant</span><span class="special">(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">move</span><span class="special">(</span><span class="identifier">x</span><span class="special">),</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">move</span><span class="special">(</span><span class="identifier">y</span><span class="special">),</span> <span class="number">5</span><span class="special">);</span>
  91 </pre>
  92 <p>
  93       To evaluate the interpolant, simply use
  94     </p>
  95 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">x</span> <span class="special">=</span> <span class="number">2.3</span><span class="special">;</span>
  96 <span class="keyword">double</span> <span class="identifier">y</span> <span class="special">=</span> <span class="identifier">interpolant</span><span class="special">(</span><span class="identifier">x</span><span class="special">);</span>
  97 </pre>
  98 <p>
  99       and to evaluate its derivative use
 100     </p>
 101 <pre class="programlisting"><span class="keyword">double</span> <span class="identifier">y</span> <span class="special">=</span> <span class="identifier">interpolant</span><span class="special">.</span><span class="identifier">prime</span><span class="special">(</span><span class="identifier">x</span><span class="special">);</span>
 102 </pre>
 103 <p>
 104       If you no longer require the interpolant, then you can get your data back:
 105     </p>
 106 <pre class="programlisting"><span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">xs</span> <span class="special">=</span> <span class="identifier">interpolant</span><span class="special">.</span><span class="identifier">return_x</span><span class="special">();</span>
 107 <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">ys</span> <span class="special">=</span> <span class="identifier">interpolant</span><span class="special">.</span><span class="identifier">return_y</span><span class="special">();</span>
 108 </pre>
 109 <p>
 110       Be aware that once you return your data, the interpolant is <span class="bold"><strong>dead</strong></span>.
 111     </p>
 112 <h4>
 113 <a name="math_toolkit.barycentric.h2"></a>
 114       <span class="phrase"><a name="math_toolkit.barycentric.caveats"></a></span><a class="link" href="barycentric.html#math_toolkit.barycentric.caveats">Caveats</a>
 115     </h4>
 116 <p>
 117       Although this algorithm is robust, it can surprise you. The main way this occurs
 118       is if the sample spacing at the endpoints is much larger than the spacing in
 119       the center. This is to be expected; all interpolants perform better in the
 120       opposite regime, where samples are clustered at the endpoints and somewhat
 121       uniformly spaced throughout the center.
 122     </p>
 123 <p>
 124       A desirable property of any interpolator <span class="emphasis"><em>f</em></span> is that for
 125       all <span class="emphasis"><em>x</em></span><sub>min</sub> &#8804; <span class="emphasis"><em>x</em></span> &#8804; <span class="emphasis"><em>x</em></span><sub>max</sub>,
 126       also <span class="emphasis"><em>y</em></span><sub>min</sub> &#8804; <span class="emphasis"><em>f</em></span>(<span class="emphasis"><em>x</em></span>)
 127       &#8804; <span class="emphasis"><em>y</em></span><sub>max</sub>.
 128     </p>
 129 <p>
 130       <span class="emphasis"><em>This property does not hold for the barycentric rational interpolator.</em></span>
 131       However, unless you deliberately try to antagonize this interpolator (by, for
 132       instance, placing the final value far from all the rest), you will probably
 133       not fall victim to this problem.
 134     </p>
 135 <p>
 136       The reference used for implementation of this algorithm is <a href="https://web.archive.org/save/_embed/http://www.mn.uio.no/math/english/people/aca/michaelf/papers/rational.pdf" target="_top">Barycentric
 137       rational interpolation with no poles and a high rate of interpolation</a>,
 138       and the evaluation of the derivative is given by <a href="http://www.ams.org/journals/mcom/1986-47-175/S0025-5718-1986-0842136-8/S0025-5718-1986-0842136-8.pdf" target="_top">Some
 139       New Aspects of Rational Interpolation</a>.
 140     </p>
 141 <h4>
 142 <a name="math_toolkit.barycentric.h3"></a>
 143       <span class="phrase"><a name="math_toolkit.barycentric.examples"></a></span><a class="link" href="barycentric.html#math_toolkit.barycentric.examples">Examples</a>
 144     </h4>
 145 <p>
 146       This example shows how to use barycentric rational interpolation, using Walter
 147       Kohn's classic paper "Solution of the Schrodinger Equation in Periodic
 148       Lattices with an Application to Metallic Lithium" In this paper, Kohn
 149       needs to repeatedly solve an ODE (the radial Schrodinger equation) given a
 150       potential which is only known at non-equally samples data.
 151     </p>
 152 <p>
 153       If he'd only had the barycentric rational interpolant of Boost.Math!
 154     </p>
 155 <p>
 156       References: Kohn, W., and N. Rostoker. "Solution of the Schrodinger equation
 157       in periodic lattices with an application to metallic lithium." Physical
 158       Review 94.5 (1954): 1111.
 159     </p>
 160 <pre class="programlisting"><span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">math</span><span class="special">/</span><span class="identifier">interpolators</span><span class="special">/</span><span class="identifier">barycentric_rational</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
 161
 162 <span class="keyword">int</span> <span class="identifier">main</span><span class="special">()</span>
 163 <span class="special">{</span>
 164     <span class="comment">// The lithium potential is given in Kohn's paper, Table I:</span>
 165     <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">r</span><span class="special">(</span><span class="number">45</span><span class="special">);</span>
 166     <span class="identifier">std</span><span class="special">::</span><span class="identifier">vector</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">mrV</span><span class="special">(</span><span class="number">45</span><span class="special">);</span>
 167
 168     <span class="comment">// We'll skip the code for filling the above vectors with data for now...</span>
 169
 170     <span class="comment">// Now we want to interpolate this potential at any r:</span>
 171     <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">barycentric_rational</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">b</span><span class="special">(</span><span class="identifier">r</span><span class="special">.</span><span class="identifier">data</span><span class="special">(),</span> <span class="identifier">mrV</span><span class="special">.</span><span class="identifier">data</span><span class="special">(),</span> <span class="identifier">r</span><span class="special">.</span><span class="identifier">size</span><span class="special">());</span>
 172
 173     <span class="keyword">for</span> <span class="special">(</span><span class="identifier">size_t</span> <span class="identifier">i</span> <span class="special">=</span> <span class="number">1</span><span class="special">;</span> <span class="identifier">i</span> <span class="special">&lt;</span> <span class="number">8</span><span class="special">;</span> <span class="special">++</span><span class="identifier">i</span><span class="special">)</span>
 174     <span class="special">{</span>
 175         <span class="keyword">double</span> <span class="identifier">r</span> <span class="special">=</span> <span class="identifier">i</span><span class="special">*</span><span class="number">0.5</span><span class="special">;</span>
 176         <span class="identifier">std</span><span class="special">::</span><span class="identifier">cout</span> <span class="special">&lt;&lt;</span>  <span class="string">"(r, V) = ("</span> <span class="special">&lt;&lt;</span> <span class="identifier">r</span> <span class="special">&lt;&lt;</span> <span class="string">", "</span> <span class="special">&lt;&lt;</span> <span class="special">-</span><span class="identifier">b</span><span class="special">(</span><span class="identifier">r</span><span class="special">)/</span><span class="identifier">r</span> <span class="special">&lt;&lt;</span> <span class="string">")\n"</span><span class="special">;</span>
 177     <span class="special">}</span>
 178 <span class="special">}</span>
 179 </pre>
 180 <p>
 181       This further example shows how to use the iterator based constructor, and then
 182       uses the function object in our root finding algorithms to locate the points
 183       where the potential achieves a specific value.
 184     </p>
 185 <pre class="programlisting"><span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">math</span><span class="special">/</span><span class="identifier">interpolators</span><span class="special">/</span><span class="identifier">barycentric_rational</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
 186 <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">range</span><span class="special">/</span><span class="identifier">adaptors</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
 187 <span class="preprocessor">#include</span> <span class="special">&lt;</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">math</span><span class="special">/</span><span class="identifier">tools</span><span class="special">/</span><span class="identifier">roots</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">&gt;</span>
 188
 189 <span class="keyword">int</span> <span class="identifier">main</span><span class="special">()</span>
 190 <span class="special">{</span>
 191     <span class="comment">// The lithium potential is given in Kohn's paper, Table I.</span>
 192     <span class="comment">// (We could equally easily use an unordered_map, a list of tuples or pairs, or a 2-dimentional array).</span>
 193     <span class="identifier">std</span><span class="special">::</span><span class="identifier">map</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">,</span> <span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">r</span><span class="special">;</span>
 194
 195     <span class="identifier">r</span><span class="special">[</span><span class="number">0.02</span><span class="special">]</span> <span class="special">=</span> <span class="number">5.727</span><span class="special">;</span>
 196     <span class="identifier">r</span><span class="special">[</span><span class="number">0.04</span><span class="special">]</span> <span class="special">=</span> <span class="number">5.544</span><span class="special">;</span>
 197     <span class="identifier">r</span><span class="special">[</span><span class="number">0.06</span><span class="special">]</span> <span class="special">=</span> <span class="number">5.450</span><span class="special">;</span>
 198     <span class="identifier">r</span><span class="special">[</span><span class="number">0.08</span><span class="special">]</span> <span class="special">=</span> <span class="number">5.351</span><span class="special">;</span>
 199     <span class="identifier">r</span><span class="special">[</span><span class="number">0.10</span><span class="special">]</span> <span class="special">=</span> <span class="number">5.253</span><span class="special">;</span>
 200     <span class="identifier">r</span><span class="special">[</span><span class="number">0.12</span><span class="special">]</span> <span class="special">=</span> <span class="number">5.157</span><span class="special">;</span>
 201     <span class="identifier">r</span><span class="special">[</span><span class="number">0.14</span><span class="special">]</span> <span class="special">=</span> <span class="number">5.058</span><span class="special">;</span>
 202     <span class="identifier">r</span><span class="special">[</span><span class="number">0.16</span><span class="special">]</span> <span class="special">=</span> <span class="number">4.960</span><span class="special">;</span>
 203     <span class="identifier">r</span><span class="special">[</span><span class="number">0.18</span><span class="special">]</span> <span class="special">=</span> <span class="number">4.862</span><span class="special">;</span>
 204     <span class="identifier">r</span><span class="special">[</span><span class="number">0.20</span><span class="special">]</span> <span class="special">=</span> <span class="number">4.762</span><span class="special">;</span>
 205     <span class="identifier">r</span><span class="special">[</span><span class="number">0.24</span><span class="special">]</span> <span class="special">=</span> <span class="number">4.563</span><span class="special">;</span>
 206     <span class="identifier">r</span><span class="special">[</span><span class="number">0.28</span><span class="special">]</span> <span class="special">=</span> <span class="number">4.360</span><span class="special">;</span>
 207     <span class="identifier">r</span><span class="special">[</span><span class="number">0.32</span><span class="special">]</span> <span class="special">=</span> <span class="number">4.1584</span><span class="special">;</span>
 208     <span class="identifier">r</span><span class="special">[</span><span class="number">0.36</span><span class="special">]</span> <span class="special">=</span> <span class="number">3.9463</span><span class="special">;</span>
 209     <span class="identifier">r</span><span class="special">[</span><span class="number">0.40</span><span class="special">]</span> <span class="special">=</span> <span class="number">3.7360</span><span class="special">;</span>
 210     <span class="identifier">r</span><span class="special">[</span><span class="number">0.44</span><span class="special">]</span> <span class="special">=</span> <span class="number">3.5429</span><span class="special">;</span>
 211     <span class="identifier">r</span><span class="special">[</span><span class="number">0.48</span><span class="special">]</span> <span class="special">=</span> <span class="number">3.3797</span><span class="special">;</span>
 212     <span class="identifier">r</span><span class="special">[</span><span class="number">0.52</span><span class="special">]</span> <span class="special">=</span> <span class="number">3.2417</span><span class="special">;</span>
 213     <span class="identifier">r</span><span class="special">[</span><span class="number">0.56</span><span class="special">]</span> <span class="special">=</span> <span class="number">3.1209</span><span class="special">;</span>
 214     <span class="identifier">r</span><span class="special">[</span><span class="number">0.60</span><span class="special">]</span> <span class="special">=</span> <span class="number">3.0138</span><span class="special">;</span>
 215     <span class="identifier">r</span><span class="special">[</span><span class="number">0.68</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.8342</span><span class="special">;</span>
 216     <span class="identifier">r</span><span class="special">[</span><span class="number">0.76</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.6881</span><span class="special">;</span>
 217     <span class="identifier">r</span><span class="special">[</span><span class="number">0.84</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.5662</span><span class="special">;</span>
 218     <span class="identifier">r</span><span class="special">[</span><span class="number">0.92</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.4242</span><span class="special">;</span>
 219     <span class="identifier">r</span><span class="special">[</span><span class="number">1.00</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.3766</span><span class="special">;</span>
 220     <span class="identifier">r</span><span class="special">[</span><span class="number">1.08</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.3058</span><span class="special">;</span>
 221     <span class="identifier">r</span><span class="special">[</span><span class="number">1.16</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.2458</span><span class="special">;</span>
 222     <span class="identifier">r</span><span class="special">[</span><span class="number">1.24</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.2035</span><span class="special">;</span>
 223     <span class="identifier">r</span><span class="special">[</span><span class="number">1.32</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.1661</span><span class="special">;</span>
 224     <span class="identifier">r</span><span class="special">[</span><span class="number">1.40</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.1350</span><span class="special">;</span>
 225     <span class="identifier">r</span><span class="special">[</span><span class="number">1.48</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.1090</span><span class="special">;</span>
 226     <span class="identifier">r</span><span class="special">[</span><span class="number">1.64</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0697</span><span class="special">;</span>
 227     <span class="identifier">r</span><span class="special">[</span><span class="number">1.80</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0466</span><span class="special">;</span>
 228     <span class="identifier">r</span><span class="special">[</span><span class="number">1.96</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0325</span><span class="special">;</span>
 229     <span class="identifier">r</span><span class="special">[</span><span class="number">2.12</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0288</span><span class="special">;</span>
 230     <span class="identifier">r</span><span class="special">[</span><span class="number">2.28</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0292</span><span class="special">;</span>
 231     <span class="identifier">r</span><span class="special">[</span><span class="number">2.44</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0228</span><span class="special">;</span>
 232     <span class="identifier">r</span><span class="special">[</span><span class="number">2.60</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0124</span><span class="special">;</span>
 233     <span class="identifier">r</span><span class="special">[</span><span class="number">2.76</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0065</span><span class="special">;</span>
 234     <span class="identifier">r</span><span class="special">[</span><span class="number">2.92</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0031</span><span class="special">;</span>
 235     <span class="identifier">r</span><span class="special">[</span><span class="number">3.08</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0015</span><span class="special">;</span>
 236     <span class="identifier">r</span><span class="special">[</span><span class="number">3.24</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0008</span><span class="special">;</span>
 237     <span class="identifier">r</span><span class="special">[</span><span class="number">3.40</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0004</span><span class="special">;</span>
 238     <span class="identifier">r</span><span class="special">[</span><span class="number">3.56</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0002</span><span class="special">;</span>
 239     <span class="identifier">r</span><span class="special">[</span><span class="number">3.72</span><span class="special">]</span> <span class="special">=</span> <span class="number">2.0001</span><span class="special">;</span>
 240
 241     <span class="comment">// Let's discover the absissa that will generate a potential of exactly 3.0,</span>
 242     <span class="comment">// start by creating 2 ranges for the x and y values:</span>
 243     <span class="keyword">auto</span> <span class="identifier">x_range</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">adaptors</span><span class="special">::</span><span class="identifier">keys</span><span class="special">(</span><span class="identifier">r</span><span class="special">);</span>
 244     <span class="keyword">auto</span> <span class="identifier">y_range</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">adaptors</span><span class="special">::</span><span class="identifier">values</span><span class="special">(</span><span class="identifier">r</span><span class="special">);</span>
 245     <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">barycentric_rational</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">b</span><span class="special">(</span><span class="identifier">x_range</span><span class="special">.</span><span class="identifier">begin</span><span class="special">(),</span> <span class="identifier">x_range</span><span class="special">.</span><span class="identifier">end</span><span class="special">(),</span> <span class="identifier">y_range</span><span class="special">.</span><span class="identifier">begin</span><span class="special">());</span>
 246     <span class="comment">//</span>
 247     <span class="comment">// We'll use a lamda expression to provide the functor to our root finder, since we want</span>
 248     <span class="comment">// the abscissa value that yields 3, not zero.  We pass the functor b by value to the</span>
 249     <span class="comment">// lambda expression since barycentric_rational is trivial to copy.</span>
 250     <span class="comment">// Here we're using simple bisection to find the root:</span>
 251     <span class="identifier">boost</span><span class="special">::</span><span class="identifier">uintmax_t</span> <span class="identifier">iterations</span> <span class="special">=</span> <span class="special">(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">numeric_limits</span><span class="special">&lt;</span><span class="identifier">boost</span><span class="special">::</span><span class="identifier">uintmax_t</span><span class="special">&gt;::</span><span class="identifier">max</span><span class="special">)();</span>
 252     <span class="keyword">double</span> <span class="identifier">abscissa_3</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">tools</span><span class="special">::</span><span class="identifier">bisect</span><span class="special">([=](</span><span class="keyword">double</span> <span class="identifier">x</span><span class="special">)</span> <span class="special">{</span> <span class="keyword">return</span> <span class="identifier">b</span><span class="special">(</span><span class="identifier">x</span><span class="special">)</span> <span class="special">-</span> <span class="number">3</span><span class="special">;</span> <span class="special">},</span> <span class="number">0.44</span><span class="special">,</span> <span class="number">1.24</span><span class="special">,</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">tools</span><span class="special">::</span><span class="identifier">eps_tolerance</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;(),</span> <span class="identifier">iterations</span><span class="special">).</span><span class="identifier">first</span><span class="special">;</span>
 253     <span class="identifier">std</span><span class="special">::</span><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Abscissa value that yields a potential of 3 = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">abscissa_3</span> <span class="special">&lt;&lt;</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">endl</span><span class="special">;</span>
 254     <span class="identifier">std</span><span class="special">::</span><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Root was found in "</span> <span class="special">&lt;&lt;</span> <span class="identifier">iterations</span> <span class="special">&lt;&lt;</span> <span class="string">" iterations."</span> <span class="special">&lt;&lt;</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">endl</span><span class="special">;</span>
 255     <span class="comment">//</span>
 256     <span class="comment">// However, we have a more efficient root finding algorithm than simple bisection:</span>
 257     <span class="identifier">iterations</span> <span class="special">=</span> <span class="special">(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">numeric_limits</span><span class="special">&lt;</span><span class="identifier">boost</span><span class="special">::</span><span class="identifier">uintmax_t</span><span class="special">&gt;::</span><span class="identifier">max</span><span class="special">)();</span>
 258     <span class="identifier">abscissa_3</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">tools</span><span class="special">::</span><span class="identifier">bracket_and_solve_root</span><span class="special">([=](</span><span class="keyword">double</span> <span class="identifier">x</span><span class="special">)</span> <span class="special">{</span> <span class="keyword">return</span> <span class="identifier">b</span><span class="special">(</span><span class="identifier">x</span><span class="special">)</span> <span class="special">-</span> <span class="number">3</span><span class="special">;</span> <span class="special">},</span> <span class="number">0.6</span><span class="special">,</span> <span class="number">1.2</span><span class="special">,</span> <span class="keyword">false</span><span class="special">,</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">tools</span><span class="special">::</span><span class="identifier">eps_tolerance</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;(),</span> <span class="identifier">iterations</span><span class="special">).</span><span class="identifier">first</span><span class="special">;</span>
 259     <span class="identifier">std</span><span class="special">::</span><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Abscissa value that yields a potential of 3 = "</span> <span class="special">&lt;&lt;</span> <span class="identifier">abscissa_3</span> <span class="special">&lt;&lt;</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">endl</span><span class="special">;</span>
 260     <span class="identifier">std</span><span class="special">::</span><span class="identifier">cout</span> <span class="special">&lt;&lt;</span> <span class="string">"Root was found in "</span> <span class="special">&lt;&lt;</span> <span class="identifier">iterations</span> <span class="special">&lt;&lt;</span> <span class="string">" iterations."</span> <span class="special">&lt;&lt;</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">endl</span><span class="special">;</span>
 261 <span class="special">}</span>
 262 </pre>
 263 <pre class="programlisting"><span class="identifier">Program</span> <span class="identifier">output</span> <span class="identifier">is</span><span class="special">:</span>
 264 </pre>
 265 <pre class="programlisting">Abscissa value that yields a potential of 3 = 0.604728
 266 Root was found in 54 iterations.
 267 Abscissa value that yields a potential of 3 = 0.604728
 268 Root was found in 10 iterations.
 269 </pre>
 270 </div>
 271 <table xmlns:rev="http://www.cs.rpi.edu/~gregod/boost/tools/doc/revision" width="100%"><tr>
 272 <td align="left"></td>
 273 <td align="right"><div class="copyright-footer">Copyright &#169; 2006-2019 Nikhar
 274       Agrawal, Anton Bikineev, Paul A. Bristow, Marco Guazzone, Christopher Kormanyos,
 275       Hubert Holin, Bruno Lalande, John Maddock, Jeremy Murphy, Matthew Pulver, Johan
 276       R&#229;de, Gautam Sewani, Benjamin Sobotta, Nicholas Thompson, Thijs van den Berg,
 277       Daryle Walker and Xiaogang Zhang<p>
 278         Distributed under the Boost Software License, Version 1.0. (See accompanying
 279         file LICENSE_1_0.txt or copy at <a href="http://www.boost.org/LICENSE_1_0.txt" target="_top">http://www.boost.org/LICENSE_1_0.txt</a>)
 280       </p>
 281 </div></td>
 282 </tr></table>
 283 <hr>
 284 <div class="spirit-nav">
 285 <a accesskey="p" href="whittaker_shannon.html"><img src="../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../interpolation.html"><img src="../../../../../doc/src/images/up.png" alt="Up"></a><a accesskey="h" href="../index.html"><img src="../../../../../doc/src/images/home.png" alt="Home"></a><a accesskey="n" href="vector_barycentric.html"><img src="../../../../../doc/src/images/next.png" alt="Next"></a>
 286 </div>
 287 </body>
 288 </html>