Approximate the definite integral of a function over a finite or infinite interval.

Framework

- Accelerate

## Overview

Quadrature provides an approximation of the definite integral of a function over a finite or infinite interval.

*Quadrature* is an historic term for determining the area under a curve. Often this was done by breaking the area into smaller shapes whose area could be easily calculated (such as rectangles), and summing these smaller areas to obtain an approximate result.

In modern terms this process is called *definite integration*. The Accelerate framework’s Quadrature functionality provides an approximation of the definite integral of a function over a finite or infinite interval, performed by evaluating the function at a series of points within the interval. The `quadrature`

function performs this calculation using any one of three algorithms described in the Integrators section below.

### The Integration Callback

Note

To avoid confusion over the word *function*, this document refers to the mathematical function that is to be integrated as the *integrand*.

To represent the integrand, use a C function of the following type, defined in `integration`

:

This function is the *integration callback*. This is a function that processes values in an input array, and produces corresponding result values in an output array; the function should not do anything else, and it must return `void`

. The input values are *x* values within the interval over which the integrand is being integrated, and the output values are the corresponding values *y = integrand(x)* at those points.

There are two other parameters: a `void*`

pointer `arg`

that you will supply at call time, and a size for the arrays (for details, consult the header file `integration`

). The pointer `arg`

is available in case you need to reference some outside object from inside your callback; however if your callback only needs the `x`

and `y`

arguments, you can ignore `arg`

.

When the integration callback has been defined, you can package it in a struct of type `quadrature`

, defined in `integration`

:

where `integration`

is the integration callback and `my`

is a `void*`

pointer value to be passed to `integration`

as the first parameter, `arg`

. If the integration callback has been written to ignore this value, just pass `NULL`

.

The struct `f`

can now be passed as the first argument to `quadrature`

, which will supply the input and output arrays, fill the input array with *x* values, and perform the integration by calling `integration`

as many times as necessary.

### Integration Options

Options that control the logic of the integration are specified as fields of a struct of type `quadrature`

, defined in `integration`

. The fields are interpreted as follows:

`integrator`

: A constant that identifies one of three available integration algorithms (see the Integrators section below):`QUADRATURE`

selects the QNG algorithm._INTEGRATE _QNG `QUADRATURE`

selects the QAG algorithm._INTEGRATE _QAG `QUADRATURE`

selects the QAGS algorithm._INTEGRATE _QAGS

`abs`

: A value of type_tolerance `double`

; requested absolute tolerance on the result.`rel`

: A value of type_tolerance `double`

; requested relative tolerance on the result.`qag`

: Number of points per subinterval. Used by the QAG integrator only; other integrators ignore this value. Can be 0, 15, 21, 31, 41, 51, 61. 0 maps to the default 21._points _per _interval `max`

: When you do not provide a workspace (see the Providing a Workspace section below), this is the maximum number of subintervals in the subdivision used by QAG and QAGS integrators._intervals `max`

is ignored by QAG and QAGS if you provide a workspace, and always ignored by QNG._intervals

### Integrators

Quadrature has three different integrators available to perform the integration, called QNG, QAG, and QAGS. Each integrator is a specific implementation of an algorithm of the same name. The algorithms are all variants of the well-known Gauss-Kronrod method, and their integrators are C ports of the corresponding routines in the QUADPACK library:

QNG is a simple nonadaptive automatic integrator using Gauss-Kronrod-Patterson quadrature coefficients. It evaluates 21, 43, or 87 points in the interval until the requested accuracy is reached.

QAG is globally adaptive – it divides the integration interval into a number of subintervals depending on the memory available and/or the

`max`

option (see the Providing a Workspace section below). In each interval it uses 0, 15, 21, 31, 41, 51, or 61 points depending on the_intervals `qag`

option._points _per _interval QAGS is globally adaptive and permits infinite bounds on the integration interval. Within each subinterval it uses 21 points if both bounds are finite, or 15 if one or both bounds are infinite. The algorithm is accelerated by Peter Wynn’s epsilon algorithm. If at least one of the interval bounds is infinite, this is equivalent to the QUADPACK QAGI routine. Otherwise, this is equivalent to the QUADPACK QAGS routine.

To select an integrator for a particular integrand, the following decision tree is recommended:

**Integration over a finite region**

If performance is not a concern and you don’t know much about the specifics of the problem, use QAGS.

Otherwise, if the integrand is smooth, use QNG – or QAG if the requested tolerance couldn't be reached with QNG.

Otherwise, if there are discontinuities or singularities of the integrand or of its derivative, and you know where they are, split the integration range at these points and analyze each subinterval.

Otherwise, if the integrand has end point singularities, use QAGS.

Otherwise, if the integrand has an oscillatory behavior of nonspecific type, and no singularities, use QAG with 61 points per interval.

Otherwise, use QAGS.

**Integration over an infinite region**

If the integrand decays rapidly toward zero, truncate the interval and use the finite interval decision tree.

Otherwise, if you are not constrained by computer time, and do not wish to analyze the problem further, use QAGS.

Otherwise, if the integrand has a non-smooth behavior in the range, and you know where it occurs, split off these regions and use the appropriate finite range routines to integrate over them. Then begin this tree again to handle the remainder of the region.

Otherwise, truncation of the interval, or application of a suitable transformation for reducing the problem to a finite range may be possible.

### Providing a Workspace

Optionally, you may allocate a memory area as workspace for quadrature. This is not a requirement and you should generally not provide a workspace unless you need precise control over memory usage. The QAG and QAGS algorithms require a workspace and will automatically allocate it if you do not provide one; the size of the automatically allocated workspace will depend on the `max`

option. The QNG algorithm does not require a workspace.

By allocating the workspace yourself, you have control over the amount of memory used; there is a tradeoff between memory usage and performance.

The QAG algorithm requires, for a given number of subintervals, at least that number times `QUADRATURE`

bytes in the workspace. If you provide more memory, the algorithm will be able to use more subintervals.

The QAGS algorithm requires, for a given number of subintervals, at least that number times `QUADRATURE`

bytes in the workspace. If you provide more memory, the algorithm will be able to use more subintervals.