Integrates a function f(x)/(x-c) in the Cauchy principal value sense.
F — User-supplied
FUNCTION to be
integrated. The form is F(X),
where
X − Independent
variable.
(Input)
F − The function
value. (Output)
F must be declared
EXTERNAL in the
calling program.
A — Lower limit of integration. (Input)
B — Upper limit of integration. (Input)
C — Singular
point. (Input)
C must not equal A or B.
RESULT — Estimate of the integral from A to B of F(X)/(X − C). (Output)
ERRABS — Absolute
accuracy desired. (Input)
Default: ERRABS = 1.e-3 for
single precision and 1.d-8 for double precision.
ERRREL — Relative
accuracy desired. (Input)
Default: ERREL =1.e-3 for
single precision and 1.d-8 for double precision.
ERREST — Estimate of the absolute value of the error. (Output)
Generic: CALL QDAWC (F, A, B, C, RESULT [,…])
Specific: The specific interface names are S_QDAWC and D_QDAWC.
Single: CALL QDAWC (F, A, B, C, ERRABS, ERRREL, RESULT, ERREST)
Double: The double precision name is DQDAWC.
The routine QDAWC
uses a globally adaptive scheme in an attempt to reduce the absolute error. This
routine computes integrals whose integrands have the special form
w(x) f(x), where
w(x) =
1/(x − c).
If c lies in the interval of integration, then the integral is
interpreted as a Cauchy principal value. A combination of modified
Clenshaw-Curtis and Gauss-Kronrod formulas are employed. In addition to the
general strategy described for the IMSL routine QDAG, this routine uses an extrapolation
procedure known as the ɛ-algorithm. The routine QDAWC
is an implementation of the subroutine QAWC,
which is fully documented by Piessens et al. (1983).
1. Workspace may be explicitly provided, if desired, by use of Q2AWC/DQ2AWC. The reference is:
CALL Q2AWC (F, A, B, C, ERRABS, ERRREL, RESULT, ERREST, MAXSUB, NEVAL, NSUBIN, ALIST, BLIST, RLIST, ELIST, IORD)
The additional arguments are as follows:
MAXSUB —
Number of subintervals allowed. (Input)
A value of 500 is used
by QDAWC.
NEVAL — Number of evaluations of F. (Output)
NSUBIN — Number of subintervals generated. (Output)
ALIST — Array of length MAXSUB containing a list of the NSUBIN left endpoints. (Output)
BLIST — Array of length MAXSUB containing a list of the NSUBIN right endpoints. (Output)
RLIST — Array of length MAXSUB containing approximations to the NSUBIN integrals over the intervals defined by ALIST, BLIST. (Output)
ELIST — Array of length MAXSUB containing the error estimates of the NSUBIN values in RLIST. (Output)
IORD —
Array of length MAXSUB. (Output)
Let
k be NSUBIN if NSUBIN.LE.(MAXSUB/2 + 2),
MAXSUB + 1
− NSUBIN otherwise.
The first k locations contain pointers to the error estimates over the
subintervals, such that ELIST(IORD(1)), …, ELIST(IORD(k))
form a decreasing sequence.
2.
Informational
errors
Type
Code
4 1 The maximum number of subintervals allowed has been reached.
3 2 Roundoff error, preventing the requested tolerance from being achieved, has been detected.
3 3 A degradation in precision has been detected.
3.
If EXACT is the
exact value, QDAWC attempts to find
RESULT such that
ABS(EXACT − RESULT) .LE. MAX(ERRABS, ERRREL * ABS(EXACT)). To specify
only a relative error, set ERRABS to zero.
Similarly, to specify only an absolute error, set ERRREL to zero.
The Cauchy principal value of
is estimated. The values of the actual and estimated error are machine dependent.
USE QDAWC_INT
USE UMACH_INT
IMPLICIT NONE
INTEGER NOUT
REAL A, ABS, ALOG, B, C, ERRABS, ERREST, ERROR, EXACT, &
F, RESULT
INTRINSIC ABS, ALOG
EXTERNAL F
! Get output unit number
CALL UMACH (2, NOUT)
! Set limits of integration and C
A = -1.0
B = 5.0
C = 0.0
! Set error tolerances
ERRABS = 0.0
CALL QDAWC (F, A, B, C, RESULT, ERRABS=ERRABS, ERREST=ERREST)
! Print results
EXACT = ALOG(125./631.)/18.
ERROR = 2*ABS(RESULT-EXACT)
WRITE (NOUT,99999) RESULT, EXACT, ERREST, ERROR
99999 FORMAT (' Computed =', F8.3, 13X, ' Exact =', F8.3, /, /, &
' Error estimate =', 1PE10.3, 6X, 'Error =', 1PE10.3)
END
!
REAL FUNCTION F (X)
REAL X
F = 1.0/(5.*X**3+6.0)
RETURN
END
Computed =
-0.090
Exact = -0.090
Error estimate =
2.022E-06 Error = 2.980E-08
Visual Numerics, Inc. PHONE: 713.784.3131 FAX:713.781.9260 |