QDAWF

Computes a Fourier integral.

Required Arguments

F — User-supplied FUNCTION to be integrated. The form is

F(X), where

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)

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 infinity of F * WEIGHT. (Output)

Optional Arguments

ERRABS — Absolute accuracy desired. (Input)

Default: ERRABS = 1.e-3 for single precision and 1.d-8 for double precision.

Default: ERRABS = 1.e-3 for single precision and 1.d-8 for double precision.

ERREST — Estimate of the absolute value of the error. (Output)

Default: ERREST = 1.e-3 for single precision and 1.d-8 for double precision.

Default: ERREST = 1.e-3 for single precision and 1.d-8 for double precision.

FORTRAN 90 Interface

Generic: CALL QDAWF (F, A, IWEIGH, OMEGA, RESULT [, …])

Specific: The specific interface names are S_QDAWF and D_QDAWF.

FORTRAN 77 Interface

Single: CALL QDAWF (F, A, IWEIGH, OMEGA, ERRABS, RESULT, ERREST)

Double: The double precision name is DQDAWF.

Description

The routine QDAWF 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. The integration interval is always semi-infinite of the form [A, ∞]. These Fourier integrals are approximated by repeated calls to the IMSL routine QDAWO followed by extrapolation. The routine QDAWF is an implementation of the subroutine QAWF, which is fully documented by Piessens et al. (1983).

Comments

1. Workspace may be explicitly provided, if desired, by use of Q2AWF/DQ2AWF. The reference is:

CALL Q2AWF (F, A, IWEIGH, OMEGA, ERRABS, RESULT, ERREST, MAXCYL, MAXSUB, MAXCBY, NEVAL, NCYCLE, RSLIST, ERLIST, IERLST, NSUBIN, WK, IWK)

The additional arguments are as follows:

MAXSUB — Maximum number of subintervals allowed. (Input)

A value of 365 is used by QDAWF.

A value of 365 is used by QDAWF.

MAXCYL — Maximum number of cycles allowed. (Input)

MAXCYL must be at least 3. QDAWF uses 50.

MAXCYL must be at least 3. QDAWF uses 50.

MAXCBY — Maximum number of Chebyshev moments allowed. (Input)

QDAWF uses 21.

QDAWF uses 21.

NEVAL — Number of evaluations of F. (Output)

NCYCLE — Number of cycles used. (Output)

RSLIST — Array of length MAXCYL containing the contributions to the integral over the interval (A + (k − 1) * C, A + k * C), for k = 1, …, NCYCLE. (Output)

C = (2 * INT(ABS(OMEGA)) + 1) * PI/ABS(OMEGA).

C = (2 * INT(ABS(OMEGA)) + 1) * PI/ABS(OMEGA).

ERLIST — Array of length MAXCYL containing the error estimates for the intervals defined in RSLIST. (Output)

IERLST — Array of length MAXCYL containing error flags for the intervals defined in RSLIST. (Output)

IERLST(K) | Meaning |
---|---|

1 | The maximum number of subdivisions (MAXSUB) has been achieved on the k-th cycle. |

2 | Roundoff error prevents the desired accuracy from being achieved on the k-th cycle. |

3 | Extremely bad integrand behavior occurs at some points of the k-th cycle. |

4 | Integration procedure does not converge (to the desired accuracy) due to roundoff in the extrapolation procedure on the k-th cycle. It is assumed that the result on this interval is the best that can be obtained. |

5 | Integral over the k-th cycle is divergent or slowly convergent. |

NSUBIN — Number of subintervals generated. (Output)

WK — Work array of length 4 * MAXSUB + 25 * MAXCBY.

IWK — Work array of length 2 * MAXSUB.

2. Informational errors

Type | Code | Description |
---|---|---|

3 | 1 | Bad integrand behavior occurred in one or more cycles. |

4 | 2 | Maximum number of cycles allowed has been reached. |

3 | 3 | Extrapolation table constructed for convergence acceleration of the series formed by the integral contributions of the cycles does not converge to the requested accuracy. |

3. If EXACT is the exact value, QDAWF attempts to find RESULT such that ABS(EXACT − RESULT) .LE. ERRABS.

Example

The value of

is estimated. The values of the actual and estimated error are machine dependent. Notice that F is coded to protect for the singularity at zero.

USE QDAWF_INT

USE UMACH_INT

USE CONST_INT

IMPLICIT NONE

INTEGER IWEIGH, NOUT

REAL A, ABS, ERRABS, ERREST, ERROR, EXACT, F, &

OMEGA, PI, RESULT

INTRINSIC ABS

EXTERNAL F

! Get output unit number

CALL UMACH (2, NOUT)

! Set lower limit of integration

A = 0.0

! Select weight W(X) = COS(PI*X/2)

IWEIGH = 1

PI = CONST('PI')

OMEGA = PI/2.0

! Set error tolerance

CALL QDAWF (F, A, IWEIGH, OMEGA, RESULT, ERREST=ERREST)

! Print results

EXACT = 1.0

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 SQRT

INTRINSIC SQRT

IF (X .GT. 0.0) THEN

F = 1.0/SQRT(X)

ELSE

F = 0.0

END IF

RETURN

END

Output

Computed = 1.000 Exact = 1.000

Error estimate = 6.267E-04 Error = 2.205E-06