mirror of
https://github.com/TheAlgorithms/C-Plus-Plus.git
synced 2026-03-21 20:31:43 +08:00
485 lines
42 KiB
HTML
485 lines
42 KiB
HTML
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" "https://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
|
||
<html xmlns="http://www.w3.org/1999/xhtml">
|
||
<head>
|
||
<meta http-equiv="Content-Type" content="text/xhtml;charset=UTF-8"/>
|
||
<meta http-equiv="X-UA-Compatible" content="IE=11"/>
|
||
<meta name="generator" content="Doxygen 1.9.3"/>
|
||
<meta name="viewport" content="width=device-width, initial-scale=1"/>
|
||
<title>Algorithms_in_C++: math/integral_approximation2.cpp File Reference</title>
|
||
<link href="../../tabs.css" rel="stylesheet" type="text/css"/>
|
||
<script type="text/javascript" src="../../jquery.js"></script>
|
||
<script type="text/javascript" src="../../dynsections.js"></script>
|
||
<link href="../../navtree.css" rel="stylesheet" type="text/css"/>
|
||
<script type="text/javascript" src="../../resize.js"></script>
|
||
<script type="text/javascript" src="../../navtreedata.js"></script>
|
||
<script type="text/javascript" src="../../navtree.js"></script>
|
||
<link href="../../search/search.css" rel="stylesheet" type="text/css"/>
|
||
<script type="text/javascript" src="../../search/searchdata.js"></script>
|
||
<script type="text/javascript" src="../../search/search.js"></script>
|
||
<script type="text/x-mathjax-config">
|
||
MathJax.Hub.Config({
|
||
extensions: ["tex2jax.js", "TeX/AMSmath.js", "TeX/AMSsymbols.js"],
|
||
jax: ["input/TeX","output/HTML-CSS"],
|
||
});
|
||
</script>
|
||
<script type="text/javascript" async="async" src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.7/MathJax.js?config=TeX-MML-AM_CHTML/MathJax.js"></script>
|
||
<link href="../../doxygen.css" rel="stylesheet" type="text/css" />
|
||
</head>
|
||
<body>
|
||
<div id="top"><!-- do not remove this div, it is closed by doxygen! -->
|
||
<div id="titlearea">
|
||
<table cellspacing="0" cellpadding="0">
|
||
<tbody>
|
||
<tr id="projectrow">
|
||
<td id="projectalign">
|
||
<div id="projectname">Algorithms_in_C++<span id="projectnumber"> 1.0.0</span>
|
||
</div>
|
||
<div id="projectbrief">Set of algorithms implemented in C++.</div>
|
||
</td>
|
||
</tr>
|
||
</tbody>
|
||
</table>
|
||
</div>
|
||
<!-- end header part -->
|
||
<!-- Generated by Doxygen 1.9.3 -->
|
||
<script type="text/javascript">
|
||
/* @license magnet:?xt=urn:btih:d3d9a9a6595521f9666a5e94cc830dab83b65699&dn=expat.txt MIT */
|
||
var searchBox = new SearchBox("searchBox", "../../search",'Search','.html');
|
||
/* @license-end */
|
||
</script>
|
||
<script type="text/javascript" src="../../menudata.js"></script>
|
||
<script type="text/javascript" src="../../menu.js"></script>
|
||
<script type="text/javascript">
|
||
/* @license magnet:?xt=urn:btih:d3d9a9a6595521f9666a5e94cc830dab83b65699&dn=expat.txt MIT */
|
||
$(function() {
|
||
initMenu('../../',true,false,'search.php','Search');
|
||
$(document).ready(function() { init_search(); });
|
||
});
|
||
/* @license-end */
|
||
</script>
|
||
<div id="main-nav"></div>
|
||
</div><!-- top -->
|
||
<div id="side-nav" class="ui-resizable side-nav-resizable">
|
||
<div id="nav-tree">
|
||
<div id="nav-tree-contents">
|
||
<div id="nav-sync" class="sync"></div>
|
||
</div>
|
||
</div>
|
||
<div id="splitbar" style="-moz-user-select:none;"
|
||
class="ui-resizable-handle">
|
||
</div>
|
||
</div>
|
||
<script type="text/javascript">
|
||
/* @license magnet:?xt=urn:btih:d3d9a9a6595521f9666a5e94cc830dab83b65699&dn=expat.txt MIT */
|
||
$(document).ready(function(){initNavTree('db/d40/integral__approximation2_8cpp.html','../../'); initResizable(); });
|
||
/* @license-end */
|
||
</script>
|
||
<div id="doc-content">
|
||
<!-- window showing the filter options -->
|
||
<div id="MSearchSelectWindow"
|
||
onmouseover="return searchBox.OnSearchSelectShow()"
|
||
onmouseout="return searchBox.OnSearchSelectHide()"
|
||
onkeydown="return searchBox.OnSearchSelectKey(event)">
|
||
</div>
|
||
|
||
<!-- iframe showing the search results (closed by default) -->
|
||
<div id="MSearchResultsWindow">
|
||
<iframe src="javascript:void(0)" frameborder="0"
|
||
name="MSearchResults" id="MSearchResults">
|
||
</iframe>
|
||
</div>
|
||
|
||
<div class="header">
|
||
<div class="summary">
|
||
<a href="#namespaces">Namespaces</a> |
|
||
<a href="#define-members">Macros</a> |
|
||
<a href="#typedef-members">Typedefs</a> |
|
||
<a href="#func-members">Functions</a> </div>
|
||
<div class="headertitle"><div class="title">integral_approximation2.cpp File Reference</div></div>
|
||
</div><!--header-->
|
||
<div class="contents">
|
||
|
||
<p><a href="https://en.wikipedia.org/wiki/Monte_Carlo_integration" target="_blank">Monte Carlo Integration</a>
|
||
<a href="#details">More...</a></p>
|
||
<div class="textblock"><code>#include <cmath></code><br />
|
||
<code>#include <cstdint></code><br />
|
||
<code>#include <ctime></code><br />
|
||
<code>#include <functional></code><br />
|
||
<code>#include <iostream></code><br />
|
||
<code>#include <random></code><br />
|
||
<code>#include <vector></code><br />
|
||
</div><div class="textblock"><div class="dynheader">
|
||
Include dependency graph for integral_approximation2.cpp:</div>
|
||
<div class="dyncontent">
|
||
<div class="center"><iframe scrolling="no" frameborder="0" src="../../d5/d5c/integral__approximation2_8cpp__incl.svg" width="607" height="112"><p><b>This browser is not able to show SVG: try Firefox, Chrome, Safari, or Opera instead.</b></p></iframe>
|
||
</div>
|
||
</div>
|
||
</div><table class="memberdecls">
|
||
<tr class="heading"><td colspan="2"><h2 class="groupheader"><a id="namespaces" name="namespaces"></a>
|
||
Namespaces</h2></td></tr>
|
||
<tr class="memitem:dd/d47/namespacemath"><td class="memItemLeft" align="right" valign="top">namespace  </td><td class="memItemRight" valign="bottom"><a class="el" href="../../dd/d47/namespacemath.html">math</a></td></tr>
|
||
<tr class="memdesc:dd/d47/namespacemath"><td class="mdescLeft"> </td><td class="mdescRight">for M_PI definition and pow() <br /></td></tr>
|
||
<tr class="separator:"><td class="memSeparator" colspan="2"> </td></tr>
|
||
<tr class="memitem:d0/da4/namespacemonte__carlo"><td class="memItemLeft" align="right" valign="top">namespace  </td><td class="memItemRight" valign="bottom"><a class="el" href="../../d0/da4/namespacemonte__carlo.html">monte_carlo</a></td></tr>
|
||
<tr class="memdesc:d0/da4/namespacemonte__carlo"><td class="mdescLeft"> </td><td class="mdescRight">Functions for the <a href="https://en.wikipedia.org/wiki/Monte_Carlo_integration" target="_blank">Monte Carlo Integration</a> implementation. <br /></td></tr>
|
||
<tr class="separator:"><td class="memSeparator" colspan="2"> </td></tr>
|
||
</table><table class="memberdecls">
|
||
<tr class="heading"><td colspan="2"><h2 class="groupheader"><a id="typedef-members" name="typedef-members"></a>
|
||
Typedefs</h2></td></tr>
|
||
<tr class="memitem:a56f3ff2501246c5e47bd2ac320a97446"><td class="memItemLeft" align="right" valign="top"><a id="a56f3ff2501246c5e47bd2ac320a97446" name="a56f3ff2501246c5e47bd2ac320a97446"></a>
|
||
using </td><td class="memItemRight" valign="bottom"><b>math::monte_carlo::Function</b> = <a class="elRef" target="_blank" href="http://en.cppreference.com/w/cpp/utility/functional/function.html">std::function</a>< double(double &)></td></tr>
|
||
<tr class="separator:a56f3ff2501246c5e47bd2ac320a97446"><td class="memSeparator" colspan="2"> </td></tr>
|
||
</table><table class="memberdecls">
|
||
<tr class="heading"><td colspan="2"><h2 class="groupheader"><a id="func-members" name="func-members"></a>
|
||
Functions</h2></td></tr>
|
||
<tr class="memitem:a71249ee535f16f8ed2e9cc8f0199a2cf"><td class="memItemLeft" align="right" valign="top"><a class="elRef" target="_blank" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector</a>< double > </td><td class="memItemRight" valign="bottom"><a class="el" href="../../db/d40/integral__approximation2_8cpp.html#a71249ee535f16f8ed2e9cc8f0199a2cf">math::monte_carlo::generate_samples</a> (const double &start_point, const Function &pdf, const uint32_t &num_samples, const uint32_t &discard=100000)</td></tr>
|
||
<tr class="memdesc:a71249ee535f16f8ed2e9cc8f0199a2cf"><td class="mdescLeft"> </td><td class="mdescRight">short-hand for std::functions used in this implementation <a href="../../db/d40/integral__approximation2_8cpp.html#a71249ee535f16f8ed2e9cc8f0199a2cf">More...</a><br /></td></tr>
|
||
<tr class="separator:a71249ee535f16f8ed2e9cc8f0199a2cf"><td class="memSeparator" colspan="2"> </td></tr>
|
||
<tr class="memitem:af7da9ba8932f1f48b9bbc2d80471af51"><td class="memItemLeft" align="right" valign="top">double </td><td class="memItemRight" valign="bottom"><a class="el" href="../../db/d40/integral__approximation2_8cpp.html#af7da9ba8932f1f48b9bbc2d80471af51">math::monte_carlo::integral_monte_carlo</a> (const double &start_point, const Function &function, const Function &pdf, const uint32_t &num_samples=1000000)</td></tr>
|
||
<tr class="memdesc:af7da9ba8932f1f48b9bbc2d80471af51"><td class="mdescLeft"> </td><td class="mdescRight">Compute an approximation of an integral using Monte Carlo integration. <a href="../../db/d40/integral__approximation2_8cpp.html#af7da9ba8932f1f48b9bbc2d80471af51">More...</a><br /></td></tr>
|
||
<tr class="separator:af7da9ba8932f1f48b9bbc2d80471af51"><td class="memSeparator" colspan="2"> </td></tr>
|
||
<tr class="memitem:aa8dca7b867074164d5f45b0f3851269d"><td class="memItemLeft" align="right" valign="top">static void </td><td class="memItemRight" valign="bottom"><a class="el" href="../../db/d40/integral__approximation2_8cpp.html#aa8dca7b867074164d5f45b0f3851269d">test</a> ()</td></tr>
|
||
<tr class="memdesc:aa8dca7b867074164d5f45b0f3851269d"><td class="mdescLeft"> </td><td class="mdescRight">Self-test implementations. <a href="../../db/d40/integral__approximation2_8cpp.html#aa8dca7b867074164d5f45b0f3851269d">More...</a><br /></td></tr>
|
||
<tr class="separator:aa8dca7b867074164d5f45b0f3851269d"><td class="memSeparator" colspan="2"> </td></tr>
|
||
<tr class="memitem:ae66f6b31b5ad750f1fe042a706a4e3d4"><td class="memItemLeft" align="right" valign="top">int </td><td class="memItemRight" valign="bottom"><a class="el" href="../../db/d40/integral__approximation2_8cpp.html#ae66f6b31b5ad750f1fe042a706a4e3d4">main</a> ()</td></tr>
|
||
<tr class="memdesc:ae66f6b31b5ad750f1fe042a706a4e3d4"><td class="mdescLeft"> </td><td class="mdescRight">Main function. <a href="../../db/d40/integral__approximation2_8cpp.html#ae66f6b31b5ad750f1fe042a706a4e3d4">More...</a><br /></td></tr>
|
||
<tr class="separator:ae66f6b31b5ad750f1fe042a706a4e3d4"><td class="memSeparator" colspan="2"> </td></tr>
|
||
</table>
|
||
<a name="details" id="details"></a><h2 class="groupheader">Detailed Description</h2>
|
||
<div class="textblock"><p ><a href="https://en.wikipedia.org/wiki/Monte_Carlo_integration" target="_blank">Monte Carlo Integration</a> </p>
|
||
<p >In mathematics, Monte Carlo integration is a technique for numerical integration using random numbers. It is a particular Monte Carlo method that numerically computes a definite integral. While other algorithms usually evaluate the integrand at a regular grid, Monte Carlo randomly chooses points at which the integrand is evaluated. This method is particularly useful for higher-dimensional integrals.</p>
|
||
<p >This implementation supports arbitrary pdfs. These pdfs are sampled using the <a href="https://en.wikipedia.org/wiki/Metropolis–Hastings_algorithm" target="_blank">Metropolis-Hastings algorithm</a>. This can be swapped out by every other sampling techniques for example the inverse method. Metropolis-Hastings was chosen because it is the most general and can also be extended for a higher dimensional sampling space.</p>
|
||
<dl class="section author"><dt>Author</dt><dd><a href="https://github.com/DerAndereDomenic" target="_blank">Domenic Zingsheim</a> </dd></dl>
|
||
</div><h2 class="groupheader">Function Documentation</h2>
|
||
<a id="a71249ee535f16f8ed2e9cc8f0199a2cf" name="a71249ee535f16f8ed2e9cc8f0199a2cf"></a>
|
||
<h2 class="memtitle"><span class="permalink"><a href="#a71249ee535f16f8ed2e9cc8f0199a2cf">◆ </a></span>generate_samples()</h2>
|
||
|
||
<div class="memitem">
|
||
<div class="memproto">
|
||
<table class="memname">
|
||
<tr>
|
||
<td class="memname"><a class="elRef" target="_blank" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector</a>< double > math::monte_carlo::generate_samples </td>
|
||
<td>(</td>
|
||
<td class="paramtype">const double & </td>
|
||
<td class="paramname"><em>start_point</em>, </td>
|
||
</tr>
|
||
<tr>
|
||
<td class="paramkey"></td>
|
||
<td></td>
|
||
<td class="paramtype">const <a class="elRef" target="_blank" href="http://en.cppreference.com/w/cpp/utility/functional/function.html">Function</a> & </td>
|
||
<td class="paramname"><em>pdf</em>, </td>
|
||
</tr>
|
||
<tr>
|
||
<td class="paramkey"></td>
|
||
<td></td>
|
||
<td class="paramtype">const uint32_t & </td>
|
||
<td class="paramname"><em>num_samples</em>, </td>
|
||
</tr>
|
||
<tr>
|
||
<td class="paramkey"></td>
|
||
<td></td>
|
||
<td class="paramtype">const uint32_t & </td>
|
||
<td class="paramname"><em>discard</em> = <code>100000</code> </td>
|
||
</tr>
|
||
<tr>
|
||
<td></td>
|
||
<td>)</td>
|
||
<td></td><td></td>
|
||
</tr>
|
||
</table>
|
||
</div><div class="memdoc">
|
||
|
||
<p>short-hand for std::functions used in this implementation </p>
|
||
<p >Generate samples according to some pdf</p>
|
||
<p >This function uses Metropolis-Hastings to generate random numbers. It generates a sequence of random numbers by using a markov chain. Therefore, we need to define a start_point and the number of samples we want to generate. Because the first samples generated by the markov chain may not be distributed according to the given pdf, one can specify how many samples should be discarded before storing samples. </p><dl class="params"><dt>Parameters</dt><dd>
|
||
<table class="params">
|
||
<tr><td class="paramname">start_point</td><td>The starting point of the markov chain </td></tr>
|
||
<tr><td class="paramname">pdf</td><td>The pdf to sample </td></tr>
|
||
<tr><td class="paramname">num_samples</td><td>The number of samples to generate </td></tr>
|
||
<tr><td class="paramname">discard</td><td>How many samples should be discarded at the start </td></tr>
|
||
</table>
|
||
</dd>
|
||
</dl>
|
||
<dl class="section return"><dt>Returns</dt><dd>A vector of size num_samples with samples distributed according to the pdf </dd></dl>
|
||
<div class="fragment"><div class="line"><span class="lineno"> 67</span> {</div>
|
||
<div class="line"><span class="lineno"> 68</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector<double></a> samples;</div>
|
||
<div class="line"><span class="lineno"> 69</span> samples.<a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/container/vector/reserve.html">reserve</a>(num_samples);</div>
|
||
<div class="line"><span class="lineno"> 70</span> </div>
|
||
<div class="line"><span class="lineno"> 71</span> <span class="keywordtype">double</span> x_t = start_point;</div>
|
||
<div class="line"><span class="lineno"> 72</span> </div>
|
||
<div class="line"><span class="lineno"> 73</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/numeric/random.html">std::default_random_engine</a> generator;</div>
|
||
<div class="line"><span class="lineno"> 74</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/numeric/random/uniform_real_distribution.html">std::uniform_real_distribution<double></a> uniform(0.0, 1.0);</div>
|
||
<div class="line"><span class="lineno"> 75</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/numeric/random/normal_distribution.html">std::normal_distribution<double></a> normal(0.0, 1.0);</div>
|
||
<div class="line"><span class="lineno"> 76</span> generator.seed(<a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/chrono/c/time.html">time</a>(<span class="keyword">nullptr</span>));</div>
|
||
<div class="line"><span class="lineno"> 77</span> </div>
|
||
<div class="line"><span class="lineno"> 78</span> <span class="keywordflow">for</span> (uint32_t t = 0; t < num_samples + discard; ++t) {</div>
|
||
<div class="line"><span class="lineno"> 79</span> <span class="comment">// Generate a new proposal according to some mutation strategy.</span></div>
|
||
<div class="line"><span class="lineno"> 80</span> <span class="comment">// This is arbitrary and can be swapped.</span></div>
|
||
<div class="line"><span class="lineno"> 81</span> <span class="keywordtype">double</span> x_dash = normal(generator) + x_t;</div>
|
||
<div class="line"><span class="lineno"> 82</span> <span class="keywordtype">double</span> acceptance_probability = <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/algorithm/min.html">std::min</a>(pdf(x_dash) / pdf(x_t), 1.0);</div>
|
||
<div class="line"><span class="lineno"> 83</span> <span class="keywordtype">double</span> u = uniform(generator);</div>
|
||
<div class="line"><span class="lineno"> 84</span> </div>
|
||
<div class="line"><span class="lineno"> 85</span> <span class="comment">// Accept "new state" according to the acceptance_probability</span></div>
|
||
<div class="line"><span class="lineno"> 86</span> <span class="keywordflow">if</span> (u <= acceptance_probability) {</div>
|
||
<div class="line"><span class="lineno"> 87</span> x_t = x_dash;</div>
|
||
<div class="line"><span class="lineno"> 88</span> }</div>
|
||
<div class="line"><span class="lineno"> 89</span> </div>
|
||
<div class="line"><span class="lineno"> 90</span> <span class="keywordflow">if</span> (t >= discard) {</div>
|
||
<div class="line"><span class="lineno"> 91</span> samples.<a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/container/vector/push_back.html">push_back</a>(x_t);</div>
|
||
<div class="line"><span class="lineno"> 92</span> }</div>
|
||
<div class="line"><span class="lineno"> 93</span> }</div>
|
||
<div class="line"><span class="lineno"> 94</span> </div>
|
||
<div class="line"><span class="lineno"> 95</span> <span class="keywordflow">return</span> samples;</div>
|
||
<div class="line"><span class="lineno"> 96</span>}</div>
|
||
<div class="ttc" id="amin_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/algorithm/min.html">std::min</a></div><div class="ttdeci">T min(T... args)</div></div>
|
||
<div class="ttc" id="anormal_distribution_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/numeric/random/normal_distribution.html">std::normal_distribution</a></div></div>
|
||
<div class="ttc" id="apush_back_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/container/vector/push_back.html">std::vector::push_back</a></div><div class="ttdeci">T push_back(T... args)</div></div>
|
||
<div class="ttc" id="arandom_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/numeric/random.html">std::default_random_engine</a></div></div>
|
||
<div class="ttc" id="areserve_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/container/vector/reserve.html">std::vector::reserve</a></div><div class="ttdeci">T reserve(T... args)</div></div>
|
||
<div class="ttc" id="atime_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/chrono/c/time.html">std::time</a></div><div class="ttdeci">T time(T... args)</div></div>
|
||
<div class="ttc" id="auniform_real_distribution_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/numeric/random/uniform_real_distribution.html">std::uniform_real_distribution</a></div></div>
|
||
<div class="ttc" id="avector_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector< double ></a></div></div>
|
||
</div><!-- fragment --><div class="dynheader">
|
||
Here is the call graph for this function:</div>
|
||
<div class="dyncontent">
|
||
<div class="center"><iframe scrolling="no" frameborder="0" src="../../db/d40/integral__approximation2_8cpp_a71249ee535f16f8ed2e9cc8f0199a2cf_cgraph.svg" width="352" height="139"><p><b>This browser is not able to show SVG: try Firefox, Chrome, Safari, or Opera instead.</b></p></iframe>
|
||
</div>
|
||
</div>
|
||
|
||
</div>
|
||
</div>
|
||
<a id="af7da9ba8932f1f48b9bbc2d80471af51" name="af7da9ba8932f1f48b9bbc2d80471af51"></a>
|
||
<h2 class="memtitle"><span class="permalink"><a href="#af7da9ba8932f1f48b9bbc2d80471af51">◆ </a></span>integral_monte_carlo()</h2>
|
||
|
||
<div class="memitem">
|
||
<div class="memproto">
|
||
<table class="memname">
|
||
<tr>
|
||
<td class="memname">double math::monte_carlo::integral_monte_carlo </td>
|
||
<td>(</td>
|
||
<td class="paramtype">const double & </td>
|
||
<td class="paramname"><em>start_point</em>, </td>
|
||
</tr>
|
||
<tr>
|
||
<td class="paramkey"></td>
|
||
<td></td>
|
||
<td class="paramtype">const <a class="elRef" target="_blank" href="http://en.cppreference.com/w/cpp/utility/functional/function.html">Function</a> & </td>
|
||
<td class="paramname"><em>function</em>, </td>
|
||
</tr>
|
||
<tr>
|
||
<td class="paramkey"></td>
|
||
<td></td>
|
||
<td class="paramtype">const <a class="elRef" target="_blank" href="http://en.cppreference.com/w/cpp/utility/functional/function.html">Function</a> & </td>
|
||
<td class="paramname"><em>pdf</em>, </td>
|
||
</tr>
|
||
<tr>
|
||
<td class="paramkey"></td>
|
||
<td></td>
|
||
<td class="paramtype">const uint32_t & </td>
|
||
<td class="paramname"><em>num_samples</em> = <code>1000000</code> </td>
|
||
</tr>
|
||
<tr>
|
||
<td></td>
|
||
<td>)</td>
|
||
<td></td><td></td>
|
||
</tr>
|
||
</table>
|
||
</div><div class="memdoc">
|
||
|
||
<p>Compute an approximation of an integral using Monte Carlo integration. </p>
|
||
<p >The integration domain [a,b] is given by the pdf. The pdf has to fulfill the following conditions: 1) for all x \in [a,b] : p(x) > 0 2) for all x \not\in [a,b] : p(x) = 0 3) \int_a^b p(x) dx = 1 </p><dl class="params"><dt>Parameters</dt><dd>
|
||
<table class="params">
|
||
<tr><td class="paramname">start_point</td><td>The start point of the Markov Chain (see generate_samples) </td></tr>
|
||
<tr><td class="paramname">function</td><td>The function to integrate </td></tr>
|
||
<tr><td class="paramname">pdf</td><td>The pdf to sample </td></tr>
|
||
<tr><td class="paramname">num_samples</td><td>The number of samples used to approximate the integral </td></tr>
|
||
</table>
|
||
</dd>
|
||
</dl>
|
||
<dl class="section return"><dt>Returns</dt><dd>The approximation of the integral according to 1/N \sum_{i}^N f(x_i) / p(x_i) </dd></dl>
|
||
<div class="fragment"><div class="line"><span class="lineno"> 114</span> {</div>
|
||
<div class="line"><span class="lineno"> 115</span> <span class="keywordtype">double</span> integral = 0.0;</div>
|
||
<div class="line"><span class="lineno"> 116</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/container/vector.html">std::vector<double></a> samples =</div>
|
||
<div class="line"><span class="lineno"> 117</span> <a class="code hl_function" href="../../db/d40/integral__approximation2_8cpp.html#a71249ee535f16f8ed2e9cc8f0199a2cf">generate_samples</a>(start_point, pdf, num_samples);</div>
|
||
<div class="line"><span class="lineno"> 118</span> </div>
|
||
<div class="line"><span class="lineno"> 119</span> <span class="keywordflow">for</span> (<span class="keywordtype">double</span> sample : samples) {</div>
|
||
<div class="line"><span class="lineno"> 120</span> integral += function(sample) / pdf(sample);</div>
|
||
<div class="line"><span class="lineno"> 121</span> }</div>
|
||
<div class="line"><span class="lineno"> 122</span> </div>
|
||
<div class="line"><span class="lineno"> 123</span> <span class="keywordflow">return</span> integral / <span class="keyword">static_cast<</span><span class="keywordtype">double</span><span class="keyword">></span>(samples.size());</div>
|
||
<div class="line"><span class="lineno"> 124</span>}</div>
|
||
<div class="ttc" id="aintegral__approximation2_8cpp_html_a71249ee535f16f8ed2e9cc8f0199a2cf"><div class="ttname"><a href="../../db/d40/integral__approximation2_8cpp.html#a71249ee535f16f8ed2e9cc8f0199a2cf">math::monte_carlo::generate_samples</a></div><div class="ttdeci">std::vector< double > generate_samples(const double &start_point, const Function &pdf, const uint32_t &num_samples, const uint32_t &discard=100000)</div><div class="ttdoc">short-hand for std::functions used in this implementation</div><div class="ttdef"><b>Definition:</b> integral_approximation2.cpp:64</div></div>
|
||
</div><!-- fragment --><div class="dynheader">
|
||
Here is the call graph for this function:</div>
|
||
<div class="dyncontent">
|
||
<div class="center"><iframe scrolling="no" frameborder="0" src="../../db/d40/integral__approximation2_8cpp_af7da9ba8932f1f48b9bbc2d80471af51_cgraph.svg" width="543" height="147"><p><b>This browser is not able to show SVG: try Firefox, Chrome, Safari, or Opera instead.</b></p></iframe>
|
||
</div>
|
||
</div>
|
||
|
||
</div>
|
||
</div>
|
||
<a id="ae66f6b31b5ad750f1fe042a706a4e3d4" name="ae66f6b31b5ad750f1fe042a706a4e3d4"></a>
|
||
<h2 class="memtitle"><span class="permalink"><a href="#ae66f6b31b5ad750f1fe042a706a4e3d4">◆ </a></span>main()</h2>
|
||
|
||
<div class="memitem">
|
||
<div class="memproto">
|
||
<table class="memname">
|
||
<tr>
|
||
<td class="memname">int main </td>
|
||
<td>(</td>
|
||
<td class="paramtype">void </td>
|
||
<td class="paramname"></td><td>)</td>
|
||
<td></td>
|
||
</tr>
|
||
</table>
|
||
</div><div class="memdoc">
|
||
|
||
<p>Main function. </p>
|
||
<dl class="section return"><dt>Returns</dt><dd>0 on exit </dd></dl>
|
||
<div class="fragment"><div class="line"><span class="lineno"> 215</span> {</div>
|
||
<div class="line"><span class="lineno"> 216</span> <a class="code hl_function" href="../../db/d40/integral__approximation2_8cpp.html#aa8dca7b867074164d5f45b0f3851269d">test</a>(); <span class="comment">// run self-test implementations</span></div>
|
||
<div class="line"><span class="lineno"> 217</span> <span class="keywordflow">return</span> 0;</div>
|
||
<div class="line"><span class="lineno"> 218</span>}</div>
|
||
<div class="ttc" id="aintegral__approximation2_8cpp_html_aa8dca7b867074164d5f45b0f3851269d"><div class="ttname"><a href="../../db/d40/integral__approximation2_8cpp.html#aa8dca7b867074164d5f45b0f3851269d">test</a></div><div class="ttdeci">static void test()</div><div class="ttdoc">Self-test implementations.</div><div class="ttdef"><b>Definition:</b> integral_approximation2.cpp:133</div></div>
|
||
</div><!-- fragment --><div class="dynheader">
|
||
Here is the call graph for this function:</div>
|
||
<div class="dyncontent">
|
||
<div class="center"><iframe scrolling="no" frameborder="0" src="../../db/d40/integral__approximation2_8cpp_ae66f6b31b5ad750f1fe042a706a4e3d4_cgraph.svg" width="274" height="190"><p><b>This browser is not able to show SVG: try Firefox, Chrome, Safari, or Opera instead.</b></p></iframe>
|
||
</div>
|
||
</div>
|
||
|
||
</div>
|
||
</div>
|
||
<a id="aa8dca7b867074164d5f45b0f3851269d" name="aa8dca7b867074164d5f45b0f3851269d"></a>
|
||
<h2 class="memtitle"><span class="permalink"><a href="#aa8dca7b867074164d5f45b0f3851269d">◆ </a></span>test()</h2>
|
||
|
||
<div class="memitem">
|
||
<div class="memproto">
|
||
<table class="mlabels">
|
||
<tr>
|
||
<td class="mlabels-left">
|
||
<table class="memname">
|
||
<tr>
|
||
<td class="memname">static void test </td>
|
||
<td>(</td>
|
||
<td class="paramname"></td><td>)</td>
|
||
<td></td>
|
||
</tr>
|
||
</table>
|
||
</td>
|
||
<td class="mlabels-right">
|
||
<span class="mlabels"><span class="mlabel">static</span></span> </td>
|
||
</tr>
|
||
</table>
|
||
</div><div class="memdoc">
|
||
|
||
<p>Self-test implementations. </p>
|
||
<dl class="section return"><dt>Returns</dt><dd>void </dd></dl>
|
||
<div class="fragment"><div class="line"><span class="lineno"> 133</span> {</div>
|
||
<div class="line"><span class="lineno"> 134</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a> << <span class="stringliteral">"Disclaimer: Because this is a randomized algorithm,"</span></div>
|
||
<div class="line"><span class="lineno"> 135</span> << <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/manip/endl.html">std::endl</a>;</div>
|
||
<div class="line"><span class="lineno"> 136</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a></div>
|
||
<div class="line"><span class="lineno"> 137</span> << <span class="stringliteral">"it may happen that singular samples deviate from the true result."</span></div>
|
||
<div class="line"><span class="lineno"> 138</span> << <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/manip/endl.html">std::endl</a></div>
|
||
<div class="line"><span class="lineno"> 139</span> << <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/manip/endl.html">std::endl</a>;</div>
|
||
<div class="line"><span class="lineno"> 140</span> ;</div>
|
||
<div class="line"><span class="lineno"> 141</span> </div>
|
||
<div class="line"><span class="lineno"> 142</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/utility/functional/function.html">math::monte_carlo::Function</a> <a class="code hl_function" href="../../d4/d18/composite__simpson__rule_8cpp.html#a4251b4df4748a0b9c43a48f61bdd2397">f</a>;</div>
|
||
<div class="line"><span class="lineno"> 143</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/utility/functional/function.html">math::monte_carlo::Function</a> pdf;</div>
|
||
<div class="line"><span class="lineno"> 144</span> <span class="keywordtype">double</span> integral = 0;</div>
|
||
<div class="line"><span class="lineno"> 145</span> <span class="keywordtype">double</span> <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/algorithm/lower_bound.html">lower_bound</a> = 0, <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/algorithm/upper_bound.html">upper_bound</a> = 0;</div>
|
||
<div class="line"><span class="lineno"> 146</span> </div>
|
||
<div class="line"><span class="lineno"> 147</span> <span class="comment">/* \int_{-2}^{2} -x^2 + 4 dx */</span></div>
|
||
<div class="line"><span class="lineno"> 148</span> <a class="code hl_function" href="../../d4/d18/composite__simpson__rule_8cpp.html#a4251b4df4748a0b9c43a48f61bdd2397">f</a> = [&](<span class="keywordtype">double</span>& x) { <span class="keywordflow">return</span> -x * x + 4.0; };</div>
|
||
<div class="line"><span class="lineno"> 149</span> </div>
|
||
<div class="line"><span class="lineno"> 150</span> <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/algorithm/lower_bound.html">lower_bound</a> = -2.0;</div>
|
||
<div class="line"><span class="lineno"> 151</span> <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/algorithm/upper_bound.html">upper_bound</a> = 2.0;</div>
|
||
<div class="line"><span class="lineno"> 152</span> pdf = [&](<span class="keywordtype">double</span>& x) {</div>
|
||
<div class="line"><span class="lineno"> 153</span> <span class="keywordflow">if</span> (x >= lower_bound && x <= -1.0) {</div>
|
||
<div class="line"><span class="lineno"> 154</span> <span class="keywordflow">return</span> 0.1;</div>
|
||
<div class="line"><span class="lineno"> 155</span> }</div>
|
||
<div class="line"><span class="lineno"> 156</span> <span class="keywordflow">if</span> (x <= upper_bound && x >= 1.0) {</div>
|
||
<div class="line"><span class="lineno"> 157</span> <span class="keywordflow">return</span> 0.1;</div>
|
||
<div class="line"><span class="lineno"> 158</span> }</div>
|
||
<div class="line"><span class="lineno"> 159</span> <span class="keywordflow">if</span> (x > -1.0 && x < 1.0) {</div>
|
||
<div class="line"><span class="lineno"> 160</span> <span class="keywordflow">return</span> 0.4;</div>
|
||
<div class="line"><span class="lineno"> 161</span> }</div>
|
||
<div class="line"><span class="lineno"> 162</span> <span class="keywordflow">return</span> 0.0;</div>
|
||
<div class="line"><span class="lineno"> 163</span> };</div>
|
||
<div class="line"><span class="lineno"> 164</span> </div>
|
||
<div class="line"><span class="lineno"> 165</span> integral = <a class="code hl_function" href="../../db/d40/integral__approximation2_8cpp.html#af7da9ba8932f1f48b9bbc2d80471af51">math::monte_carlo::integral_monte_carlo</a>(</div>
|
||
<div class="line"><span class="lineno"> 166</span> (upper_bound - lower_bound) / 2.0, f, pdf);</div>
|
||
<div class="line"><span class="lineno"> 167</span> </div>
|
||
<div class="line"><span class="lineno"> 168</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a> << <span class="stringliteral">"This number should be close to 10.666666: "</span> << integral</div>
|
||
<div class="line"><span class="lineno"> 169</span> << <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/manip/endl.html">std::endl</a>;</div>
|
||
<div class="line"><span class="lineno"> 170</span> </div>
|
||
<div class="line"><span class="lineno"> 171</span> <span class="comment">/* \int_{0}^{1} e^x dx */</span></div>
|
||
<div class="line"><span class="lineno"> 172</span> <a class="code hl_function" href="../../d4/d18/composite__simpson__rule_8cpp.html#a4251b4df4748a0b9c43a48f61bdd2397">f</a> = [&](<span class="keywordtype">double</span>& x) { <span class="keywordflow">return</span> <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/numeric/math/exp.html">std::exp</a>(x); };</div>
|
||
<div class="line"><span class="lineno"> 173</span> </div>
|
||
<div class="line"><span class="lineno"> 174</span> <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/algorithm/lower_bound.html">lower_bound</a> = 0.0;</div>
|
||
<div class="line"><span class="lineno"> 175</span> <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/algorithm/upper_bound.html">upper_bound</a> = 1.0;</div>
|
||
<div class="line"><span class="lineno"> 176</span> pdf = [&](<span class="keywordtype">double</span>& x) {</div>
|
||
<div class="line"><span class="lineno"> 177</span> <span class="keywordflow">if</span> (x >= lower_bound && x <= 0.2) {</div>
|
||
<div class="line"><span class="lineno"> 178</span> <span class="keywordflow">return</span> 0.1;</div>
|
||
<div class="line"><span class="lineno"> 179</span> }</div>
|
||
<div class="line"><span class="lineno"> 180</span> <span class="keywordflow">if</span> (x > 0.2 && x <= 0.4) {</div>
|
||
<div class="line"><span class="lineno"> 181</span> <span class="keywordflow">return</span> 0.4;</div>
|
||
<div class="line"><span class="lineno"> 182</span> }</div>
|
||
<div class="line"><span class="lineno"> 183</span> <span class="keywordflow">if</span> (x > 0.4 && x < upper_bound) {</div>
|
||
<div class="line"><span class="lineno"> 184</span> <span class="keywordflow">return</span> 1.5;</div>
|
||
<div class="line"><span class="lineno"> 185</span> }</div>
|
||
<div class="line"><span class="lineno"> 186</span> <span class="keywordflow">return</span> 0.0;</div>
|
||
<div class="line"><span class="lineno"> 187</span> };</div>
|
||
<div class="line"><span class="lineno"> 188</span> </div>
|
||
<div class="line"><span class="lineno"> 189</span> integral = <a class="code hl_function" href="../../db/d40/integral__approximation2_8cpp.html#af7da9ba8932f1f48b9bbc2d80471af51">math::monte_carlo::integral_monte_carlo</a>(</div>
|
||
<div class="line"><span class="lineno"> 190</span> (upper_bound - lower_bound) / 2.0, f, pdf);</div>
|
||
<div class="line"><span class="lineno"> 191</span> </div>
|
||
<div class="line"><span class="lineno"> 192</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a> << <span class="stringliteral">"This number should be close to 1.7182818: "</span> << integral</div>
|
||
<div class="line"><span class="lineno"> 193</span> << <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/manip/endl.html">std::endl</a>;</div>
|
||
<div class="line"><span class="lineno"> 194</span> </div>
|
||
<div class="line"><span class="lineno"> 195</span> <span class="comment">/* \int_{-\infty}^{\infty} sinc(x) dx, sinc(x) = sin(pi * x) / (pi * x)</span></div>
|
||
<div class="line"><span class="lineno"> 196</span><span class="comment"> This is a difficult integral because of its infinite domain.</span></div>
|
||
<div class="line"><span class="lineno"> 197</span><span class="comment"> Therefore, it may deviate largely from the expected result.</span></div>
|
||
<div class="line"><span class="lineno"> 198</span><span class="comment"> */</span></div>
|
||
<div class="line"><span class="lineno"> 199</span> <a class="code hl_function" href="../../d4/d18/composite__simpson__rule_8cpp.html#a4251b4df4748a0b9c43a48f61bdd2397">f</a> = [&](<span class="keywordtype">double</span>& x) { <span class="keywordflow">return</span> <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/numeric/math/sin.html">std::sin</a>(M_PI * x) / (M_PI * x); };</div>
|
||
<div class="line"><span class="lineno"> 200</span> </div>
|
||
<div class="line"><span class="lineno"> 201</span> pdf = [&](<span class="keywordtype">double</span>& x) {</div>
|
||
<div class="line"><span class="lineno"> 202</span> <span class="keywordflow">return</span> 1.0 / <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/numeric/math/sqrt.html">std::sqrt</a>(2.0 * M_PI) * <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/numeric/math/exp.html">std::exp</a>(-x * x / 2.0);</div>
|
||
<div class="line"><span class="lineno"> 203</span> };</div>
|
||
<div class="line"><span class="lineno"> 204</span> </div>
|
||
<div class="line"><span class="lineno"> 205</span> integral = <a class="code hl_function" href="../../db/d40/integral__approximation2_8cpp.html#af7da9ba8932f1f48b9bbc2d80471af51">math::monte_carlo::integral_monte_carlo</a>(0.0, f, pdf, 10000000);</div>
|
||
<div class="line"><span class="lineno"> 206</span> </div>
|
||
<div class="line"><span class="lineno"> 207</span> <a class="code hl_classRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a> << <span class="stringliteral">"This number should be close to 1.0: "</span> << integral</div>
|
||
<div class="line"><span class="lineno"> 208</span> << <a class="code hl_functionRef" target="_blank" href="http://en.cppreference.com/w/cpp/io/manip/endl.html">std::endl</a>;</div>
|
||
<div class="line"><span class="lineno"> 209</span>}</div>
|
||
<div class="ttc" id="abasic_ostream_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/io/basic_ostream.html">std::cout</a></div></div>
|
||
<div class="ttc" id="acomposite__simpson__rule_8cpp_html_a4251b4df4748a0b9c43a48f61bdd2397"><div class="ttname"><a href="../../d4/d18/composite__simpson__rule_8cpp.html#a4251b4df4748a0b9c43a48f61bdd2397">numerical_methods::simpson_method::f</a></div><div class="ttdeci">double f(double x)</div><div class="ttdoc">A function f(x) that will be used to test the method.</div><div class="ttdef"><b>Definition:</b> composite_simpson_rule.cpp:113</div></div>
|
||
<div class="ttc" id="aendl_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/io/manip/endl.html">std::endl</a></div><div class="ttdeci">T endl(T... args)</div></div>
|
||
<div class="ttc" id="aexp_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/numeric/math/exp.html">std::exp</a></div><div class="ttdeci">T exp(T... args)</div></div>
|
||
<div class="ttc" id="afunction_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/utility/functional/function.html">std::function</a></div></div>
|
||
<div class="ttc" id="aintegral__approximation2_8cpp_html_af7da9ba8932f1f48b9bbc2d80471af51"><div class="ttname"><a href="../../db/d40/integral__approximation2_8cpp.html#af7da9ba8932f1f48b9bbc2d80471af51">math::monte_carlo::integral_monte_carlo</a></div><div class="ttdeci">double integral_monte_carlo(const double &start_point, const Function &function, const Function &pdf, const uint32_t &num_samples=1000000)</div><div class="ttdoc">Compute an approximation of an integral using Monte Carlo integration.</div><div class="ttdef"><b>Definition:</b> integral_approximation2.cpp:112</div></div>
|
||
<div class="ttc" id="alower_bound_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/algorithm/lower_bound.html">std::lower_bound</a></div><div class="ttdeci">T lower_bound(T... args)</div></div>
|
||
<div class="ttc" id="asin_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/numeric/math/sin.html">std::sin</a></div><div class="ttdeci">T sin(T... args)</div></div>
|
||
<div class="ttc" id="asqrt_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/numeric/math/sqrt.html">std::sqrt</a></div><div class="ttdeci">T sqrt(T... args)</div></div>
|
||
<div class="ttc" id="aupper_bound_html"><div class="ttname"><a href="http://en.cppreference.com/w/cpp/algorithm/upper_bound.html">std::upper_bound</a></div><div class="ttdeci">T upper_bound(T... args)</div></div>
|
||
</div><!-- fragment --><div class="dynheader">
|
||
Here is the call graph for this function:</div>
|
||
<div class="dyncontent">
|
||
<div class="center"><iframe scrolling="no" frameborder="0" src="../../db/d40/integral__approximation2_8cpp_aa8dca7b867074164d5f45b0f3851269d_cgraph.svg" width="175" height="190"><p><b>This browser is not able to show SVG: try Firefox, Chrome, Safari, or Opera instead.</b></p></iframe>
|
||
</div>
|
||
</div>
|
||
|
||
</div>
|
||
</div>
|
||
</div><!-- contents -->
|
||
</div><!-- doc-content -->
|
||
<!-- start footer part -->
|
||
<div id="nav-path" class="navpath"><!-- id is needed for treeview function! -->
|
||
<ul>
|
||
<li class="navelem"><a class="el" href="../../dir_296d53ceaeaa7e099814a6def439fe8a.html">math</a></li><li class="navelem"><a class="el" href="../../db/d40/integral__approximation2_8cpp.html">integral_approximation2.cpp</a></li>
|
||
<li class="footer">Generated by <a href="https://www.doxygen.org/index.html"><img class="footer" src="../../doxygen.svg" width="104" height="31" alt="doxygen"/></a> 1.9.3 </li>
|
||
</ul>
|
||
</div>
|
||
</body>
|
||
</html>
|