Integrates a function containing a sine or a cosine.
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)
IWEIGH — Type of weight function used. (Input)
IWEIGH Weight
1 COS(OMEGA * X)
2 SIN(OMEGA * X)
OMEGA — Parameter in the weight function. (Input)
RESULT — Estimate of the integral from A to B of F * WEIGHT. (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: ERRREL = 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 QDAWO (F, A, B, IWEIGH, OMEGA, RESULT [,…])
Specific: The specific interface names are S_QDAWO and D_QDAWO.
Single: CALL QDAWO (F, A, B, IWEIGH, OMEGA, ERRABS, ERRREL, RESULT, ERREST)
Double: The double precision name is DQDAWO.
The routine QDAWO 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) is either cos ωx or sin ωx. Depending on the length of the subinterval in relation to the size of ω, either a modified Clenshaw-Curtis procedure or a Gauss-Kronrod 7/15 rule is employed to approximate the integral on a subinterval. In addition to the general strategy described for the IMSL routine QDAG, this subroutine uses an extrapolation procedure known as the ɛ-algorithm. The routine QDAWO is an implementation of the subroutine QAWO, which is fully documented by Piessens et al. (1983).
1. Workspace may be explicitly provided, if desired, by use of Q2AWO/DQ2AWO. The reference is:
CALL Q2AWO (F, A, B, IWEIGH, OMEGA, ERRABS, ERRREL, RESULT, ERREST, MAXSUB, MAXCBY, NEVAL, NSUBIN, ALIST, BLIST, RLIST, ELIST, IORD, NNLOG, WK)
The additional arguments are as follows:
MAXSUB —
Maximum number of subintervals allowed. (Input)
A value of 390
is used by QDAWO.
MAXCBY — Upper bound on the number of Chebyshev moments which can be stored. That is, for the intervals of lengths ABS(B − A) * 2**(−L), L = 0, 1, …, MAXCBY − 2, MAXCBY.GE.1. The routine QDAWO uses 21. (Input)
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. 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. (Output)
NNLOG — Array of length MAXSUB containing the subdivision levels of the subintervals, i.e. NNLOG(I) = L means that the subinterval numbered I is of length ABS(B − A) * (1− L). (Output)
WK — Array of length 25 * MAXCBY. (Workspace)
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 4 Roundoff error in the extrapolation table, preventing the requested tolerances from being achieved, has been detected.
4 5 Integral is probably divergent or slowly convergent.
3.
If EXACT is the
exact value, QDAWO 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 value of
is estimated. The values of the actual and estimated error are machine dependent. Notice that the log function is coded to protect for the singularity at zero.
USE QDAWO_INT
USE UMACH_INT
USE CONST_INT
IMPLICIT NONE
INTEGER IWEIGH, NOUT
REAL A, ABS, B, ERRABS, ERREST, ERROR, &
EXACT, F, OMEGA, PI, RESULT
INTRINSIC ABS
EXTERNAL F
! Get output unit number
CALL UMACH (2, NOUT)
! Set limits of integration
A = 0.0
B = 1.0
! Weight function = sin(10.*pi*x)
IWEIGH = 2
PI = CONST('PI')
OMEGA = 10.*PI
! Set error tolerances
ERRABS = 0.0
CALL QDAWO (F, A, B, IWEIGH, OMEGA, RESULT, ERRABS=ERRABS, &
ERREST=ERREST)
! Print results
EXACT = -0.1281316
ERROR = 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
REAL ALOG
INTRINSIC ALOG
IF (X .EQ. 0.) THEN
F = 0.0
ELSE
F = ALOG(X)
END IF
RETURN
END
Computed =
-0.128
Exact = -0.128
Error estimate =
7.504E-05 Error = 5.260E-06
Visual Numerics, Inc. PHONE: 713.784.3131 FAX:713.781.9260 |