Chapter 4: Quadrature

.p>.CMCH4.DOC!INT_FCN_FOURIER;int_fcn_fourier

Computes a Fourier sine or cosine transform.

Synopsis

#include <imsl.h>

float imsl_f_int_fcn_fourier (float fcn(), float a, Imsl_quad weight, float omega, ¼, 0)

The type double function is imsl_d_int_fcn_fourier.

Required Arguments

float fcn (float x)   (Input)
User-supplied function to be integrated.

float a   (Input)
Lower limit of integration. The upper limit of integration is ¥.

 

Imsl_quad weight and float omega   (Input)
These two parameters are used to describe the trigonometric weight. The parameter weight can take on the two values described below, and the parameter omega = ¥ specifies the frequency of the trigonometric weighting function.

weight

Integration Weight

IMSL_COS

cos (wx)

IMSL_SIN

sin (wx)

Return Value

The return value is

if weight = IMSL_COS. If weight = IMSL_SIN, then the cosine factor is replaced with a sine factor. If no value can be computed, NaN is returned.

Synopsis with Optional Arguments

#include <imsl.h>

float imsl_f_int_fcn_fourier (float fcn(), float a, Imsl_quad weight, float omega,
IMSL_ERR_ABS, float err_abs,
IMSL_ERR_EST, float *err_est,
IMSL_MAX_SUBINTER, int max_subinter,
IMSL_MAX_CYCLES, int max_cycles,
IMSL_MAX_MOMENTS, int max_moments,
IMSL_N_CYCLES, int *n_cycles,
IMSL_N_EVALS, int *n_evals,
IMSL_FCN_W_DATA, float fcn(), void *data,
0)

Optional Arguments

IMSL_ERR_ABS, float err_abs   (Input)
Absolute accuracy desired.
Default:

where ɛ is the machine precision

IMSL_ERR_EST, float *err_est   (Output)
Address to store an estimate of the absolute value of the error.

IMSL_MAX_SUBINTER, int max_subinter   (Input)
Number of subintervals allowed.
Default: max_subinter = 500

IMSL_MAX_CYCLES, int max_cycles   (Input)
Number of cycles allowed.
Default: max_subinter = 50

IMSL_MAX_MOMENTS, int max_moments   (Input)
Number of subintervals allowed in the partition of each cycle.
Default: max_moments = 21

IMSL_N_CYCLES, int *n_cycles   (Output)
Address to store the number of cycles generated.

IMSL_N_EVALS, int *n_evals   (Output)
Address to store the number of evaluations of fcn.

IMSL_FCN_W_DATA, float fcn (float x, void *data), void *data (Input)
User supplied function to be integrated, which also accepts a pointer to data that is supplied by the user.  data is a pointer to the data to be passed to the user-supplied function.  See the Introduction, Passing Data to User-Supplied Functions at the beginning of this manual for more details.

Description

The function imsl_f_int_fcn_fourier is a special-purpose integrator that uses a globally adaptive scheme to reduce the absolute error. It 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 function imsl_f_int_fcn_trig followed by extrapolation.

The function imsl_f_int_fcn_fourier is based on the subroutine QAWF by Piessens et al. (1983).

Examples

Example 1

The value of

is computed. Notice that the integrand is coded to protect for the singularity at zero.

#include <math.h>
#include <imsl.h>

float           fcn(float x);
 
main()
{
    float       q, exact, omega;

    omega = imsl_f_constant("pi",0) / 2.;
                                /* Evaluate the integral */
    q = imsl_f_int_fcn_fourier (fcn, 0.0,
                                IMSL_COS, omega,
                                0);
                                /* Print the result and the */
                                /* exact answer */
    exact = 1.0;
    printf("integral  = %10.3f\nexact     = %10.3f\n", q, exact);
}

float fcn(float x)
{  
    return  (x==0.) ? 0. : 1./sqrt(x);
}

Output

integral  =      1.000
exact     =      1.000

Example 2

The value of

is again computed. The values of the actual and estimated error are printed as well. Note that these numbers are machine dependent. Furthermore, the error estimate is usually pessimistic. That is, the actual error is usually smaller than the error estimate,
as is the case in this example.The number of function evaluations also are printed. Notice that the integrand is coded to protect for the singularity at zero.

#include <math.h>
#include <imsl.h>

float           fcn(float x);
 
main()
{
    int         n_evals;
    float       q, exact, omega, err_est, exact_err;

    omega = imsl_f_constant("pi",0) / 2.0;
                                /* Evaluate the integral */
    q = imsl_f_int_fcn_fourier (fcn, 0.0,
                                IMSL_COS, omega,
                                IMSL_ERR_EST, &err_est,
                                IMSL_N_EVALS, &n_evals,
                                0);
                                /* Print the result and the */
                                /* exact answer */
    exact = 1.;
    exact_err = fabs(exact - q);
    printf("integral  = %10.3f\nexact     = %10.3f\n", q, exact);
    printf("error estimate   = %e\nexact error      = %e\n", err_est,
             exact_err);
    printf("The number of function evaluations  =  %d\n", n_evals);
}

float fcn(float x)
{  
    return  (x==0.) ? 0. : 1./sqrt(x);
}

Output

integral  =      1.000
exact     =      1.000
error estimate   = 1.803637e-04
exact error      = 1.013279e-06
The number of function evaluations  =  405

Warning Errors

IMSL_BAD_INTEGRAND_BEHAVIOR         Bad integrand behavior occurred in one or more cycles.

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

Fatal Errors

IMSL_MAX_CYCLES                                     Maximum number of cycles allowed has been reached.


Visual Numerics, Inc.
Visual Numerics - Developers of IMSL and PV-WAVE
http://www.vni.com/
PHONE: 713.784.3131
FAX:713.781.9260