csharpfftfsharpintegrationinterpolationlinear-algebramathdifferentiationmatrixnumericsrandomregressionstatisticsmathnet
You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
246 lines
21 KiB
246 lines
21 KiB
<!DOCTYPE html>
|
|
<html lang="en">
|
|
|
|
<head>
|
|
<meta charset="utf-8">
|
|
<title>Numerical Integration
|
|
</title>
|
|
<meta name="viewport" content="width=device-width, initial-scale=1, shrink-to-fit=no">
|
|
<meta name="author" content="Christoph Ruegg, Marcus Cuda, Jurgen Van Gael">
|
|
|
|
<link rel="stylesheet" id="theme_link" href="https://cdnjs.cloudflare.com/ajax/libs/bootswatch/4.6.0/materia/bootstrap.min.css">
|
|
<script src="https://code.jquery.com/jquery-3.4.1.min.js"></script>
|
|
<script src="https://cdn.jsdelivr.net/npm/bootstrap@4.6.0/dist/js/bootstrap.bundle.min.js" integrity="sha384-Piv4xVNRyMGpqkS2by6br4gNJ7DXjqk09RmUpJ8jgGtD7zP9yug3goQfGII0yAns" crossorigin="anonymous"></script>
|
|
|
|
<script type="text/javascript" async src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.7/MathJax.js?config=TeX-MML-AM_CHTML"></script>
|
|
|
|
<link rel="shortcut icon" type="image/x-icon" href="/favicon.ico">
|
|
<link type="text/css" rel="stylesheet" href="https://numerics.mathdotnet.com/content/navbar-fixed-left.css" />
|
|
<link type="text/css" rel="stylesheet" href="https://numerics.mathdotnet.com/content/fsdocs-default.css" />
|
|
<link type="text/css" rel="stylesheet" href="https://numerics.mathdotnet.com/content/fsdocs-custom.css" />
|
|
<script type="text/javascript" src="https://numerics.mathdotnet.com/content/fsdocs-tips.js"></script>
|
|
<!-- HTML5 shim, for IE6-8 support of HTML5 elements -->
|
|
<!--[if lt IE 9]>
|
|
<script src="http://html5shim.googlecode.com/svn/trunk/html5.js"></script>
|
|
<![endif]-->
|
|
<!-- BEGIN SEARCH BOX: this adds support for the search box -->
|
|
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/JavaScript-autoComplete/1.0.4/auto-complete.css" />
|
|
<!-- END SEARCH BOX: this adds support for the search box -->
|
|
|
|
</head>
|
|
|
|
<body>
|
|
<nav class="navbar navbar-expand-md navbar-light bg-secondary fixed-left" id="fsdocs-nav">
|
|
<button class="navbar-toggler" type="button" data-toggle="collapse" data-target="#navbarsExampleDefault" aria-controls="navbarsExampleDefault" aria-expanded="false" aria-label="Toggle navigation">
|
|
<span class="navbar-toggler-icon"></span>
|
|
</button>
|
|
<div class="collapse navbar-collapse navbar-nav-scroll" id="navbarsExampleDefault">
|
|
<a href="https://numerics.mathdotnet.com/"><img id="fsdocs-logo" src="/logo.png" /></a>
|
|
<!-- BEGIN SEARCH BOX: this adds support for the search box -->
|
|
<div id="header">
|
|
<div class="searchbox" id="fsdocs-searchbox">
|
|
<label for="search-by">
|
|
<i class="fas fa-search"></i>
|
|
</label>
|
|
<input data-search-input="" id="search-by" type="search" placeholder="Search..." />
|
|
<span data-search-clear="">
|
|
<i class="fas fa-times"></i>
|
|
</span>
|
|
</div>
|
|
</div>
|
|
<!-- END SEARCH BOX: this adds support for the search box -->
|
|
<ul class="navbar-nav">
|
|
<li class="nav-header">Math.NET Numerics</li>
|
|
<li class="nav-item"><a class="nav-link" href="Packages.html">NuGet & Binaries</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="ReleaseNotes.html">Release Notes</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="https://github.com/mathnet/mathnet-numerics/blob/master/LICENSE.md">MIT License</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Compatibility.html">Platform Support</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="https://numerics.mathdotnet.com/api/">Class Reference</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="https://github.com/mathnet/mathnet-numerics/issues">Issues & Bugs</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Users.html">Who is using Math.NET?</a></li>
|
|
|
|
<li class="nav-header">Contributing</li>
|
|
<li class="nav-item"><a class="nav-link" href="Contributors.html">Contributors</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Contributing.html">Contributing</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Build.html">Build & Tools</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="https://github.com/mathnet/mathnet-numerics/discussions/categories/ideas">Your Ideas</a></li>
|
|
|
|
<li class="nav-header">Getting Help</li>
|
|
<li class="nav-item"><a class="nav-link" href="https://discuss.mathdotnet.com/c/numerics">Discuss</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="https://stackoverflow.com/questions/tagged/mathdotnet">Stack Overflow</a></li>
|
|
|
|
<li class="nav-header">Getting Started</li>
|
|
<l class="nav-item"i><a class="nav-link" href="/">Getting started</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Constants.html">Constants</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Matrix.html">Matrices and Vectors</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Euclid.html">Euclid & Number Theory</a></li>
|
|
<li class="nav-item">Combinatorics</li>
|
|
|
|
<li class="nav-header">Evaluation</li>
|
|
<li class="nav-item"><a class="nav-link" href="Functions.html">Special Functions</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Integration.html">Integration</a></li>
|
|
|
|
<li class="nav-header">Statistics/Probability</li>
|
|
<li class="nav-item"><a class="nav-link" href="DescriptiveStatistics.html">Descriptive Statistics</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Probability.html">Probability Distributions</a></li>
|
|
|
|
<li class="nav-header">Generation</li>
|
|
<li class="nav-item"><a class="nav-link" href="Generate.html">Generating Data</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="Random.html">Random Numbers</a></li>
|
|
|
|
<li class="nav-header">Solving Equations</li>
|
|
<li class="nav-item"><a class="nav-link" href="LinearEquations.html">Linear Equation Systems</a></li>
|
|
|
|
<li class="nav-header">Optimization</li>
|
|
<li class="nav-item"><a class="nav-link" href="Distance.html">Distance Metrics</a></li>
|
|
|
|
<li class="nav-header">Curve Fitting</li>
|
|
<li class="nav-item"><a class="nav-link" href="Regression.html">Regression</a></li>
|
|
|
|
<li class="nav-header">Native Providers</li>
|
|
<li class="nav-item"><a class="nav-link" href="MKL.html">Intel MKL</a></li>
|
|
|
|
<li class="nav-header">Working Together</li>
|
|
<li class="nav-item"><a class="nav-link" href="CSV.html">Delimited Text Files (CSV)</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="MatrixMarket.html">NIST MatrixMarket</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="MatlabFiles.html">MATLAB</a></li>
|
|
<li class="nav-item"><a class="nav-link" href="IFSharpNotebook.html">IF# Notebook</a></li>
|
|
</ul>
|
|
</div>
|
|
</nav>
|
|
<div class="container">
|
|
<div class="masthead">
|
|
<h3 class="muted">
|
|
<a href="https://numerics.mathdotnet.com">Math.NET Numerics</a> |
|
|
<a href="https://www.mathdotnet.com">Math.NET Project</a> |
|
|
<a href="https://github.com/mathnet/mathnet-numerics">GitHub</a>
|
|
</h3>
|
|
</div>
|
|
<hr />
|
|
<div class="container" id="fsdocs-content">
|
|
<h1><a name="Numerical-Integration" class="anchor" href="#Numerical-Integration">Numerical Integration</a></h1>
|
|
<p>The following double precision numerical integration or quadrature rules are supported in Math.NET Numerics under the <code>MathNet.Numerics.Integration</code> namespace. Unless stated otherwise, the examples below evaluate the integral <span class="math">\(\int_0^{10} x^2 \, dx = \frac{1000}{3} \approx 333.\overline{3}\)</span>.</p>
|
|
<h2><a name="Simpson-s-Rule" class="anchor" href="#Simpson-s-Rule">Simpson's Rule</a></h2>
|
|
<table class="pre"><tr><td class="snippet"><pre class="fssnip highlighted"><code lang="csharp"><span class="c">// Composite approximation with 4 partitions</span>
|
|
<span class="k">double</span> composite <span class="o">=</span> SimpsonRule.IntegrateComposite(x <span class="o">=</span><span class="o">></span> x <span class="o">*</span> x, <span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">4</span>);
|
|
|
|
<span class="c">// Approximate value using IntegrateComposite with 4 partitions is: 333.33333333333337</span>
|
|
Console.WriteLine(<span class="s">"Approximate value using IntegrateComposite with 4 partitions is: "</span> <span class="o">+</span> composite);
|
|
|
|
<span class="c">// Three point approximation</span>
|
|
<span class="k">double</span> threePoint <span class="o">=</span> SimpsonRule.IntegrateThreePoint(x <span class="o">=</span><span class="o">></span> x <span class="o">*</span> x, <span class="n">0.0</span>, <span class="n">10.0</span>);
|
|
|
|
<span class="c">// Approximate value using IntegrateThreePoint is: 333.333333333333</span>
|
|
Console.WriteLine(<span class="s">"Approximate value using IntegrateThreePoint is: "</span> <span class="o">+</span> threePoint);
|
|
</code></pre></td></tr></table>
|
|
<h2><a name="Newton-Cotes-Trapezium-Rule" class="anchor" href="#Newton-Cotes-Trapezium-Rule">Newton Cotes Trapezium Rule</a></h2>
|
|
<table class="pre"><tr><td class="snippet"><pre class="fssnip highlighted"><code lang="csharp"><span class="c">// Adaptive approximation with a relative error of 1e-5</span>
|
|
<span class="k">double</span> adaptive <span class="o">=</span> NewtonCotesTrapeziumRule.IntegrateAdaptive(x <span class="o">=</span><span class="o">></span> x <span class="o">*</span> x, <span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">1</span>e<span class="n">-5</span>);
|
|
|
|
<span class="c">// Approximate value of the integral using IntegrateAdaptive with a relative error of 1e-5 is: 333.333969116211</span>
|
|
Console.WriteLine(<span class="s">"Approximate value using IntegrateAdaptive with a relative error of 1e-5: "</span> <span class="o">+</span> adaptive);
|
|
|
|
<span class="c">// Composite approximation with 15 partitions</span>
|
|
<span class="k">double</span> composite <span class="o">=</span> NewtonCotesTrapeziumRule.IntegrateComposite(x <span class="o">=</span><span class="o">></span> x <span class="o">*</span> x, <span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">15</span>);
|
|
|
|
<span class="c">//Approximate value of the integral using IntegrateComposite with 15 partitions is: 334.074074074074</span>
|
|
Console.WriteLine(<span class="s">"Approximate value using IntegrateComposite with 15 partitions is: "</span> <span class="o">+</span> composite);
|
|
|
|
<span class="c">// Two point approximation</span>
|
|
<span class="k">double</span> twoPoint <span class="o">=</span> NewtonCotesTrapeziumRule.IntegrateTwoPoint(x <span class="o">=</span><span class="o">></span> x <span class="o">*</span> x, <span class="n">0.0</span>, <span class="n">10.0</span>);
|
|
|
|
<span class="c">//Approximate value using IntegrateTwoPoint is: 500</span>
|
|
Console.WriteLine(<span class="s">"Approximate value using IntegrateTwoPoint is: "</span> <span class="o">+</span> twoPoint);
|
|
</code></pre></td></tr></table>
|
|
<h2><a name="Double-Exponential-Transformation" class="anchor" href="#Double-Exponential-Transformation">Double-Exponential Transformation</a></h2>
|
|
<p>The Double-Exponential Transformation is suited for integration of smooth functions with no discontinuities, derivative discontinuities, and poles inside the interval.</p>
|
|
<table class="pre"><tr><td class="snippet"><pre class="fssnip highlighted"><code lang="csharp"><span class="c">// Approximate using a relative error of 1e-5.</span>
|
|
<span class="k">double</span> integrate <span class="o">=</span> DoubleExponentialTransformation.Integrate(x <span class="o">=</span><span class="o">></span> x <span class="o">*</span> x, <span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">1</span>e<span class="n">-5</span>);
|
|
|
|
<span class="c">// Approximate value using a relative error of 1e-5 is: 333.333333333332</span>
|
|
Console.WriteLine(<span class="s">"Approximate value using a relative error of 1e-5 is: "</span> <span class="o">+</span> integrate);
|
|
</code></pre></td></tr></table>
|
|
<h2><a name="Gauss-Legendre-Rule" class="anchor" href="#Gauss-Legendre-Rule">Gauss-Legendre Rule</a></h2>
|
|
<p>A fixed-order Gauss-Legendre integration routine is provided for fast integration of smooth functions with known polynomial order. The N-point Gauss-Legendre rule is exact for polynomials of order <span class="math">\(2N-1\)</span> or less. For example, these rules are useful when integrating basis functions to form mass matrices for the Galerkin method <a href="https://www.gnu.org/software/gsl/">[GSL]</a>.</p>
|
|
<p>The basic idea of Gauss-Legendre integration is to approximate the integral of a function <span class="math">\(f(x)\)</span> using <span class="math">\(N\)</span> Weights <span class="math">\(w_i\)</span> and abscissas (or nodes) <span class="math">\(x_i\)</span>.</p>
|
|
<p><span class="math">\[\int_a^b f(x) \, dx \approx \sum_{i = 0}^{N - 1} w_i f(x_i)\]</span></p>
|
|
<p>This algorithm calculates the abscissas and weights for a given order and integration interval. For efficiency, pre-computed abscissas and weights for the orders <span class="math">\(N = 2 - 20, \, 32, \, 64, \, 96, 100, \, 128, \, 256, \, 512, \, 1024\)</span> are used. Otherwise, they are calculated on the fly using Newton's method. For more information on the algorithm see <a href="https://www.holoborodko.com/pavel/numerical-methods/numerical-integration/">[Holoborodko, Pavel] </a>.</p>
|
|
<h3><a name="Abscissas-and-Weights" class="anchor" href="#Abscissas-and-Weights">Abscissas and Weights</a></h3>
|
|
<p>We'll first use the abscissas and weights to approximate an integral using a 5-point Gauss-Legendre rule</p>
|
|
<table class="pre"><tr><td class="snippet"><pre class="fssnip highlighted"><code lang="csharp"><span class="c">// Create a 5-point Gauss-Legendre rule over the integration interval [0, 10]</span>
|
|
GaussLegendreRule rule <span class="o">=</span> <span class="k">new</span> GaussLegendreRule(<span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">5</span>);
|
|
|
|
<span class="k">double</span> sum <span class="o">=</span> <span class="n">0</span>; <span class="c">// Will hold the approximate value of the integral</span>
|
|
<span class="k">for</span> (<span class="k">int</span> i <span class="o">=</span> <span class="n">0</span>; i <span class="o"><</span> rule.Order; i+<span class="o">+</span>) <span class="c">// rule.Order = 5</span>
|
|
{
|
|
<span class="c">// Access the ith abscissa and weight</span>
|
|
sum <span class="o">+</span><span class="o">=</span> rule.GetWeight(i) <span class="o">*</span> rule.GetAbscissa(i) <span class="o">*</span> rule.GetAbscissa(i);
|
|
}
|
|
|
|
<span class="c">// Approximate value is: 333.333333333333</span>
|
|
Console.WriteLine(<span class="s">"Approximate value is: "</span> <span class="o">+</span> sum);
|
|
</code></pre></td></tr></table>
|
|
<p>If you prefer direct access to the abscissas and weights, as opposed to using the methods</p>
|
|
<ul>
|
|
<li><code>double GetAbscissa(int i)</code></li>
|
|
<li><code>double GetWeight(int i)</code></li>
|
|
</ul>
|
|
<p>then use the properties <code>Abscissas</code> and <code>Weights</code></p>
|
|
<table class="pre"><tr><td class="snippet"><pre class="fssnip highlighted"><code lang="csharp"><span class="c">// Create a 5-point Gauss-Legendre rule over the integration interval [0, 10]</span>
|
|
GaussLegendreRule rule <span class="o">=</span> <span class="k">new</span> GaussLegendreRule(<span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">5</span>);
|
|
|
|
<span class="k">double</span>[] x <span class="o">=</span> rule.Abscissas; <span class="c">// Creates a clone and returns array of abscissas</span>
|
|
<span class="k">double</span>[] w <span class="o">=</span> rule.Weights; <span class="c">// Creates a clone and returns array of weights</span>
|
|
|
|
<span class="k">double</span> sum <span class="o">=</span> <span class="n">0</span>; <span class="c">// Will hold the approximate value of the integral</span>
|
|
<span class="k">for</span> (<span class="k">int</span> i <span class="o">=</span> <span class="n">0</span>; i <span class="o"><</span> rule.Order; i+<span class="o">+</span>) <span class="c">// rule.Order = 5</span>
|
|
{
|
|
<span class="c">// Access the ith abscissa and weight</span>
|
|
sum <span class="o">+</span><span class="o">=</span> w[i] <span class="o">*</span> x[i] <span class="o">*</span> x[i];
|
|
}
|
|
|
|
<span class="c">// Approximate value is: 333.333333333333</span>
|
|
Console.WriteLine(<span class="s">"Approximate value is: "</span> <span class="o">+</span> sum);;
|
|
</code></pre></td></tr></table>
|
|
<p>In addition to obtaining the abscissas and weights, the order and integration interval can be obtained</p>
|
|
<table class="pre"><tr><td class="snippet"><pre class="fssnip highlighted"><code lang="csharp"><span class="c">// Create a 5-point Gauss-Legendre rule over the integration interval [0, 10]</span>
|
|
GaussLegendreRule rule <span class="o">=</span> <span class="k">new</span> GaussLegendreRule(<span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">5</span>);
|
|
|
|
<span class="c">// The order of the rule is: 5</span>
|
|
Console.WriteLine(<span class="s">"The order of the rule is: "</span> <span class="o">+</span> rule.Order);
|
|
|
|
<span class="c">// The lower integral bound is 0</span>
|
|
Console.WriteLine(<span class="s">"The lower integral bound is: "</span> <span class="o">+</span> rule.IntervalBegin);
|
|
|
|
<span class="c">// The upper integral bound is 10</span>
|
|
Console.WriteLine(<span class="s">"The upper integral bound is: "</span> <span class="o">+</span> rule.IntervalEnd);
|
|
</code></pre></td></tr></table>
|
|
<h3><a name="Integrate-Method" class="anchor" href="#Integrate-Method">Integrate Method</a></h3>
|
|
<p>For convenience, we provide an overloaded static method <code>double Integrate(...)</code> which preforms 1D and 2D integration of a function. The first parameter to the method is a delegate of type <code>Func<double, double></code> or <code>Func<double, double, double></code> for 1D and 2D integration respectively. So for example</p>
|
|
<table class="pre"><tr><td class="snippet"><pre class="fssnip highlighted"><code lang="csharp"><span class="c">// 1D integration using a 5-point Gauss-Legendre rule over the integration interval [0, 10]</span>
|
|
<span class="k">double</span> integrate<span class="n">1</span>D <span class="o">=</span> GaussLegendreRule.Integrate(x <span class="o">=</span><span class="o">></span> x <span class="o">*</span> x, <span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">5</span>);
|
|
|
|
<span class="c">// Approximate value of the 1D integral is: 333.333333333333</span>
|
|
Console.WriteLine(<span class="s">"Approximate value of the 1D integral is: "</span> <span class="o">+</span> integrate<span class="n">1</span>D);
|
|
|
|
<span class="c">// 2D integration using a 5-point Gauss-Legendre rule over the integration interval [0, 10] X [1, 2]</span>
|
|
<span class="k">double</span> integrate<span class="n">2</span>D <span class="o">=</span> GaussLegendreRule.Integrate((x, y) <span class="o">=</span><span class="o">></span> (x <span class="o">*</span> x) <span class="o">*</span> (y <span class="o">*</span> y), <span class="n">0.0</span>, <span class="n">10.0</span>, <span class="n">1.0</span>, <span class="n">2.0</span>, <span class="n">5</span>);
|
|
|
|
<span class="c">// Approximate value of the 2D integral is: 777.777777777778</span>
|
|
Console.WriteLine(<span class="s">"Approximate value of the 2D integral is: "</span> <span class="o">+</span> integrate<span class="n">2</span>D);
|
|
</code></pre></td></tr></table>
|
|
<p>where we used <span class="math">\(\int_0^{10}\int_1^2 x^2 y^2 \,dydx = \frac{7000}{9} \approx 777.\overline{7}\)</span> for the 2D integral example.</p>
|
|
|
|
|
|
</div>
|
|
<!-- BEGIN SEARCH BOX: this adds support for the search box -->
|
|
<link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/JavaScript-autoComplete/1.0.4/auto-complete.css" />
|
|
<script type="text/javascript">var fsdocs_search_baseurl = 'https://numerics.mathdotnet.com/';</script>
|
|
<script type="text/javascript" src="https://cdnjs.cloudflare.com/ajax/libs/lunr.js/2.3.8/lunr.min.js"></script>
|
|
<script type="text/javascript" src="https://cdnjs.cloudflare.com/ajax/libs/JavaScript-autoComplete/1.0.4/auto-complete.min.js"></script>
|
|
<script type="text/javascript" src="https://numerics.mathdotnet.com/content/fsdocs-search.js"></script>
|
|
<!-- END SEARCH BOX: this adds support for the search box -->
|
|
</div>
|
|
</body>
|
|
|
|
</html>
|
|
|