next up previous index
Next: D114 Adaptive Multidimensional Up: CERNLIB Previous: D110 Gaussian Quadrature

D113 Adaptive Complex Integration Along a Line Segment

Routine ID: D113
Author(s): G.A. ErskineLibrary: MATHLIB
Submitter: K.S. KölbigSubmitted: 07.12.1970
Language: FortranRevised: 15.03.1993

Function subprograms CGAUSS and WGAUSS compute, to an attempted specified accuracy, the value of the complex integral

I=∫ABf(z)dz.

The path of integration is the directed line segment AB in the complex z-plane. The function f(z) must be single-valued on this segment.

The double-precision version WGAUSS is available only on computers which support a COMPLEX*16 Fortran data type.

Structure:

FUNCTION subprograms
User Entry Names: CGAUSS, WGAUSS
Files Referenced: Unit 6
External References: MTLMTR (N002), ABEND (Z035), user-supplied FUNCTION subprogram

Usage:

In any arithmetic expression, CGAUSS(F,A,B,EPS) or WGAUSS(F,A,B,EPS)

has the approximate value of the integral I.

F
Name of a user-supplied FUNCTION subprogram, declared EXTERNAL in the calling program. This subroutine must set F(Z)= f(Z) .
A,B
End-points of integration interval.
EPS
Accuracy parameter (see Accuracy).
CGAUSS is of type COMPLEX, WGAUSS is of type COMPLEX*16, and the arguments F, A, B, and Z (in F) have the same type as the function name. EPS is of type REAL for CGAUSS and of type DOUBLE PRECISION for WGAUSS.

Method:

For any line segment [a,b] we define g8(a,b) and g16(a,b) to be the 8-point and 16-point Gaussian quadrature approximations to

abf(z)dz

and define

r(a,b) ={|g16(a,b) - g8(a,b)|1+|g16(a,b)|}.

Then, with G = CGAUSS or WGAUSS,

G =∑i=1kg16(zi-1,zi),

where, starting with z0=A and finishing with zk=B , the subdivision points zi (i=1,2,...) are given by

zi= zi-1+ λ(B-zi-1),

with λ equal to the first member of the sequence 1,1/2,1/4,... for which r(zi-1,zi) < EPS . If, at any stage in the process of subdivision, the ratio

q=|{zi-zi-1B-A}|

is so small that 1+0.005q is indistinguishable from 1 to machine accuracy, an error exit occurs with the function value set equal to zero.

Accuracy:

Unless there is severe cancellation of positive and negative values of f(z) over the interval [A,B] , the argument EPS may be considered as specifying a bound on the relative error of I in the case |I| > 1, and a bound on the absolute error in the case |I|<1. More precisely, if k is the number of sub-intervals contributing to the approximation (see Method), and if

Iabs= ∫AB|f(z)|dz,

then the relation

{|G - I|Iabs+k}< EPS

will nearly always be true, provided the routine terminates without printing an error message. For functions f having no singularities in the closed interval [A,B]

the accuracy will usually be much higher than this.

Error handling:

Error D113.1: The requested accuracy (see Method) cannot be obtained. The function value is set equal to zero, and a message is written on Unit 6, unless subroutine MTLSET (N002) has been called.

Notes:

Values of the function f(z) at the end-points of the line segment A and B are not required. The subprogram may therefore be used when these values are undefined.

D114



next up previous index
Next: D114 Adaptive Multidimensional Up: CERNLIB Previous: D110 Gaussian Quadrature


Janne Saarela
Mon Apr 3 15:06:23 METDST 1995