dqawfe.f
SUBROUTINE DQAWFE (F, A, OMEGA, INTEGR, EPSABS, LIMLST, LIMIT,
+ MAXP1, RESULT, ABSERR, NEVAL, IER, RSLST, ERLST, IERLST, LST,
+ ALIST, BLIST, RLIST, ELIST, IORD, NNLOG, CHEBMO)
C***BEGIN PROLOGUE DQAWFE
C***PURPOSE The routine calculates an approximation result to a
C given Fourier integral
C I = Integral of F(X)*W(X) over (A,INFINITY)
C where W(X)=COS(OMEGA*X) or W(X)=SIN(OMEGA*X),
C hopefully satisfying following claim for accuracy
C ABS(I-RESULT).LE.EPSABS.
C***LIBRARY SLATEC (QUADPACK)
C***CATEGORY H2A3A1
C***TYPE DOUBLE PRECISION (QAWFE-S, DQAWFE-D)
C***KEYWORDS AUTOMATIC INTEGRATOR, CONVERGENCE ACCELERATION,
C FOURIER INTEGRALS, INTEGRATION BETWEEN ZEROS, QUADPACK,
C QUADRATURE, SPECIAL-PURPOSE INTEGRAL
C***AUTHOR Piessens, Robert
C Applied Mathematics and Programming Division
C K. U. Leuven
C de Doncker, Elise
C Applied Mathematics and Programming Division
C K. U. Leuven
C***DESCRIPTION
C
C Computation of Fourier integrals
C Standard fortran subroutine
C Double precision version
C
C PARAMETERS
C ON ENTRY
C F - Double precision
C Function subprogram defining the integrand
C Function F(X). The actual name for F needs to
C be declared E X T E R N A L in the driver program.
C
C A - Double precision
C Lower limit of integration
C
C OMEGA - Double precision
C Parameter in the WEIGHT function
C
C INTEGR - Integer
C Indicates which WEIGHT function is used
C INTEGR = 1 W(X) = COS(OMEGA*X)
C INTEGR = 2 W(X) = SIN(OMEGA*X)
C If INTEGR.NE.1.AND.INTEGR.NE.2, the routine will
C end with IER = 6.
C
C EPSABS - Double precision
C absolute accuracy requested, EPSABS.GT.0
C If EPSABS.LE.0, the routine will end with IER = 6.
C
C LIMLST - Integer
C LIMLST gives an upper bound on the number of
C cycles, LIMLST.GE.1.
C If LIMLST.LT.3, the routine will end with IER = 6.
C
C LIMIT - Integer
C Gives an upper bound on the number of subintervals
C allowed in the partition of each cycle, LIMIT.GE.1
C each cycle, LIMIT.GE.1.
C
C MAXP1 - Integer
C Gives an upper bound on the number of
C Chebyshev moments which can be stored, I.E.
C for the intervals of lengths ABS(B-A)*2**(-L),
C L=0,1, ..., MAXP1-2, MAXP1.GE.1
C
C ON RETURN
C RESULT - Double precision
C Approximation to the integral X
C
C ABSERR - Double precision
C Estimate of the modulus of the absolute error,
C which should equal or exceed ABS(I-RESULT)
C
C NEVAL - Integer
C Number of integrand evaluations
C
C IER - IER = 0 Normal and reliable termination of
C the routine. It is assumed that the
C requested accuracy has been achieved.
C IER.GT.0 Abnormal termination of the routine. The
C estimates for integral and error are less
C reliable. It is assumed that the requested
C accuracy has not been achieved.
C ERROR MESSAGES
C If OMEGA.NE.0
C IER = 1 Maximum number of cycles allowed
C Has been achieved., i.e. of subintervals
C (A+(K-1)C,A+KC) where
C C = (2*INT(ABS(OMEGA))+1)*PI/ABS(OMEGA),
C for K = 1, 2, ..., LST.
C One can allow more cycles by increasing
C the value of LIMLST (and taking the
C according dimension adjustments into
C account).
C Examine the array IWORK which contains
C the error flags on the cycles, in order to
C look for eventual local integration
C difficulties. If the position of a local
C difficulty can be determined (e.g.
C SINGULARITY, DISCONTINUITY within the
C interval) one will probably gain from
C splitting up the interval at this point
C and calling appropriate integrators on
C the subranges.
C = 4 The extrapolation table constructed for
C convergence acceleration of the series
C formed by the integral contributions over
C the cycles, does not converge to within
C the requested accuracy. As in the case of
C IER = 1, it is advised to examine the
C array IWORK which contains the error
C flags on the cycles.
C = 6 The input is invalid because
C (INTEGR.NE.1 AND INTEGR.NE.2) or
C EPSABS.LE.0 or LIMLST.LT.3.
C RESULT, ABSERR, NEVAL, LST are set
C to zero.
C = 7 Bad integrand behaviour occurs within one
C or more of the cycles. Location and type
C of the difficulty involved can be
C determined from the vector IERLST. Here
C LST is the number of cycles actually
C needed (see below).
C IERLST(K) = 1 The maximum number of
C subdivisions (= LIMIT) has
C been achieved on the K th
C cycle.
C = 2 Occurrence of roundoff error
C is detected and prevents the
C tolerance imposed on the
C K th cycle, from being
C achieved.
C = 3 Extremely bad integrand
C behaviour occurs at some
C points of the K th cycle.
C = 4 The integration procedure
C over the K th cycle does
C not converge (to within the
C required accuracy) due to
C roundoff in the
C extrapolation procedure
C invoked on this cycle. It
C is assumed that the result
C on this interval is the
C best which can be obtained.
C = 5 The integral over the K th
C cycle is probably divergent
C or slowly convergent. It
C must be noted that
C divergence can occur with
C any other value of
C IERLST(K).
C If OMEGA = 0 and INTEGR = 1,
C The integral is calculated by means of DQAGIE
C and IER = IERLST(1) (with meaning as described
C for IERLST(K), K = 1).
C
C RSLST - Double precision
C Vector of dimension at least LIMLST
C RSLST(K) contains the integral contribution
C over the interval (A+(K-1)C,A+KC) where
C C = (2*INT(ABS(OMEGA))+1)*PI/ABS(OMEGA),
C K = 1, 2, ..., LST.
C Note that, if OMEGA = 0, RSLST(1) contains
C the value of the integral over (A,INFINITY).
C
C ERLST - Double precision
C Vector of dimension at least LIMLST
C ERLST(K) contains the error estimate corresponding
C with RSLST(K).
C
C IERLST - Integer
C Vector of dimension at least LIMLST
C IERLST(K) contains the error flag corresponding
C with RSLST(K). For the meaning of the local error
C flags see description of output parameter IER.
C
C LST - Integer
C Number of subintervals needed for the integration
C If OMEGA = 0 then LST is set to 1.
C
C ALIST, BLIST, RLIST, ELIST - Double precision
C vector of dimension at least LIMIT,
C
C IORD, NNLOG - Integer
C Vector of dimension at least LIMIT, providing
C space for the quantities needed in the subdivision
C process of each cycle
C
C CHEBMO - Double precision
C Array of dimension at least (MAXP1,25), providing
C space for the Chebyshev moments needed within the
C cycles
C
C***REFERENCES (NONE)
C***ROUTINES CALLED D1MACH, DQAGIE, DQAWOE, DQELG
C***REVISION HISTORY (YYMMDD)
C 800101 DATE WRITTEN
C 890531 Changed all specific intrinsics to generic. (WRB)
C 890831 Modified array declarations. (WRB)
C 891009 Removed unreferenced variable. (WRB)
C 891009 REVISION DATE from Version 3.2
C 891214 Prologue converted to Version 4.0 format. (BAB)
C***END PROLOGUE DQAWFE