hide long namesshow long names
hide short namesshow short names
Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

NAG Toolbox: nag_quad_2d_fin (d01da)

Purpose

nag_quad_2d_fin (d01da) attempts to evaluate a double integral to a specified absolute accuracy by repeated applications of the method described by Patterson (1968) and Patterson (1969).

Syntax

[ans, npts, ifail] = d01da(ya, yb, phi1, phi2, f, absacc)
[ans, npts, ifail] = nag_quad_2d_fin(ya, yb, phi1, phi2, f, absacc)

Description

nag_quad_2d_fin (d01da) attempts to evaluate a definite integral of the form
bφ2(y)
I = f(x,y)dxdy
aφ1(y)
I= ab ϕ1(y) ϕ2(y) f(x,y) dx dy
where aa and bb are constants and φ1(y)ϕ1(y) and φ2(y)ϕ2(y) are functions of the variable yy.
The integral is evaluated by expressing it as
b φ2(y)
I = F(y)dy,   where  F(y) = f(x,y)dx.
a φ1(y)
I=abF(y)dy,   where  F(y)= ϕ1(y) ϕ2(y) f(x,y)dx.
Both the outer integral II and the inner integrals F(y)F(y) are evaluated by the method, described by
Patterson (1968) and Patterson (1969), of the optimum addition of points to Gauss quadrature formulae.
This method uses a family of interlacing common point formulae. Beginning with the 33-point Gauss rule, formulae using 77, 1515, 3131, 6363, 127127 and finally 255255 points are derived. Each new formula contains all the points of the earlier formulae so that no function evaluations are wasted. Each integral is evaluated by applying these formulae successively until two results are obtained which differ by less than the specified absolute accuracy.

References

Patterson T N L (1968) On some Gauss and Lobatto based integration formulae Math. Comput. 22 877–881
Patterson T N L (1969) The optimum addition of points to quadrature formulae, errata Math. Comput. 23 892

Parameters

Compulsory Input Parameters

1:     ya – double scalar
aa, the lower limit of the integral.
2:     yb – double scalar
bb, the upper limit of the integral. It is not necessary that a < ba<b.
3:     phi1 – function handle or string containing name of m-file
phi1 must return the lower limit of the inner integral for a given value of yy.
[result] = phi1(y)

Input Parameters

1:     y – double scalar
The value of yy for which the lower limit must be evaluated.

Output Parameters

1:     result – double scalar
The result of the function.
4:     phi2 – function handle or string containing name of m-file
phi2 must return the upper limit of the inner integral for a given value of yy.
[result] = phi2(y)

Input Parameters

1:     y – double scalar
The value of yy for which the upper limit must be evaluated.

Output Parameters

1:     result – double scalar
The result of the function.
5:     f – function handle or string containing name of m-file
f must return the value of the integrand ff at a given point.
[result] = f(x, y)

Input Parameters

1:     x – double scalar
2:     y – double scalar
The coordinates of the point (x,y)(x,y) at which the integrand ff must be evaluated.

Output Parameters

1:     result – double scalar
The result of the function.
6:     absacc – double scalar
The absolute accuracy requested.

Optional Input Parameters

None.

Input Parameters Omitted from the MATLAB Interface

None.

Output Parameters

1:     ans – double scalar
The estimated value of the integral.
2:     npts – int64int32nag_int scalar
The total number of function evaluations.
3:     ifail – int64int32nag_int scalar
ifail = 0ifail=0 unless the function detects an error (see [Error Indicators and Warnings]).

Error Indicators and Warnings

Note: nag_quad_2d_fin (d01da) may return useful information for one or more of the following detected errors or warnings.
Errors or warnings detected by the function:

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

W ifail = 1ifail=1
This indicates that 255255 points have been used in the outer integral and convergence has not been obtained. All the inner integrals have, however, converged. In this case ans may still contain an approximate estimate of the integral.
W ifail = 10 × nifail=10×n
This indicates that the outer integral has converged but nn inner integrals have failed to converge with the use of 255255 points. In this case ans may still contain an approximate estimate of the integral, but its reliability will decrease as nn increases.
W ifail = 10 × n + 1ifail=10×n+1
This indicates that both the outer integral and nn of the inner integrals have not converged. ans may still contain an approximate estimate of the integral, but its reliability will decrease as nn increases.

Accuracy

The absolute accuracy is specified by the variable absacc. If, on exit, ifail = 0ifail=0 then the result is most likely correct to this accuracy. Even if ifail is nonzero on exit, it is still possible that the calculated result could differ from the true value by less than the given accuracy.

Further Comments

The time taken by nag_quad_2d_fin (d01da) depends upon the complexity of the integrand and the accuracy requested.
With Patterson's method accidental convergence may occasionally occur, when two estimates of an integral agree to within the requested accuracy, but both estimates differ considerably from the true result. This could occur in either the outer integral or in one or more of the inner integrals.
If it occurs in the outer integral then apparent convergence is likely to be obtained with considerably fewer integrand evaluations than may be expected. If it occurs in an inner integral, the incorrect value could make the function F(y)F(y) appear to be badly behaved, in which case a very large number of pivots may be needed for the overall evaluation of the integral. Thus both unexpectedly small and unexpectedly large numbers of integrand evaluations should be considered as indicating possible trouble. If accidental convergence is suspected, the integral may be recomputed, requesting better accuracy; if the new request is more stringent than the degree of accidental agreement (which is of course unknown), improved results should be obtained. This is only possible when the accidental agreement is not better than machine accuracy. It should be noted that the function requests the same accuracy for the inner integrals as for the outer integral. In practice it has been found that in the vast majority of cases this has proved to be adequate for the overall result of the double integral to be accurate to within the specified value.
The function is not well-suited to non-smooth integrands, i.e., integrands having some kind of analytic discontinuity (such as a discontinuous or infinite partial derivative of some low-order) in, on the boundary of, or near, the region of integration. Warning: such singularities may be induced by incautiously presenting an apparently smooth interval over the positive quadrant of the unit circle, RR 
I = R(x + y)dx dy.
I=R(x+y)dx dy.
This may be presented to nag_quad_2d_fin (d01da) as
1 sqrt(1y2) 1
I = dy(x + y)dx = ((1/2)(1y2) + y×sqrt(1y2))dy
0 0 0
I=01 dy 01-y2 (x+y)dx=01 (12(1-y2)+y1-y2) dy
but here the outer integral has an induced square-root singularity stemming from the way the region has been presented to nag_quad_2d_fin (d01da). This situation should be avoided by re-casting the problem. For the example given, the use of polar coordinates would avoid the difficulty:
1 π/2
I = drr2(cosυ + sinυ)dυ.
0 0
I=01dr0π2r2(cosυ+sinυ)dυ.

Example

function nag_quad_2d_fin_example
ya = 0;
yb = 1;
absacc = 1e-06;
[ans, npts, ifail] = nag_quad_2d_fin(ya, yb, @phi1, @phi2, @f, absacc)

function result = phi1(y)
  result=0;
function result = phi2(y)
  result = sqrt(1-y^2);
function result = f(x,y)
  result=x+y;
 

ans =

    0.6667


npts =

                  189


ifail =

                    0


function d01da_example
ya = 0;
yb = 1;
absacc = 1e-06;
[ans, npts, ifail] = d01da(ya, yb, @phi1, @phi2, @f, absacc)

function result = phi1(y)
  result=0;
function result = phi2(y)
  result = sqrt(1-y^2);
function result = f(x,y)
  result=x+y;
 

ans =

    0.6667


npts =

                  189


ifail =

                    0



PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2013