cquad is a new doubly-adaptive general-purpose quadrature
routine which can handle most types of singularities,
non-numerical function values such as
as well as some divergent integrals. It generally requires more
function evaluations than the integration routines in
quadpack, yet fails less often for difficult integrands.
The underlying algorithm uses a doubly-adaptive scheme in which Clenshaw-Curtis quadrature rules of increasing degree are used to compute the integral in each interval. The L_2-norm of the difference between the underlying interpolatory polynomials of two successive rules is used as an error estimate. The interval is subdivided if the difference between two successive rules is too large or a rule of maximum degree has been reached.
This function allocates a workspace sufficient to hold the data for n intervals. The number n is not the maximum number of intervals that will be evaluated. If the workspace is full, intervals with smaller error estimates will be discarded. A minimum of 3 intervals is required and or most functions, a workspace of size 100 is sufficient.
This function frees the memory associated with the workspace w.
This function computes the integral of f over (a,b) within the desired absolute and relative error limits, epsabs and epsrel using the cquad algorithm. The function returns the final approximation, result, an estimate of the absolute error, abserr, and the number of function evaluations required, nevals.
The cquad algorithm divides the integration region into subintervals, and in each iteration, the subinterval with the largest estimated error is processed. The algorithm uses Clenshaw-Curits quadrature rules of degree 4, 8, 16 and 32 over 5, 9, 17 and 33 nodes respectively. Each interval is initialized with the lowest-degree rule. When an interval is processed, the next-higher degree rule is evaluated and an error estimate is computed based on the L_2-norm of the difference between the underlying interpolating polynomials of both rules. If the highest-degree rule has already been used, or the interpolatory polynomials differ significantly, the interval is bisected.
The subintervals and their results are stored in the memory provided by workspace. If the error estimate or the number of function evaluations is not needed, the pointers abserr and nevals can be set to