NAG FL Interface
d01arf (dim1_indef)
1
Purpose
d01arf computes definite and indefinite integrals over a finite range to a specified relative or absolute accuracy, using the method described in
Patterson (1968).
2
Specification
Fortran Interface
Subroutine d01arf ( 
a, b, f, relacc, absacc, maxrul, iparm, acc, ans, n, alpha, ifail) 
Integer, Intent (In) 
:: 
maxrul, iparm 
Integer, Intent (Inout) 
:: 
ifail 
Integer, Intent (Out) 
:: 
n 
Real (Kind=nag_wp), External 
:: 
f 
Real (Kind=nag_wp), Intent (In) 
:: 
a, b, relacc, absacc 
Real (Kind=nag_wp), Intent (Inout) 
:: 
alpha(390) 
Real (Kind=nag_wp), Intent (Out) 
:: 
acc, ans 

C Header Interface
#include <nag.h>
void 
d01arf_ (const double *a, const double *b, double (NAG_CALL *f)(const double *x), const double *relacc, const double *absacc, const Integer *maxrul, const Integer *iparm, double *acc, double *ans, Integer *n, double alpha[], Integer *ifail) 

C++ Header Interface
#include <nag.h> extern "C" {
void 
d01arf_ (const double &a, const double &b, double (NAG_CALL *f)(const double &x), const double &relacc, const double &absacc, const Integer &maxrul, const Integer &iparm, double &acc, double &ans, Integer &n, double alpha[], Integer &ifail) 
}

The routine may be called by the names d01arf or nagf_quad_dim1_indef.
3
Description
d01arf evaluates definite and indefinite integrals of the form:
using the method described in
Patterson (1968).
3.1
Definite Integrals
In this case
d01arf must be called with
${\mathbf{iparm}}=0$. By linear transformation the integral is changed to
where
and is then approximated by an
$n$point quadrature rule
where
${w}_{k}$ are the weights and
${x}_{k}$ are the abscissae.
The routine uses a family of nine interlacing rules based on the optimal extension of the threepoint Gauss rule. These rules use $1$, $3$, $7$, $15$, $31$, $63$, $127$, $255$ and $511$ points and have respective polynomial integrating degrees $1$, $5$, $11$, $23$, $47$, $95$, $191$, $383$ and $767$. Each rule has the property that the next in sequence includes all the points of its predecessor and has the greatest possible increase in integrating degree.
The integration method is based on the successive application of these rules until the absolute value of the difference of two successive results differs by not more than
absacc, or relatively by not more than
relacc. The result of the last rule used is taken as the value of the integral (
ans), and the absolute difference of the results of the last two rules used is taken as an estimate of the absolute error (
acc). Due to their interlacing form no integrand evaluations are wasted in passing from one rule to the next.
3.2
Indefinite Integrals
Suppose the value of the integral
is required for a number of subintervals
$\left[c,d\right]$, all of which lie in an interval
$\left[a,b\right]$.
In this case d01arf should first be called with the argument ${\mathbf{iparm}}=1$ and the interval set to $\left[a,b\right]$. The routine then calculates the integral over $\left[a,b\right]$ and the Legendre expansion of the integrand, using the same integrand values. If the routine is subsequently called with ${\mathbf{iparm}}=2$ and the interval set to $\left[c,d\right]$, the integral over $\left[c,d\right]$ is calculated by analytical integration of the Legendre expansion, without further evaluations of the integrand.
For the interval
$\left[1,1\right]$ the expansion takes the form
where
${P}_{i}\left(x\right)$ is the order
$i$ Legendre polynomial. Assuming that the integral over the full range
$\left[1,1\right]$ was evaluated to the required accuracy using an
$n$point rule, then the coefficients
are evaluated by that same rule, up to
The accuracy for indefinite integration should be of the same order as that obtained for the definite integral over the full range. The indefinite integrals will be exact when
$F\left(x\right)$ is a polynomial of degree
$\text{}\le m$.
4
References
Patterson T N L (1968) The Optimum addition of points to quadrature formulae Math. Comput. 22 847–856
5
Arguments

1:
$\mathbf{a}$ – Real (Kind=nag_wp)
Input

On entry: $a$, the lower limit of integration.

2:
$\mathbf{b}$ – Real (Kind=nag_wp)
Input

On entry: $b$, the upper limit of integration. It is not necessary that $a<b$.

3:
$\mathbf{f}$ – real (Kind=nag_wp) Function, supplied by the user.
External Procedure

f must return the value of the integrand at a given point.
The specification of
f is:
Fortran Interface
Real (Kind=nag_wp) 
:: 
f 
Real (Kind=nag_wp), Intent (In) 
:: 
x 

C Header Interface
double 
f_ (const double *x) 

C++ Header Interface
#include <nag.h> extern "C" {
double 
f_ (const double &x) 
}


1:
$\mathbf{x}$ – Real (Kind=nag_wp)
Input

On entry: the point in $\left[a,b\right]$ at which the integrand $f$ must be evaluated.
f must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which
d01arf is called. Arguments denoted as
Input must
not be changed by this procedure.
Note: f should not return floatingpoint NaN (Not a Number) or infinity values, since these are not handled by
d01arf. If your code inadvertently
does return any NaNs or infinities,
d01arf is likely to produce unexpected results.
If
${\mathbf{iparm}}=2$,
f is not called.

4:
$\mathbf{relacc}$ – Real (Kind=nag_wp)
Input

On entry: the relative accuracy required. If convergence according to absolute accuracy is required,
relacc should be set to zero (but see also
Section 7). If
${\mathbf{relacc}}<0.0$, its absolute value is used.
If
${\mathbf{iparm}}=2$,
relacc is not used.

5:
$\mathbf{absacc}$ – Real (Kind=nag_wp)
Input

On entry: the absolute accuracy required. If convergence according to relative accuracy is required,
absacc should be set to zero (but see also
Section 7). If
${\mathbf{absacc}}<0.0$, its absolute value is used.
If
${\mathbf{iparm}}=2$,
absacc is not used.

6:
$\mathbf{maxrul}$ – Integer
Input

On entry: the maximum number of successive rules that may be used.
Constraint:
$1\le {\mathbf{maxrul}}\le 9$. If
maxrul is outside these limits, the value
$9$ is assumed.
If
${\mathbf{iparm}}=2$,
maxrul is not used.

7:
$\mathbf{iparm}$ – Integer
Input

On entry: indicates the task to be performed by the routine.
 ${\mathbf{iparm}}=0$
 Only the definite integral over $\left[a,b\right]$ is evaluated.
 ${\mathbf{iparm}}=1$
 As well as the definite integral, the expansion of the integrand in Legendre polynomials over $\left[a,b\right]$ is calculated, using the same values of the integrand as used to compute the integral. The expansion coefficients, and some other quantities, are returned in alpha for later use in computing indefinite integrals.
 ${\mathbf{iparm}}=2$
 $f\left(t\right)$ is integrated analytically over $\left[a,b\right]$ using the previously computed expansion, stored in alpha. No further evaluations of the integrand are required. The routine must previously have been called with ${\mathbf{iparm}}=1$ and the interval $\left[a,b\right]$ must lie within that specified for the previous call. In this case only the arguments a, b, iparm, ans, alpha and ifail are used.
Constraint:
${\mathbf{iparm}}=0$, $1$ or $2$.

8:
$\mathbf{acc}$ – Real (Kind=nag_wp)
Output

On exit: if
${\mathbf{iparm}}=0$ or
$1$,
acc contains the absolute value of the difference between the last two successive estimates of the integral. This may be used as a measure of the accuracy actually achieved.
If
${\mathbf{iparm}}=2$,
acc is not used.

9:
$\mathbf{ans}$ – Real (Kind=nag_wp)
Output

On exit: the estimated value of the integral.

10:
$\mathbf{n}$ – Integer
Output

On exit: when
${\mathbf{iparm}}=0$ or
$1$,
n contains the number of integrand evaluations used in the calculation of the integral.
If
${\mathbf{iparm}}=2$,
n is not used.

11:
$\mathbf{alpha}\left(390\right)$ – Real (Kind=nag_wp) array
Input/Output

On entry: if
${\mathbf{iparm}}=2$,
alpha must contain the coefficients of the Legendre expansions of the integrand, as returned by a previous call of
d01arf with
${\mathbf{iparm}}=1$ and a range containing the present range.
If
${\mathbf{iparm}}=0$ or
$1$,
alpha need not be set on entry.
On exit: if
${\mathbf{iparm}}=1$, the first
$m$ elements of
alpha hold the coefficients of the Legendre expansion of the integrand, and the value of
$m$ is stored in
${\mathbf{alpha}}\left(390\right)$.
alpha must not be changed between a call with
${\mathbf{iparm}}=1$ and subsequent calls with
${\mathbf{iparm}}=2$.
If
${\mathbf{iparm}}=2$, the first
$m$ elements of
alpha are unchanged on exit.

12:
$\mathbf{ifail}$ – Integer
Input/Output

On entry:
ifail must be set to
$0$,
$1$ or
$1$ to set behaviour on detection of an error; these values have no effect when no error is detected.
A value of $0$ causes the printing of an error message and program execution will be halted; otherwise program execution continues. A value of $1$ means that an error message is printed while a value of $1$ means that it is not.
If halting is not appropriate, the value
$1$ or
$1$ is recommended. If message printing is undesirable, then the value
$1$ is recommended. Otherwise, the value
$1$ is recommended since useful values can be provided in some output arguments even when
${\mathbf{ifail}}\ne {\mathbf{0}}$ on exit.
When the value $\mathbf{1}$ or $\mathbf{1}$ is used it is essential to test the value of ifail on exit.
On exit:
${\mathbf{ifail}}={\mathbf{0}}$ unless the routine detects an error or a warning has been flagged (see
Section 6).
6
Error Indicators and Warnings
If on entry
${\mathbf{ifail}}=0$ or
$1$, explanatory error messages are output on the current error message unit (as defined by
x04aaf).
Errors or warnings detected by the routine:
Note: in some cases d01arf may return useful information.
 ${\mathbf{ifail}}=1$

All
maxrul rules have been used but the integral has not converged to the accuracy requested.
${\mathbf{maxrul}}=\u2329\mathit{\text{value}}\u232a$.
Failure in convergence during the run with ${\mathbf{iparm}}=1$ in which the Legendre expansion was created.
 ${\mathbf{ifail}}=2$

On entry, ${\mathbf{iparm}}=\u2329\mathit{\text{value}}\u232a$.
Constraint: ${\mathbf{iparm}}=0$, $1$ or $2$.
 ${\mathbf{ifail}}=3$

The routine was called with ${\mathbf{iparm}}=2$, but a previous call with ${\mathbf{iparm}}=1$ has been omitted, or made with zero integration interval.
 ${\mathbf{ifail}}=4$

The routine was called with ${\mathbf{iparm}}=2$, but the interval for indefinite integration is not contained within the interval specified in the previous call with ${\mathbf{iparm}}=1$.
 ${\mathbf{ifail}}=99$
An unexpected error has been triggered by this routine. Please
contact
NAG.
See
Section 7 in the Introduction to the NAG Library FL Interface for further information.
 ${\mathbf{ifail}}=399$
Your licence key may have expired or may not have been installed correctly.
See
Section 8 in the Introduction to the NAG Library FL Interface for further information.
 ${\mathbf{ifail}}=999$
Dynamic memory allocation failed.
See
Section 9 in the Introduction to the NAG Library FL Interface for further information.
7
Accuracy
The relative or absolute accuracy required is specified by you in the variables
relacc or
absacc.
d01arf will terminate whenever either the relative accuracy specified by
relacc or the absolute accuracy specified by
absacc is reached. One or other of these criteria may be ‘forced’ by setting the argument for the other to zero. If both
relacc and
absacc are specified as zero, then the routine uses the value
$10.0\times \left(\mathit{machineprecision}\right)$ for
relacc.
If on exit ${\mathbf{ifail}}={\mathbf{0}}$, then it is likely that the result is correct to one or other of these accuracies. If on exit ${\mathbf{ifail}}={\mathbf{1}}$, then it is likely that neither of the requested accuracies has been reached.
When you have no prior idea of the magnitude of the integral, it is possible that an unreasonable accuracy may be requested, e.g., a relative accuracy for an integral which turns out to be zero, or a small absolute accuracy for an integral which turns out to be very large. Even if failure is reported in such a case, the value of the integral may still be satisfactory. The device of setting the other ‘unused’ accuracy argument to a small positive value (e.g., ${10}^{9}$ for an implementation of $11$digit precision) rather than zero, may prevent excessive calculation in such a situation.
To avoid spurious convergence, it is recommended that relative accuracies larger than about ${10}^{3}$ be avoided.
8
Parallelism and Performance
d01arf is not threaded in any implementation.
The time taken by d01arf depends on the complexity of the integrand and the accuracy required.
This routine uses the Patterson method over the whole integration interval and should therefore be suitable for well behaved functions. However, for very irregular functions it would be more efficient to submit the differently behaved regions separately for integration.
10
Example
This example evaluates the following integrals

(i)Definite integral only $\left({\mathbf{iparm}}=0\right)$ for

(ii)Definite integral together with expansion coefficients $\left({\mathbf{iparm}}=1\right)$ for

(iii)Indefinite integral using previous expansion $\left({\mathbf{iparm}}=2\right)$ for
10.1
Program Text
10.2
Program Data
None.
10.3
Program Results