<html lang="en"> <head> <title>Functions of One Variable - GNU Octave</title> <meta http-equiv="Content-Type" content="text/html"> <meta name="description" content="GNU Octave"> <meta name="generator" content="makeinfo 4.13"> <link title="Top" rel="start" href="index.html#Top"> <link rel="up" href="Numerical-Integration.html#Numerical-Integration" title="Numerical Integration"> <link rel="next" href="Orthogonal-Collocation.html#Orthogonal-Collocation" title="Orthogonal Collocation"> <link href="http://www.gnu.org/software/texinfo/" rel="generator-home" title="Texinfo Homepage"> <meta http-equiv="Content-Style-Type" content="text/css"> <style type="text/css"><!-- pre.display { font-family:inherit } pre.format { font-family:inherit } pre.smalldisplay { font-family:inherit; font-size:smaller } pre.smallformat { font-family:inherit; font-size:smaller } pre.smallexample { font-size:smaller } pre.smalllisp { font-size:smaller } span.sc { font-variant:small-caps } span.roman { font-family:serif; font-weight:normal; } span.sansserif { font-family:sans-serif; font-weight:normal; } --></style> </head> <body> <div class="node"> <a name="Functions-of-One-Variable"></a> <p> Next: <a rel="next" accesskey="n" href="Orthogonal-Collocation.html#Orthogonal-Collocation">Orthogonal Collocation</a>, Up: <a rel="up" accesskey="u" href="Numerical-Integration.html#Numerical-Integration">Numerical Integration</a> <hr> </div> <h3 class="section">23.1 Functions of One Variable</h3> <p>Octave supports five different algorithms for computing the integral of a function f over the interval from a to b. These are <dl> <dt><code>quad</code><dd>Numerical integration based on Gaussian quadrature. <br><dt><code>quadv</code><dd>Numerical integration using an adaptive vectorized Simpson's rule. <br><dt><code>quadl</code><dd>Numerical integration using an adaptive Lobatto rule. <br><dt><code>quadgk</code><dd>Numerical integration using an adaptive Gauss-Konrod rule. <br><dt><code>quadcc</code><dd>Numerical integration using adaptive Clenshaw-Curtis rules. <br><dt><code>trapz, cumtrapz</code><dd>Numerical integration of data using the trapezoidal method. </dl> <p class="noindent">The best quadrature algorithm to use depends on the integrand. If you have empirical data, rather than a function, the choice is <code>trapz</code> or <code>cumtrapz</code>. If you are uncertain about the characteristics of the integrand, <code>quadcc</code> will be the most robust as it can handle discontinuities, singularities, oscillatory functions, and infinite intervals. When the integrand is smooth <code>quadgk</code> may be the fastest of the algorithms. <p><table summary=""><tr align="left"><th valign="top" width="5%"></th><th valign="top" width="15%">Function </th><th valign="top" width="80%">Characteristics <br></th></tr><tr align="left"><td valign="top" width="5%"></td><td valign="top" width="15%">quad </td><td valign="top" width="80%">Low accuracy with nonsmooth integrands <br></td></tr><tr align="left"><td valign="top" width="5%"></td><td valign="top" width="15%">quadv </td><td valign="top" width="80%">Medium accuracy with smooth integrands <br></td></tr><tr align="left"><td valign="top" width="5%"></td><td valign="top" width="15%">quadl </td><td valign="top" width="80%">Medium accuracy with smooth integrands. Slower than quadgk. <br></td></tr><tr align="left"><td valign="top" width="5%"></td><td valign="top" width="15%">quadgk </td><td valign="top" width="80%">Medium accuracy (1e^-6–1e^-9) with smooth integrands. <br></td></tr><tr align="left"><td valign="top" width="5%"></td><td valign="top" width="15%"></td><td valign="top" width="80%">Handles oscillatory functions and infinite bounds <br></td></tr><tr align="left"><td valign="top" width="5%"></td><td valign="top" width="15%">quadcc </td><td valign="top" width="80%">Low to High accuracy with nonsmooth/smooth integrands <br></td></tr><tr align="left"><td valign="top" width="5%"></td><td valign="top" width="15%"></td><td valign="top" width="80%">Handles oscillatory functions, singularities, and infinite bounds <br></td></tr></table> <p>Here is an example of using <code>quad</code> to integrate the function <pre class="example"> <var>f</var>(<var>x</var>) = <var>x</var> * sin (1/<var>x</var>) * sqrt (abs (1 - <var>x</var>)) </pre> <p class="noindent">from <var>x</var> = 0 to <var>x</var> = 3. <p>This is a fairly difficult integration (plot the function over the range of integration to see why). <p>The first step is to define the function: <pre class="example"> function y = f (x) y = x .* sin (1./x) .* sqrt (abs (1 - x)); endfunction </pre> <p>Note the use of the `dot' forms of the operators. This is not necessary for the <code>quad</code> integrator, but is required by the other integrators. In any case, it makes it much easier to generate a set of points for plotting because it is possible to call the function with a vector argument to produce a vector result. <p>The second step is to call quad with the limits of integration: <pre class="example"> [q, ier, nfun, err] = quad ("f", 0, 3) ⇒ 1.9819 ⇒ 1 ⇒ 5061 ⇒ 1.1522e-07 </pre> <p>Although <code>quad</code> returns a nonzero value for <var>ier</var>, the result is reasonably accurate (to see why, examine what happens to the result if you move the lower bound to 0.1, then 0.01, then 0.001, etc.). <p>The function "f" can be the string name of a function, a function handle, or an inline function. These options make it quite easy to do integration without having to fully define a function in an m-file. For example: <pre class="example"> # Verify integral (x^3) = x^4/4 f = inline ("x.^3"); quadgk (f, 0, 1) ⇒ 0.25000 # Verify gamma function = (n-1)! for n = 4 f = @(x) x.^3 .* exp (-x); quadcc (f, 0, Inf) ⇒ 6.0000 </pre> <!-- quad src/DLD-FUNCTIONS/quad.cc --> <p><a name="doc_002dquad"></a> <div class="defun"> — Loadable Function: <var>q</var> = <b>quad</b> (<var>f, a, b</var>)<var><a name="index-quad-2303"></a></var><br> — Loadable Function: <var>q</var> = <b>quad</b> (<var>f, a, b, tol</var>)<var><a name="index-quad-2304"></a></var><br> — Loadable Function: <var>q</var> = <b>quad</b> (<var>f, a, b, tol, sing</var>)<var><a name="index-quad-2305"></a></var><br> — Loadable Function: [<var>q</var>, <var>ier</var>, <var>nfun</var>, <var>err</var>] = <b>quad</b> (<var><small class="dots">...</small></var>)<var><a name="index-quad-2306"></a></var><br> <blockquote><p>Numerically evaluate the integral of <var>f</var> from <var>a</var> to <var>b</var> using Fortran routines from <span class="sc">quadpack</span><!-- /@w -->. <var>f</var> is a function handle, inline function, or a string containing the name of the function to evaluate. The function must have the form <code>y = f (x)</code> where <var>y</var> and <var>x</var> are scalars. <p><var>a</var> and <var>b</var> are the lower and upper limits of integration. Either or both may be infinite. <p>The optional argument <var>tol</var> is a vector that specifies the desired accuracy of the result. The first element of the vector is the desired absolute tolerance, and the second element is the desired relative tolerance. To choose a relative test only, set the absolute tolerance to zero. To choose an absolute test only, set the relative tolerance to zero. Both tolerances default to <code>sqrt(eps)</code> or approximately 1.5e^-8. <p>The optional argument <var>sing</var> is a vector of values at which the integrand is known to be singular. <p>The result of the integration is returned in <var>q</var>. <var>ier</var> contains an integer error code (0 indicates a successful integration). <var>nfun</var> indicates the number of function evaluations that were made, and <var>err</var> contains an estimate of the error in the solution. <p>The function <code>quad_options</code> can set other optional parameters for <code>quad</code>. <p>Note: because <code>quad</code> is written in Fortran it cannot be called recursively. This prevents its use in integrating over more than one variable by routines <code>dblquad</code> and <code>triplequad</code>. <!-- Texinfo @sp should work but in practice produces ugly results for HTML. --> <!-- A simple blank line produces the correct behavior. --> <!-- @sp 1 --> <p class="noindent"><strong>See also:</strong> <a href="doc_002dquad_005foptions.html#doc_002dquad_005foptions">quad_options</a>, <a href="doc_002dquadv.html#doc_002dquadv">quadv</a>, <a href="doc_002dquadl.html#doc_002dquadl">quadl</a>, <a href="doc_002dquadgk.html#doc_002dquadgk">quadgk</a>, <a href="doc_002dquadcc.html#doc_002dquadcc">quadcc</a>, <a href="doc_002dtrapz.html#doc_002dtrapz">trapz</a>, <a href="doc_002ddblquad.html#doc_002ddblquad">dblquad</a>, <a href="doc_002dtriplequad.html#doc_002dtriplequad">triplequad</a>. </p></blockquote></div> <!-- quad_options src/DLD-FUNCTIONS/quad.cc --> <p><a name="doc_002dquad_005foptions"></a> <div class="defun"> — Loadable Function: <b>quad_options</b> ()<var><a name="index-quad_005foptions-2307"></a></var><br> — Loadable Function: val = <b>quad_options</b> (<var>opt</var>)<var><a name="index-quad_005foptions-2308"></a></var><br> — Loadable Function: <b>quad_options</b> (<var>opt, val</var>)<var><a name="index-quad_005foptions-2309"></a></var><br> <blockquote><p>Query or set options for the function <code>quad</code>. When called with no arguments, the names of all available options and their current values are displayed. Given one argument, return the value of the corresponding option. When called with two arguments, <code>quad_options</code> set the option <var>opt</var> to value <var>val</var>. <p>Options include <dl> <dt><code>"absolute tolerance"</code><dd>Absolute tolerance; may be zero for pure relative error test. <br><dt><code>"relative tolerance"</code><dd>Non-negative relative tolerance. If the absolute tolerance is zero, the relative tolerance must be greater than or equal to <code>max (50*eps, 0.5e-28)</code><!-- /@w -->. <br><dt><code>"single precision absolute tolerance"</code><dd>Absolute tolerance for single precision; may be zero for pure relative error test. <br><dt><code>"single precision relative tolerance"</code><dd>Non-negative relative tolerance for single precision. If the absolute tolerance is zero, the relative tolerance must be greater than or equal to <code>max (50*eps, 0.5e-28)</code><!-- /@w -->. </dl> </p></blockquote></div> <!-- quadv scripts/general/quadv.m --> <p><a name="doc_002dquadv"></a> <div class="defun"> — Function File: <var>q</var> = <b>quadv</b> (<var>f, a, b</var>)<var><a name="index-quadv-2310"></a></var><br> — Function File: <var>q</var> = <b>quadv</b> (<var>f, a, b, tol</var>)<var><a name="index-quadv-2311"></a></var><br> — Function File: <var>q</var> = <b>quadv</b> (<var>f, a, b, tol, trace</var>)<var><a name="index-quadv-2312"></a></var><br> — Function File: <var>q</var> = <b>quadv</b> (<var>f, a, b, tol, trace, p1, p2, <small class="dots">...</small></var>)<var><a name="index-quadv-2313"></a></var><br> — Function File: [<var>q</var>, <var>nfun</var>] = <b>quadv</b> (<var><small class="dots">...</small></var>)<var><a name="index-quadv-2314"></a></var><br> <blockquote> <p>Numerically evaluate the integral of <var>f</var> from <var>a</var> to <var>b</var> using an adaptive Simpson's rule. <var>f</var> is a function handle, inline function, or string containing the name of the function to evaluate. <code>quadv</code> is a vectorized version of <code>quad</code> and the function defined by <var>f</var> must accept a scalar or vector as input and return a scalar, vector, or array as output. <p><var>a</var> and <var>b</var> are the lower and upper limits of integration. Both limits must be finite. <p>The optional argument <var>tol</var> defines the tolerance used to stop the adaptation procedure. The default value is 1e^-6. <p>The algorithm used by <code>quadv</code> involves recursively subdividing the integration interval and applying Simpson's rule on each subinterval. If <var>trace</var> is true then after computing each of these partial integrals display: (1) the total number of function evaluations, (2) the left end of the subinterval, (3) the length of the subinterval, (4) the approximation of the integral over the subinterval. <p>Additional arguments <var>p1</var>, etc., are passed directly to the function <var>f</var>. To use default values for <var>tol</var> and <var>trace</var>, one may pass empty matrices ([]). <p>The result of the integration is returned in <var>q</var>. <var>nfun</var> indicates the number of function evaluations that were made. <p>Note: <code>quadv</code> is written in Octave's scripting language and can be used recursively in <code>dblquad</code> and <code>triplequad</code>, unlike the similar <code>quad</code> function. <!-- Texinfo @sp should work but in practice produces ugly results for HTML. --> <!-- A simple blank line produces the correct behavior. --> <!-- @sp 1 --> <p class="noindent"><strong>See also:</strong> <a href="doc_002dquad.html#doc_002dquad">quad</a>, <a href="doc_002dquadl.html#doc_002dquadl">quadl</a>, <a href="doc_002dquadgk.html#doc_002dquadgk">quadgk</a>, <a href="doc_002dquadcc.html#doc_002dquadcc">quadcc</a>, <a href="doc_002dtrapz.html#doc_002dtrapz">trapz</a>, <a href="doc_002ddblquad.html#doc_002ddblquad">dblquad</a>, <a href="doc_002dtriplequad.html#doc_002dtriplequad">triplequad</a>. </p></blockquote></div> <!-- quadl scripts/general/quadl.m --> <p><a name="doc_002dquadl"></a> <div class="defun"> — Function File: <var>q</var> = <b>quadl</b> (<var>f, a, b</var>)<var><a name="index-quadl-2315"></a></var><br> — Function File: <var>q</var> = <b>quadl</b> (<var>f, a, b, tol</var>)<var><a name="index-quadl-2316"></a></var><br> — Function File: <var>q</var> = <b>quadl</b> (<var>f, a, b, tol, trace</var>)<var><a name="index-quadl-2317"></a></var><br> — Function File: <var>q</var> = <b>quadl</b> (<var>f, a, b, tol, trace, p1, p2, <small class="dots">...</small></var>)<var><a name="index-quadl-2318"></a></var><br> <blockquote> <p>Numerically evaluate the integral of <var>f</var> from <var>a</var> to <var>b</var> using an adaptive Lobatto rule. <var>f</var> is a function handle, inline function, or string containing the name of the function to evaluate. The function <var>f</var> must be vectorized and return a vector of output values if given a vector of input values. <p><var>a</var> and <var>b</var> are the lower and upper limits of integration. Both limits must be finite. <p>The optional argument <var>tol</var> defines the relative tolerance with which to perform the integration. The default value is <code>eps</code>. <p>The algorithm used by <code>quadl</code> involves recursively subdividing the integration interval. If <var>trace</var> is defined then for each subinterval display: (1) the left end of the subinterval, (2) the length of the subinterval, (3) the approximation of the integral over the subinterval. <p>Additional arguments <var>p1</var>, etc., are passed directly to the function <var>f</var>. To use default values for <var>tol</var> and <var>trace</var>, one may pass empty matrices ([]). <p>Reference: W. Gander and W. Gautschi, <cite>Adaptive Quadrature - Revisited</cite>, BIT Vol. 40, No. 1, March 2000, pp. 84–101. <a href="http://www.inf.ethz.ch/personal/gander/">http://www.inf.ethz.ch/personal/gander/</a> <!-- Texinfo @sp should work but in practice produces ugly results for HTML. --> <!-- A simple blank line produces the correct behavior. --> <!-- @sp 1 --> <p class="noindent"><strong>See also:</strong> <a href="doc_002dquad.html#doc_002dquad">quad</a>, <a href="doc_002dquadv.html#doc_002dquadv">quadv</a>, <a href="doc_002dquadgk.html#doc_002dquadgk">quadgk</a>, <a href="doc_002dquadcc.html#doc_002dquadcc">quadcc</a>, <a href="doc_002dtrapz.html#doc_002dtrapz">trapz</a>, <a href="doc_002ddblquad.html#doc_002ddblquad">dblquad</a>, <a href="doc_002dtriplequad.html#doc_002dtriplequad">triplequad</a>. </p></blockquote></div> <!-- quadgk scripts/general/quadgk.m --> <p><a name="doc_002dquadgk"></a> <div class="defun"> — Function File: <var>q</var> = <b>quadgk</b> (<var>f, a, b</var>)<var><a name="index-quadgk-2319"></a></var><br> — Function File: <var>q</var> = <b>quadgk</b> (<var>f, a, b, abstol</var>)<var><a name="index-quadgk-2320"></a></var><br> — Function File: <var>q</var> = <b>quadgk</b> (<var>f, a, b, abstol, trace</var>)<var><a name="index-quadgk-2321"></a></var><br> — Function File: <var>q</var> = <b>quadgk</b> (<var>f, a, b, prop, val, <small class="dots">...</small></var>)<var><a name="index-quadgk-2322"></a></var><br> — Function File: [<var>q</var>, <var>err</var>] = <b>quadgk</b> (<var><small class="dots">...</small></var>)<var><a name="index-quadgk-2323"></a></var><br> <blockquote> <p>Numerically evaluate the integral of <var>f</var> from <var>a</var> to <var>b</var> using adaptive Gauss-Konrod quadrature. <var>f</var> is a function handle, inline function, or string containing the name of the function to evaluate. The formulation is based on a proposal by L.F. Shampine, <cite>"Vectorized adaptive quadrature in </cite><span class="sc">matlab</span><cite>", Journal of Computational and Applied Mathematics, pp131-140, Vol 211, Issue 2, Feb 2008</cite> where all function evaluations at an iteration are calculated with a single call to <var>f</var>. Therefore, the function <var>f</var> must be vectorized and must accept a vector of input values <var>x</var> and return an output vector representing the function evaluations at the given values of <var>x</var>. <p><var>a</var> and <var>b</var> are the lower and upper limits of integration. Either or both limits may be infinite or contain weak end singularities. Variable transformation will be used to treat any infinite intervals and weaken the singularities. For example: <pre class="example"> quadgk (@(x) 1 ./ (sqrt (x) .* (x + 1)), 0, Inf) </pre> <p class="noindent">Note that the formulation of the integrand uses the element-by-element operator <code>./</code> and all user functions to <code>quadgk</code> should do the same. <p>The optional argument <var>tol</var> defines the absolute tolerance used to stop the integration procedure. The default value is 1e^-10. <p>The algorithm used by <code>quadgk</code> involves subdividing the integration interval and evaluating each subinterval. If <var>trace</var> is true then after computing each of these partial integrals display: (1) the number of subintervals at this step, (2) the current estimate of the error <var>err</var>, (3) the current estimate for the integral <var>q</var>. <p>Alternatively, properties of <code>quadgk</code> can be passed to the function as pairs <code>"</code><var>prop</var><code>", </code><var>val</var>. Valid properties are <dl> <dt><code>AbsTol</code><dd>Define the absolute error tolerance for the quadrature. The default absolute tolerance is 1e-10. <br><dt><code>RelTol</code><dd>Define the relative error tolerance for the quadrature. The default relative tolerance is 1e-5. <br><dt><code>MaxIntervalCount</code><dd><code>quadgk</code> initially subdivides the interval on which to perform the quadrature into 10 intervals. Subintervals that have an unacceptable error are subdivided and re-evaluated. If the number of subintervals exceeds 650 subintervals at any point then a poor convergence is signaled and the current estimate of the integral is returned. The property 'MaxIntervalCount' can be used to alter the number of subintervals that can exist before exiting. <br><dt><code>WayPoints</code><dd>Discontinuities in the first derivative of the function to integrate can be flagged with the <code>"WayPoints"</code> property. This forces the ends of a subinterval to fall on the breakpoints of the function and can result in significantly improved estimation of the error in the integral, faster computation, or both. For example, <pre class="example"> quadgk (@(x) abs (1 - x.^2), 0, 2, "Waypoints", 1) </pre> <p class="noindent">signals the breakpoint in the integrand at <var>x</var><code> = 1</code>. <br><dt><code>Trace</code><dd>If logically true <code>quadgk</code> prints information on the convergence of the quadrature at each iteration. </dl> <p>If any of <var>a</var>, <var>b</var>, or <var>waypoints</var> is complex then the quadrature is treated as a contour integral along a piecewise continuous path defined by the above. In this case the integral is assumed to have no edge singularities. For example, <pre class="example"> quadgk (@(z) log (z), 1+1i, 1+1i, "WayPoints", [1-1i, -1,-1i, -1+1i]) </pre> <p class="noindent">integrates <code>log (z)</code> along the square defined by <code>[1+1i, 1-1i, -1-1i, -1+1i]</code> <p>The result of the integration is returned in <var>q</var>. <var>err</var> is an approximate bound on the error in the integral <code>abs (</code><var>q</var><code> - </code><var>I</var><code>)</code>, where <var>I</var> is the exact value of the integral. <!-- Texinfo @sp should work but in practice produces ugly results for HTML. --> <!-- A simple blank line produces the correct behavior. --> <!-- @sp 1 --> <p class="noindent"><strong>See also:</strong> <a href="doc_002dquad.html#doc_002dquad">quad</a>, <a href="doc_002dquadv.html#doc_002dquadv">quadv</a>, <a href="doc_002dquadl.html#doc_002dquadl">quadl</a>, <a href="doc_002dquadcc.html#doc_002dquadcc">quadcc</a>, <a href="doc_002dtrapz.html#doc_002dtrapz">trapz</a>, <a href="doc_002ddblquad.html#doc_002ddblquad">dblquad</a>, <a href="doc_002dtriplequad.html#doc_002dtriplequad">triplequad</a>. </p></blockquote></div> <!-- quadcc src/DLD-FUNCTIONS/quadcc.cc --> <p><a name="doc_002dquadcc"></a> <div class="defun"> — Function File: <var>q</var> = <b>quadcc</b> (<var>f, a, b</var>)<var><a name="index-quadcc-2324"></a></var><br> — Function File: <var>q</var> = <b>quadcc</b> (<var>f, a, b, tol</var>)<var><a name="index-quadcc-2325"></a></var><br> — Function File: <var>q</var> = <b>quadcc</b> (<var>f, a, b, tol, sing</var>)<var><a name="index-quadcc-2326"></a></var><br> — Function File: [<var>q</var>, <var>err</var>, <var>nr_points</var>] = <b>quadcc</b> (<var><small class="dots">...</small></var>)<var><a name="index-quadcc-2327"></a></var><br> <blockquote><p>Numerically evaluate the integral of <var>f</var> from <var>a</var> to <var>b</var> using the doubly-adaptive Clenshaw-Curtis quadrature described by P. Gonnet in <cite>Increasing the Reliability of Adaptive Quadrature Using Explicit Interpolants</cite>. <var>f</var> is a function handle, inline function, or string containing the name of the function to evaluate. The function <var>f</var> must be vectorized and must return a vector of output values if given a vector of input values. For example, <pre class="example"> f = @(x) x .* sin (1./x) .* sqrt (abs (1 - x)); </pre> <p class="noindent">which uses the element-by-element `dot' form for all operators. <p><var>a</var> and <var>b</var> are the lower and upper limits of integration. Either or both limits may be infinite. <code>quadcc</code> handles an inifinite limit by substituting the variable of integration with <code>x=tan(pi/2*u)</code>. <p>The optional argument <var>tol</var> defines the relative tolerance used to stop the integration procedure. The default value is 1e^-6. <p>The optional argument <var>sing</var> contains a list of points where the integrand has known singularities, or discontinuities in any of its derivatives, inside the integration interval. For the example above, which has a discontinuity at x=1, the call to <code>quadcc</code> would be as follows <pre class="example"> int = quadcc (f, a, b, 1.0e-6, [ 1 ]); </pre> <p>The result of the integration is returned in <var>q</var>. <var>err</var> is an estimate of the absolute integration error and <var>nr_points</var> is the number of points at which the integrand was evaluated. If the adaptive integration did not converge, the value of <var>err</var> will be larger than the requested tolerance. Therefore, it is recommended to verify this value for difficult integrands. <p><code>quadcc</code> is capable of dealing with non-numeric values of the integrand such as <code>NaN</code> or <code>Inf</code>. If the integral diverges, and <code>quadcc</code> detects this, then a warning is issued and <code>Inf</code> or <code>-Inf</code> is returned. <p>Note: <code>quadcc</code> is a general purpose quadrature algorithm and, as such, may be less efficient for a smooth or otherwise well-behaved integrand than other methods such as <code>quadgk</code>. <p>The algorithm uses Clenshaw-Curtis quadrature rules of increasing degree in each interval and bisects the interval if either the function does not appear to be smooth or a rule of maximum degree has been reached. The error estimate is computed from the L2-norm of the difference between two successive interpolations of the integrand over the nodes of the respective quadrature rules. <p>Reference: P. Gonnet, <cite>Increasing the Reliability of Adaptive Quadrature Using Explicit Interpolants</cite>, ACM Transactions on Mathematical Software, Vol. 37, Issue 3, Article No. 3, 2010. <!-- Texinfo @sp should work but in practice produces ugly results for HTML. --> <!-- A simple blank line produces the correct behavior. --> <!-- @sp 1 --> <p class="noindent"><strong>See also:</strong> <a href="doc_002dquad.html#doc_002dquad">quad</a>, <a href="doc_002dquadv.html#doc_002dquadv">quadv</a>, <a href="doc_002dquadl.html#doc_002dquadl">quadl</a>, <a href="doc_002dquadgk.html#doc_002dquadgk">quadgk</a>, <a href="doc_002dtrapz.html#doc_002dtrapz">trapz</a>, <a href="doc_002ddblquad.html#doc_002ddblquad">dblquad</a>, <a href="doc_002dtriplequad.html#doc_002dtriplequad">triplequad</a>. </p></blockquote></div> <p>Sometimes one does not have the function, but only the raw (x, y) points from which to perform an integration. This can occur when collecting data in an experiment. The <code>trapz</code> function can integrate these values as shown in the following example where "data" has been collected on the cosine function over the range [0, pi/2). <pre class="example"> x = 0:0.1:pi/2; # Uniformly spaced points y = cos (x); trapz (x, y) ⇒ 0.99666 </pre> <p>The answer is reasonably close to the exact value of 1. Ordinary quadrature is sensitive to the characteristics of the integrand. Empirical integration depends not just on the integrand, but also on the particular points chosen to represent the function. Repeating the example above with the sine function over the range [0, pi/2) produces far inferior results. <pre class="example"> x = 0:0.1:pi/2; # Uniformly spaced points y = sin (x); trapz (x, y) ⇒ 0.92849 </pre> <p>However, a slightly different choice of data points can change the result significantly. The same integration, with the same number of points, but spaced differently produces a more accurate answer. <pre class="example"> x = linspace (0, pi/2, 16); # Uniformly spaced, but including endpoint y = sin (x); trapz (x, y) ⇒ 0.99909 </pre> <p>In general there may be no way of knowing the best distribution of points ahead of time. Or the points may come from an experiment where there is no freedom to select the best distribution. In any case, one must remain aware of this issue when using <code>trapz</code>. <!-- trapz scripts/general/trapz.m --> <p><a name="doc_002dtrapz"></a> <div class="defun"> — Function File: <var>q</var> = <b>trapz</b> (<var>y</var>)<var><a name="index-trapz-2328"></a></var><br> — Function File: <var>q</var> = <b>trapz</b> (<var>x, y</var>)<var><a name="index-trapz-2329"></a></var><br> — Function File: <var>q</var> = <b>trapz</b> (<var><small class="dots">...</small>, dim</var>)<var><a name="index-trapz-2330"></a></var><br> <blockquote> <p>Numerically evaluate the integral of points <var>y</var> using the trapezoidal method. <code>trapz (</code><var>y</var><code>)</code><!-- /@w --> computes the integral of <var>y</var> along the first non-singleton dimension. When the argument <var>x</var> is omitted an equally spaced <var>x</var> vector with unit spacing (1) is assumed. <code>trapz (</code><var>x</var><code>, </code><var>y</var><code>)</code> evaluates the integral with respect to the spacing in <var>x</var> and the values in <var>y</var>. This is useful if the points in <var>y</var> have been sampled unevenly. If the optional <var>dim</var> argument is given, operate along this dimension. <p>If <var>x</var> is not specified then unit spacing will be used. To scale the integral to the correct value you must multiply by the actual spacing value (deltaX). As an example, the integral of x^3 over the range [0, 1] is x^4/4 or 0.25. The following code uses <code>trapz</code> to calculate the integral in three different ways. <pre class="example"> x = 0:0.1:1; y = x.^3; q = trapz (y) ⇒ q = 2.525 # No scaling q * 0.1 ⇒ q = 0.2525 # Approximation to integral by scaling trapz (x, y) ⇒ q = 0.2525 # Same result by specifying <var>x</var> </pre> <!-- Texinfo @sp should work but in practice produces ugly results for HTML. --> <!-- A simple blank line produces the correct behavior. --> <!-- @sp 1 --> <p class="noindent"><strong>See also:</strong> <a href="doc_002dcumtrapz.html#doc_002dcumtrapz">cumtrapz</a>. </p></blockquote></div> <!-- cumtrapz scripts/general/cumtrapz.m --> <p><a name="doc_002dcumtrapz"></a> <div class="defun"> — Function File: <var>q</var> = <b>cumtrapz</b> (<var>y</var>)<var><a name="index-cumtrapz-2331"></a></var><br> — Function File: <var>q</var> = <b>cumtrapz</b> (<var>x, y</var>)<var><a name="index-cumtrapz-2332"></a></var><br> — Function File: <var>q</var> = <b>cumtrapz</b> (<var><small class="dots">...</small>, dim</var>)<var><a name="index-cumtrapz-2333"></a></var><br> <blockquote> <p>Cumulative numerical integration of points <var>y</var> using the trapezoidal method. <code>cumtrapz (</code><var>y</var><code>)</code><!-- /@w --> computes the cumulative integral of <var>y</var> along the first non-singleton dimension. Where <code>trapz</code> reports only the overall integral sum, <code>cumtrapz</code> reports the current partial sum value at each point of <var>y</var>. When the argument <var>x</var> is omitted an equally spaced <var>x</var> vector with unit spacing (1) is assumed. <code>cumtrapz (</code><var>x</var><code>, </code><var>y</var><code>)</code> evaluates the integral with respect to the spacing in <var>x</var> and the values in <var>y</var>. This is useful if the points in <var>y</var> have been sampled unevenly. If the optional <var>dim</var> argument is given, operate along this dimension. <p>If <var>x</var> is not specified then unit spacing will be used. To scale the integral to the correct value you must multiply by the actual spacing value (deltaX). <!-- Texinfo @sp should work but in practice produces ugly results for HTML. --> <!-- A simple blank line produces the correct behavior. --> <!-- @sp 1 --> <p class="noindent"><strong>See also:</strong> <a href="doc_002dtrapz.html#doc_002dtrapz">trapz</a>, <a href="doc_002dcumsum.html#doc_002dcumsum">cumsum</a>. </p></blockquote></div> </body></html>