Chapter 4: Quadrature > int_fcn_sing

int_fcn_sing

Integrates a function, which may have endpoint singularities, using a globally adaptive scheme based on Gauss-Kronrod rules.

Synopsis

#include <imsl.h>

float imsl_f_int_fcn_sing (float fcn(), float a, float b,, 0)

The type double function is imsl_d_int_fcn_sing.

Required Arguments           

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

float a   (Input)
Lower limit of integration.

float b   (Input)
Upper limit of integration.

Return Value

An estimate of

If no value can be computed, NaN is returned.

Synopsis with Optional Arguments

#include <imsl.h>

float imsl_f_int_fcn_sing (float fcn(), float a, float b,

IMSL_ERR_ABS, float err_abs,

IMSL_ERR_REL, float err_rel,

IMSL_ERR_EST, float *err_est,

IMSL_MAX_SUBINTER, int max_subinter,

IMSL_N_SUBINTER, int *n_subinter,

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_REL, float err_rel   (Input)
Relative 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_N_SUBINTER, int *n_subinter   (Output)
Address to store the number of subintervals 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

This function is designed to handle functions with endpoint singularities. However, the performance on functions that are well-behaved at the endpoints is also quite good.

The function imsl_f_int_fcn_sing is a general-purpose integrator that uses a globally adaptive scheme in order to reduce the absolute error. It subdivides the interval [a, b] and uses a 21-point Gauss-Kronrod rule to estimate the integral over each subinterval. The error for each subinterval is estimated by comparison with the 10-point Gauss quadrature rule. The subinterval with the largest estimated error is then bisected, and the same procedure is applied to both halves. The bisection process is continued until either the error criterion is satisfied, roundoff error is detected, the subintervals become too small, or the maximum number of subintervals allowed is reached. This function uses an extrapolation procedure known as the ɛ-algorithm.

On some platforms, imsl_f_int_fcn_sing can evaluate the user-supplied function fcn in parallel. This is done only if the function imsl_omp_options is called to flag user-defined functions as thread-safe. A function is thread-safe if there are no dependencies between calls. Such dependencies are usually the result of writing to global or static variables.

The function imsl_f_int_fcn_sing is based on the subroutine QAGS by Piessens et al. (1983).

Examples

Example 1

The value of

is estimated.

#include <math.h>

#include <imsl.h>

 

float           fcn(float x);

 

int main()

{

    float       q, exact;

 

    imsl_omp_options(IMSL_SET_FUNCTIONS_THREAD_SAFE, 1, 0);

                                /* Evaluate the integral */

    q = imsl_f_int_fcn_sing (fcn, 0.0, 1.0, 0);

                                /* Print the result and */

                                /*the exact answer */

    exact = -4.0;

    printf("integral  = %10.3f\nexact     = %10.3f\n", q, exact);

}

 

float fcn(float x)

{  

    return log(x)/sqrt(x);

}

Output

integral  =     -4.000

exact     =     -4.000

Example 2

The value of

is again estimated. The values of the actual and estimated errors are printed as well. Note that these numbers are machine dependent. Furthermore, usually the error estimate is pessimistic. That is, the actual error is usually smaller than the error estimate as is in this example.

#include <math.h>

#include <imsl.h>

 

float           fcn(float x);

 

int main()

{

    float       q, exact, err_est, exact_err;

 

    imsl_omp_options(IMSL_SET_FUNCTIONS_THREAD_SAFE, 1, 0);

                                /* Evaluate the integral */

    q = imsl_f_int_fcn_sing (fcn, 0.0, 1.0,

                             IMSL_ERR_EST, &err_est,

                             0);

                                /* Print the result and */

                                /* the exact answer */

    exact = -4.0;

    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);

}

 

float fcn(float x)

{  

    return log(x)/sqrt(x);

}

Output

integral  =     -4.000

exact     =     -4.000

error estimate   = 2.708435e-004

exact error      = 2.241135e-005

Warning Errors

 

IMSL_ROUNDOFF_CONTAMINATION

Roundoff error, preventing the requested tolerance from being achieved, has been detected.

IMSL_PRECISION_DEGRADATION

A degradation in precision has been detected.

IMSL_EXTRAPOLATION_ROUNDOFF

Roundoff error in the extrapolation table, preventing the requested tolerance from being achieved, has been detected.

Fatal Errors

 

IMSL_DIVERGENT

Integral is probably divergent or slowly convergent.

IMSL_PRECISION_DEGRADATION

Integral is probably divergent or slowly convergent.

IMSL_MAX_SUBINTERVALS

The maximum number of subintervals allowed has been reached.


RW_logo.jpg
Contact Support