2
0
mirror of https://github.com/boostorg/math.git synced 2026-01-22 05:22:15 +00:00
Files
math/doc/html/math_toolkit/double_exponential/de_exp_sinh.html
2024-08-06 13:18:32 +01:00

118 lines
16 KiB
HTML

<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
<title>exp_sinh</title>
<link rel="stylesheet" href="../../math.css" type="text/css">
<meta name="generator" content="DocBook XSL Stylesheets V1.79.1">
<link rel="home" href="../../index.html" title="Math Toolkit 4.2.1">
<link rel="up" href="../double_exponential.html" title="Double-exponential quadrature">
<link rel="prev" href="de_sinh_sinh.html" title="sinh_sinh">
<link rel="next" href="de_tol.html" title="Setting the Termination Condition for Integration">
<meta name="viewport" content="width=device-width, initial-scale=1">
</head>
<body bgcolor="white" text="black" link="#0000FF" vlink="#840084" alink="#0000FF">
<table cellpadding="2" width="100%"><tr>
<td valign="top"><img alt="Boost C++ Libraries" width="277" height="86" src="../../../../../../boost.png"></td>
<td align="center"><a href="../../../../../../index.html">Home</a></td>
<td align="center"><a href="../../../../../../libs/libraries.htm">Libraries</a></td>
<td align="center"><a href="http://www.boost.org/users/people.html">People</a></td>
<td align="center"><a href="http://www.boost.org/users/faq.html">FAQ</a></td>
<td align="center"><a href="../../../../../../more/index.htm">More</a></td>
</tr></table>
<hr>
<div class="spirit-nav">
<a accesskey="p" href="de_sinh_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_tol.html"><img src="../../../../../../doc/src/images/next.png" alt="Next"></a>
</div>
<div class="section">
<div class="titlepage"><div><div><h3 class="title">
<a name="math_toolkit.double_exponential.de_exp_sinh"></a><a class="link" href="de_exp_sinh.html" title="exp_sinh">exp_sinh</a>
</h3></div></div></div>
<pre class="programlisting"><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>
<span class="keyword">class</span> <span class="identifier">exp_sinh</span>
<span class="special">{</span>
<span class="keyword">public</span><span class="special">:</span>
<span class="identifier">exp_sinh</span><span class="special">(</span><span class="identifier">size_t</span> <span class="identifier">max_refinements</span> <span class="special">=</span> <span class="number">9</span><span class="special">);</span>
<span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">F</span><span class="special">&gt;</span>
<span class="keyword">auto</span> <span class="identifier">integrate</span><span class="special">(</span><span class="keyword">const</span> <span class="identifier">F</span> <span class="identifier">f</span><span class="special">,</span> <span class="identifier">Real</span> <span class="identifier">a</span><span class="special">,</span> <span class="identifier">Real</span> <span class="identifier">b</span><span class="special">,</span>
<span class="identifier">Real</span> <span class="identifier">tol</span> <span class="special">=</span> <span class="identifier">sqrt</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">Real</span><span class="special">&gt;::</span><span class="identifier">epsilon</span><span class="special">()),</span>
<span class="identifier">Real</span><span class="special">*</span> <span class="identifier">error</span> <span class="special">=</span> <span class="keyword">nullptr</span><span class="special">,</span>
<span class="identifier">Real</span><span class="special">*</span> <span class="identifier">L1</span> <span class="special">=</span> <span class="keyword">nullptr</span><span class="special">,</span>
<span class="identifier">size_t</span><span class="special">*</span> <span class="identifier">levels</span> <span class="special">=</span> <span class="keyword">nullptr</span><span class="special">)-&gt;</span><span class="keyword">decltype</span><span class="special">(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">declval</span><span class="special">&lt;</span><span class="identifier">F</span><span class="special">&gt;()(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">declval</span><span class="special">&lt;</span><span class="identifier">Real</span><span class="special">&gt;()))</span> <span class="keyword">const</span><span class="special">;</span>
<span class="keyword">template</span><span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">F</span><span class="special">&gt;</span>
<span class="keyword">auto</span> <span class="identifier">integrate</span><span class="special">(</span><span class="keyword">const</span> <span class="identifier">F</span> <span class="identifier">f</span><span class="special">,</span>
<span class="identifier">Real</span> <span class="identifier">tol</span> <span class="special">=</span> <span class="identifier">sqrt</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">Real</span><span class="special">&gt;::</span><span class="identifier">epsilon</span><span class="special">()),</span>
<span class="identifier">Real</span><span class="special">*</span> <span class="identifier">error</span> <span class="special">=</span> <span class="keyword">nullptr</span><span class="special">,</span>
<span class="identifier">Real</span><span class="special">*</span> <span class="identifier">L1</span> <span class="special">=</span> <span class="keyword">nullptr</span><span class="special">,</span>
<span class="identifier">size_t</span><span class="special">*</span> <span class="identifier">levels</span> <span class="special">=</span> <span class="keyword">nullptr</span><span class="special">)-&gt;</span><span class="keyword">decltype</span><span class="special">(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">declval</span><span class="special">&lt;</span><span class="identifier">F</span><span class="special">&gt;()(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">declval</span><span class="special">&lt;</span><span class="identifier">Real</span><span class="special">&gt;()))</span> <span class="keyword">const</span><span class="special">;</span>
<span class="special">};</span>
</pre>
<p>
For half-infinite intervals, the <code class="computeroutput"><span class="identifier">exp</span><span class="special">-</span><span class="identifier">sinh</span></code> quadrature
is provided:
</p>
<pre class="programlisting"><span class="identifier">exp_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>
<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">exp</span><span class="special">(-</span><span class="number">3</span><span class="special">*</span><span class="identifier">x</span><span class="special">);</span> <span class="special">};</span>
<span class="keyword">double</span> <span class="identifier">termination</span> <span class="special">=</span> <span class="identifier">sqrt</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="keyword">double</span><span class="special">&gt;::</span><span class="identifier">epsilon</span><span class="special">());</span>
<span class="keyword">double</span> <span class="identifier">error</span><span class="special">;</span>
<span class="keyword">double</span> <span class="identifier">L1</span><span class="special">;</span>
<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="identifier">termination</span><span class="special">,</span> <span class="special">&amp;</span><span class="identifier">error</span><span class="special">,</span> <span class="special">&amp;</span><span class="identifier">L1</span><span class="special">);</span>
</pre>
<p>
The native integration range of this integrator is (0, ∞), but we also support
/(a, ∞), (-∞, 0)/ and /(-∞, b)/ via argument transformations.
</p>
<p>
Endpoint singularities and complex-valued integrands are supported by <code class="computeroutput"><span class="identifier">exp</span><span class="special">-</span><span class="identifier">sinh</span></code>.
</p>
<p>
For example, the modified Bessel function K can be represented via:
</p>
<div class="blockquote"><blockquote class="blockquote"><p>
<span class="inlinemediaobject"><img src="../../../equations/complex_bessel_k_integral.svg"></span>
</p></blockquote></div>
<p>
Which we can code up as:
</p>
<pre class="programlisting"><span class="keyword">template</span> <span class="special">&lt;</span><span class="keyword">class</span> <span class="identifier">Complex</span><span class="special">&gt;</span>
<span class="identifier">Complex</span> <span class="identifier">bessel_K</span><span class="special">(</span><span class="identifier">Complex</span> <span class="identifier">alpha</span><span class="special">,</span> <span class="identifier">Complex</span> <span class="identifier">z</span><span class="special">)</span>
<span class="special">{</span>
<span class="keyword">typedef</span> <span class="keyword">typename</span> <span class="identifier">Complex</span><span class="special">::</span><span class="identifier">value_type</span> <span class="identifier">value_type</span><span class="special">;</span>
<span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">cosh</span><span class="special">;</span> <span class="keyword">using</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">exp</span><span class="special">;</span>
<span class="keyword">auto</span> <span class="identifier">f</span> <span class="special">=</span> <span class="special">[&amp;,</span> <span class="identifier">alpha</span><span class="special">,</span> <span class="identifier">z</span><span class="special">](</span><span class="identifier">value_type</span> <span class="identifier">t</span><span class="special">)</span>
<span class="special">{</span>
<span class="identifier">value_type</span> <span class="identifier">ct</span> <span class="special">=</span> <span class="identifier">cosh</span><span class="special">(</span><span class="identifier">t</span><span class="special">);</span>
<span class="keyword">if</span> <span class="special">(</span><span class="identifier">ct</span> <span class="special">&gt;</span> <span class="identifier">log</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">value_type</span><span class="special">&gt;::</span><span class="identifier">max</span><span class="special">()))</span>
<span class="keyword">return</span> <span class="identifier">Complex</span><span class="special">(</span><span class="number">0</span><span class="special">);</span>
<span class="keyword">return</span> <span class="identifier">exp</span><span class="special">(-</span><span class="identifier">z</span> <span class="special">*</span> <span class="identifier">ct</span><span class="special">)</span> <span class="special">*</span> <span class="identifier">cosh</span><span class="special">(</span><span class="identifier">alpha</span> <span class="special">*</span> <span class="identifier">t</span><span class="special">);</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">quadrature</span><span class="special">::</span><span class="identifier">exp_sinh</span><span class="special">&lt;</span><span class="identifier">value_type</span><span class="special">&gt;</span> <span class="identifier">integrator</span><span class="special">;</span>
<span class="keyword">return</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="special">}</span>
</pre>
<p>
The only wrinkle in the above code is the need to check for large <code class="computeroutput"><span class="identifier">cosh</span><span class="special">(</span><span class="identifier">t</span><span class="special">)</span></code> in which case we assume that <code class="computeroutput"><span class="identifier">exp</span><span class="special">(-</span><span class="identifier">x</span>
<span class="identifier">cosh</span><span class="special">(</span><span class="identifier">t</span><span class="special">))</span></code> tends
to zero faster than <code class="computeroutput"><span class="identifier">cosh</span><span class="special">(</span><span class="identifier">alpha</span> <span class="identifier">x</span><span class="special">)</span></code> tends
to infinity and return <code class="computeroutput"><span class="number">0</span></code>. Without
that check we end up with <code class="computeroutput"><span class="number">0</span> <span class="special">*</span> <span class="identifier">Infinity</span></code>
as the result (a NaN).
</p>
</div>
<div class="copyright-footer">Copyright © 2006-2021 Nikhar Agrawal, Anton Bikineev, Matthew Borland,
Paul A. Bristow, Marco Guazzone, Christopher Kormanyos, Hubert Holin, Bruno
Lalande, John Maddock, Evan Miller, Jeremy Murphy, Matthew Pulver, Johan Råde,
Gautam Sewani, Benjamin Sobotta, Nicholas Thompson, Thijs van den Berg, Daryle
Walker and Xiaogang Zhang<p>
Distributed under the Boost Software License, Version 1.0. (See accompanying
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>)
</p>
</div>
<hr>
<div class="spirit-nav">
<a accesskey="p" href="de_sinh_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_tol.html"><img src="../../../../../../doc/src/images/next.png" alt="Next"></a>
</div>
</body>
</html>