libs/math/doc/html/math_toolkit/double_exponential/de_tanh_sinh_2_arg.html

   1 <html>
   2 <head>
   3 <meta http-equiv="Content-Type" content="text/html; charset=US-ASCII">
   4 <title>Handling functions with large features near an endpoint with tanh-sinh quadrature</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="../double_exponential.html" title="Double-exponential quadrature">
   9 <link rel="prev" href="de_tanh_sinh.html" title="tanh_sinh">
  10 <link rel="next" href="de_sinh_sinh.html" title="sinh_sinh">
  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="de_tanh_sinh.html"><img src="../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../double_exponential.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="de_sinh_sinh.html"><img src="../../../../../../doc/src/images/next.png" alt="Next"></a>
  24 </div>
  25 <div class="section">
  26 <div class="titlepage"><div><div><h3 class="title">
  27 <a name="math_toolkit.double_exponential.de_tanh_sinh_2_arg"></a><a class="link" href="de_tanh_sinh_2_arg.html" title="Handling functions with large features near an endpoint with tanh-sinh quadrature">Handling
  28       functions with large features near an endpoint with tanh-sinh quadrature</a>
  29 </h3></div></div></div>
  30 <p>
  31         Tanh-sinh quadrature has a unique feature which makes it well suited to handling
  32         integrals with either singularities or large features of interest near one
  33         or both endpoints, it turns out that when we calculate and store the abscissa
  34         values at which we will be evaluating the function, we can equally well calculate
  35         the difference between the abscissa value and the nearest endpoint. This
  36         makes it possible to perform quadrature arbitrarily close to an endpoint,
  37         without suffering cancellation error. Note however, that we never actually
  38         reach the endpoint, so any endpoint singularity will always be excluded from
  39         the quadrature.
  40       </p>
  41 <p>
  42         The tanh_sinh integration routine will use this additional information internally
  43         when performing range transformation, so that for example, if one end of
  44         the range is zero (or infinite) then our transformations will get arbitrarily
  45         close to the endpoint without precision loss.
  46       </p>
  47 <p>
  48         However, there are some integrals which may have all of their area near
  49         <span class="emphasis"><em>both</em></span> endpoints, or else near the non-zero endpoint,
  50         and in that situation there is a very real risk of loss of precision. For
  51         example:
  52       </p>
  53 <pre class="programlisting"><span class="identifier">tanh_sinh</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">integrator</span><span class="special">;</span>
  54 <span class="keyword">auto</span> <span class="identifier">f</span> <span class="special">=</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">sqrt</span><span class="special">(</span><span class="identifier">x</span> <span class="special">/</span> <span class="special">(</span><span class="number">1</span> <span class="special">-</span> <span class="identifier">x</span> <span class="special">*</span> <span class="identifier">x</span><span class="special">);</span> <span class="special">};</span>
  55 <span class="keyword">double</span> <span class="identifier">Q</span> <span class="special">=</span> <span class="identifier">integrator</span><span class="special">.</span><span class="identifier">integrate</span><span class="special">(</span><span class="identifier">f</span><span class="special">,</span> <span class="number">0.0</span><span class="special">,</span> <span class="number">1.0</span><span class="special">);</span>
  56 </pre>
  57 <p>
  58         Results in very low accuracy as all the area of the integral is near 1, and
  59         the <code class="computeroutput"><span class="number">1</span> <span class="special">-</span>
  60         <span class="identifier">x</span> <span class="special">*</span> <span class="identifier">x</span></code> term suffers from cancellation error
  61         here.
  62       </p>
  63 <p>
  64         However, both of tanh_sinh's integration routines will automatically handle
  65         2 argument functors: in this case the first argument is the abscissa value
  66         as before, while the second is the distance to the nearest endpoint, ie
  67         <code class="computeroutput"><span class="identifier">a</span> <span class="special">-</span>
  68         <span class="identifier">x</span></code> or <code class="computeroutput"><span class="identifier">b</span>
  69         <span class="special">-</span> <span class="identifier">x</span></code>
  70         if we are integrating over (a,b). You can always differentiate between these
  71         2 cases because the second argument will be negative if we are nearer to
  72         the left endpoint.
  73       </p>
  74 <p>
  75         Knowing this, we can rewrite our lambda expression to take advantage of this
  76         additional information:
  77       </p>
  78 <pre class="programlisting"><span class="identifier">tanh_sinh</span><span class="special">&lt;</span><span class="keyword">double</span><span class="special">&gt;</span> <span class="identifier">integrator</span><span class="special">;</span>
  79 <span class="keyword">auto</span> <span class="identifier">f</span> <span class="special">=</span> <span class="special">[](</span><span class="keyword">double</span> <span class="identifier">x</span><span class="special">,</span> <span class="keyword">double</span> <span class="identifier">xc</span><span class="special">)</span> <span class="special">{</span> <span class="keyword">return</span> <span class="identifier">x</span> <span class="special">&lt;=</span> <span class="number">0.5</span> <span class="special">?</span> <span class="identifier">sqrt</span><span class="special">(</span><span class="identifier">x</span><span class="special">)</span> <span class="special">/</span> <span class="identifier">sqrt</span><span class="special">(</span><span class="number">1</span> <span class="special">-</span> <span class="identifier">x</span> <span class="special">*</span> <span class="identifier">x</span><span class="special">)</span> <span class="special">:</span> <span class="identifier">sqrt</span><span class="special">(</span><span class="identifier">x</span> <span class="special">/</span> <span class="special">((</span><span class="identifier">x</span> <span class="special">+</span> <span class="number">1</span><span class="special">)</span> <span class="special">*</span> <span class="special">(</span><span class="identifier">xc</span><span class="special">)));</span> <span class="special">};</span>
  80 <span class="keyword">double</span> <span class="identifier">Q</span> <span class="special">=</span> <span class="identifier">integrator</span><span class="special">.</span><span class="identifier">integrate</span><span class="special">(</span><span class="identifier">f</span><span class="special">,</span> <span class="number">0.0</span><span class="special">,</span> <span class="number">1.0</span><span class="special">);</span>
  81 </pre>
  82 <p>
  83         Not only is this form accurate to full machine-precision, but it converges
  84         to the result faster as well.
  85       </p>
  86 </div>
  87 <table xmlns:rev="http://www.cs.rpi.edu/~gregod/boost/tools/doc/revision" width="100%"><tr>
  88 <td align="left"></td>
  89 <td align="right"><div class="copyright-footer">Copyright &#169; 2006-2019 Nikhar
  90       Agrawal, Anton Bikineev, Paul A. Bristow, Marco Guazzone, Christopher Kormanyos,
  91       Hubert Holin, Bruno Lalande, John Maddock, Jeremy Murphy, Matthew Pulver, Johan
  92       R&#229;de, Gautam Sewani, Benjamin Sobotta, Nicholas Thompson, Thijs van den Berg,
  93       Daryle Walker and Xiaogang Zhang<p>
  94         Distributed under the Boost Software License, Version 1.0. (See accompanying
  95         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>)
  96       </p>
  97 </div></td>
  98 </tr></table>
  99 <hr>
 100 <div class="spirit-nav">
 101 <a accesskey="p" href="de_tanh_sinh.html"><img src="../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../double_exponential.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="de_sinh_sinh.html"><img src="../../../../../../doc/src/images/next.png" alt="Next"></a>
 102 </div>
 103 </body>
 104 </html>