quadgk

Numerically evaluate integral — Gauss-Kronrod quadrature

Syntax

q = quadgk(fun,a,b)

[q,errbnd] = quadgk(fun,a,b)

[___] = quadgk(fun,a,b,Name,Value)

Description

q = quadgk(fun,a,b) integrates the function handle fun from a to b using high-order global adaptive quadrature and default error tolerances.

[q,errbnd] = quadgk(fun,a,b) additionally returns an approximate upper bound on the absolute error |q - I|, where I is the exact value of the integral.

[___] = quadgk(fun,a,b,Name,Value) specifies additional options with one or more name-value pair arguments using either of the previous output argument combinations. For example, specify 'Waypoints' followed by a vector of real or complex numbers to indicate specific points for the integrator to use.

Examples

Integral with Singularity at Endpoint

Evaluate the integral

$$\mathit{q}={\int }_0^1{\mathit{e}}^{\mathit{x}}\ln \left(\mathit{x}\right)\mathrm{dx}.$$

This integral has a singularity at the point $\mathit{x}=0$ because $\ln \left(0\right)$ diverges to $-\infty$.

Create an anonymous function for the integrand. The log function calculates $\ln \left(\mathit{x}\right)$.

f = @(x) exp(x).*log(x);

Integrate f from 0 to 1.

q = quadgk(f,0,1)

q = 
-1.3179

Complex Contour Integration

Integrate a complex function around a pole by specifying a contour.

Evaluate the complex contour integral

$$\mathit{q}=\oint \frac{\mathrm{dz}}{2\mathit{z}-1}.$$

The integrand has a simple pole at $\mathit{z}=1/2$, so use a rectangular contour that encloses that point. The contour starts and ends at $\mathit{x}=1$ on the real number line. Use the 'Waypoints' name-value pair to specify the piecewise segments in the contour.

f = @(z) 1./(2.*z-1);
contour_segments = [1+1i 0+1i 0-1i 1-1i];
q = quadgk(f,1,1,'Waypoints',contour_segments)

q = 
-0.0000 + 3.1416i

Examine Absolute Error

Use quadgk to evaluate an oscillatory integrand that is difficult to evaluate.

Evaluate the integral

$$\mathit{Q}={\int }_0^{\pi }\sin \left(20000\pi \mathit{x}\right)\mathrm{dx}.$$

The integrand oscillates very quickly, so it is difficult to evaluate. Use quadgk to evaluate the integral, and specify two outputs to examine how close the error tolerances are to being met.

fun = @(x) sin(2e4*pi*x);
[Q,errbnd] = quadgk(fun,0,pi)

Warning: Reached the limit on the maximum number of intervals in use. Approximate bound on error is  5.7e-01. The integral may not exist, or it may be difficult to approximate numerically. Increase MaxIntervalCount to 1272 to enable QUADGK to continue for another iteration.

Q = 
-0.0082

errbnd = 
0.5723

The warning message indicates how to adjust MaxIntervalCount to allow for another iteration in the solution process.

Solve the integral again, but specify MaxIntervalCount as 1e5. With many more intervals, quadgk is able to meet the absolute error tolerance for the problem (1e-10 for double precision).

[Q,errbnd] = quadgk(fun,0,pi,'MaxIntervalCount',1e5)

Q = 
1.6656e-06

errbnd = 
2.6323e-12

Input Arguments

fun — Integrand

function handle

Integrand, specified as a function handle that defines the function to be integrated from a to b.

For scalar-valued problems, the function y = fun(x) must accept a vector argument x and return a vector result y, where y is the integrand evaluated at each element of x. This requirement generally means that fun must use array operators (.^, .*, …) instead of matrix operators (^, *, …).

Parameterizing Functions explains how to provide additional parameters to the function fun, if necessary.

Example: q = quadgk(@(x) exp(1-x.^2),a,b) integrates an anonymous function handle.

Example: q = quadgk(@myFun,a,b) integrates the function myFun, which is saved as a file.

Data Types: function_handle

a, b — Integration limits (as separate arguments)

scalars

Integration limits, specified as separate arguments of real or complex scalars. The limits a and b can be -Inf or Inf. If both are finite, they can be complex. If at least one is complex, the integral is approximated over a straight line path from a to b in the complex plane.

Example: quadgk(fun,0,1) integrates fun from 0 to 1.

Data Types: single | double
Complex Number Support: Yes

Name-Value Arguments

Specify optional pairs of arguments as Name1=Value1,...,NameN=ValueN, where Name is the argument name and Value is the corresponding value. Name-value arguments must appear after other arguments, but the order of the pairs does not matter.

Before R2021a, use commas to separate each name and value, and enclose Name in quotes.

Example: q = quadgk(fun,a,b,'Waypoints',[0.1 1.1 2.1]) uses the 'Waypoints' option to specify a few points of interest where the integrand should be evaluated.

AbsTol — Absolute error tolerance

1e-10 (double), 1e-5 (single) (default) | nonnegative real number

Absolute error tolerance, specified as the comma-separated pair consisting of 'AbsTol' and a nonnegative real number. quadgk uses the absolute error tolerance to limit an estimate of the absolute error, |q – I|, where q is the computed value of the integral and I is the (unknown) exact value. quadgk might provide more decimal places of precision if you decrease the absolute error tolerance.

quadgk attempts to satisfy

errbnd <= max(AbsTol,RelTol*abs(q))

This relation is absolute error control when |q| is sufficiently small and relative error control when |q| is larger. For pure absolute error control, use 'AbsTol' > 0 and 'RelTol'= 0. For pure relative error control use 'RelTol' > 0 and 'AbsTol' = 0. Except when using pure absolute error control, the minimum relative tolerance is 'RelTol' >= 100*eps(class(q)).

Example: quadgk(fun,a,b,'AbsTol',1e-12) sets the absolute error tolerance to approximately 12 decimal places of accuracy.

Example: quadgk(fun,a,b,'AbsTol',tol,'RelTol',0) uses a purely absolute error control, requiring that errbnd <= tol.

Data Types: single | double

RelTol — Relative error tolerance

1e-6 (double), 1e-4 (single) (default) | nonnegative real number

Relative error tolerance, specified as the comma-separated pair consisting of 'RelTol' and a nonnegative real number. quadgk uses the relative error tolerance to limit an estimate of the relative error, |q - I|/|I|, where q is the computed value of the integral and I is the (unknown) exact value. quadgk might provide more significant digits of precision if you decrease the relative error tolerance.

quadgk attempts to satisfy

errbnd <= max(AbsTol,RelTol*abs(q))

This relation is absolute error control when |q| is sufficiently small and relative error control when |q| is larger. For pure absolute error control, use 'AbsTol' > 0 and 'RelTol'= 0. For pure relative error control use 'RelTol' > 0 and 'AbsTol' = 0. Except when using pure absolute error control, the minimum relative tolerance is 'RelTol' >= 100*eps(class(q)).

Example: quadgk(fun,a,b,'RelTol',1e-9) sets the relative error tolerance to approximately 9 significant digits.

Example: quadgk(fun,a,b,'AbsTol',0,'RelTol',tol) uses a purely relative error tolerance, requiring that errbnd <= |I|*tol.

Data Types: single | double

Waypoints — Integration waypoints

vector

Integration waypoints, specified as the comma-separated pair consisting of 'Waypoints' and a vector of real or complex numbers. Use waypoints to indicate points in the integration interval that you would like the integrator to use in the initial mesh:

• Add more evaluation points near interesting features of the function, such as a local extrema.

• Integrate efficiently across discontinuities of the integrand by specifying the locations of the discontinuities.

• Perform complex contour integrations by specifying complex numbers as waypoints. If xmin, xmax, or any entry of the waypoints vector is complex, then the integration is performed over a sequence of straight line paths in the complex plane. In this case, all of the integration limits and waypoints must be finite.

Do not use waypoints to specify singularities. Instead, split the interval and add the results of separate integrations with the singularities at the endpoints.

Example: 'Waypoints',[1+1i,1-1i] specifies two complex waypoints along the interval of integration.

Data Types: single | double
Complex Number Support: Yes

MaxIntervalCount — Maximum number of intervals allowed

650 (default) | scalar

Maximum number of intervals allowed, specified as a scalar. This option limits the number of intervals that quadgk uses at any one time after the first iteration. A warning is issued if quadgk returns early because of this limit. Routinely increasing this value is not recommended, but it may be appropriate when errbnd is small enough that the desired accuracy has nearly been achieved.

Example: quadgk(fun,a,b,'MaxIntervalCount',700)

Output Arguments

q — Value of integral

scalar

Value of integral, returned as a scalar.

errbnd — Approximate upper bound on absolute error

scalar

Approximate upper bound on absolute error, returned as a scalar. The approximate upper bound on absolute error in the integration is errbnd = |q – I|, where q is the computed value of the integral and I is the (unknown) exact value. quadgk attempts to satisfy

errbnd <= max(AbsTol,RelTol*abs(q))

Specify this output argument to see how well the integration meets the AbsTol and RelTol error tolerances. In cases where errbnd is close to the desired value, you might be able to reach the desired value by increasing the value of MaxIntervalCount.

Tips

• quadgk and integral use essentially the same integration method. You should generally use integral rather than quadgk. However, you can use quadgk to:

  • Monitor solution accuracy with the errbnd output argument.

  • Specify a large value for MaxIntervalCount when integral warns about reaching the maximum number of intervals.

• quadgk can integrate functions that are singular at finite endpoints if the singularities are not too strong. For example, it can integrate functions that behave at an endpoint c like log|x-c| or |x-c|p for p >= -1/2. If the function is singular at points inside the integration limits [a b], then write the integral as a sum of integrals over subintervals with the singular points as endpoints, compute them with quadgk, and add the results.

• If the interval is infinite, $$\left[a,\infty \right)$$, then for the integral of fun(x) to exist, fun(x) must decay as x approaches infinity, and quadgk requires it to decay rapidly.

References

[1] Shampine, L.F. "Vectorized Adaptive Quadrature in MATLAB^®." Journal of Computational and Applied Mathematics. Vol. 211, 2008, pp.131–140.

Extended Capabilities

C/C++ Code Generation
Generate C and C++ code using MATLAB® Coder™.

Version History

Introduced in R2007b