The AdaptiveIntegrator class provides the most
robust, reliable, and in most cases the fastest numerical integration algorithm. It uses one of the
GaussKronrodIntegrator derived integration rules to approximate the integral over subintervals. In each
step, the interval with the largest estimated error is divided in two, and new approximations for these two
subintervals are calculated. This process continues until the total estimated error over all intervals is within the
desired tolerance.
The IntegrationRule property
lets you specify which integration rule to use to approximate integrals on a subinterval. The default is the 21-point
Gauss-Kronrod rule for normal integrands, and the 15-point Gauss-Kronrod rule when singularities are expected. For
oscillating target functions, the higher order rules will tend to give better results. This value is not limited to
Gauss-Kronrod integration rules. Any integration rule that provides an estimate for the error can be used.
AdaptiveIntegrator integrator = new AdaptiveIntegrator();
integrator.RelativeTolerance = 1e-5;
integrator.ConvergenceCriterion = ConvergenceCriterion.WithinRelativeTolerance;
double result = integrator.Integrate(Math.Sin, 0, 2);
Console.WriteLine("sin(x) on [0,2]");
Console.WriteLine("Romberg integrator:");
Console.WriteLine(" Value: {0}", result);
Console.WriteLine(" Status: {0}", integrator.Status);
Console.WriteLine(" Estimated error: {0}", integrator.EstimatedError);
Console.WriteLine(" Iterations: {0}", integrator.IterationsNeeded);
Console.WriteLine(" Function evaluations: {0}", integrator.EvaluationsNeeded);
Console.WriteLine(" Order: {0}", integrator.Order);Several enhancements make this the integration algorithm of choice, especially for 'difficult' integrands.
An extra property, UseExtrapolation, lets you specify whether
special precautions should be taken to deal with singularities. If UseExtrapolation is
true, an extrapolation algorithm is used to approximate the integral near a singularity. Because this
option requires extra resources, it is recommended only for target functions that actually have singularities inside
or at the boundaries of the integration interval.
The SetSingularities method
lets you specify the points inside the integration interval where the integrand has a singularity or a discontinuity.
The value of this property is a Double array containing the problem points. Extrapolation is used once again to
approximate the integral near the singularity or singularities. Specifying the problem points in the target function
can lead to speed-ups of several orders of magnitude.
Func<double, double> f = x => x * x * x
* Math.Log(Math.Abs((x * x - 1) * (x * x - 2)));
integrator.AbsoluteTolerance = 1e-8;
integrator.ConvergenceCriterion = ConvergenceCriterion.WithinAbsoluteTolerance;
integrator.SetSingularities(1, Math.Sqrt(2));
result = integrator.Integrate(f, 0, 3);
Console.WriteLine("Singularities inside the integration interval");
Console.WriteLine(" Value: {0}", integrator.Result);
Console.WriteLine(" Exact value: 52.740748383471444998");
Console.WriteLine(" Status: {0}", integrator.Status);
Console.WriteLine(" Estimated error: {0}", integrator.EstimatedError);
Console.WriteLine(" Iterations: {0}", integrator.IterationsNeeded);
Console.WriteLine(" Function evaluations: {0}", integrator.EvaluationsNeeded);
If the only discontinuities are on a boundary of the integration interval, you do not need to specify the
singularities. You have to make sure, however, that UseExtrapolation is set to true.
Func<double, double> f2 = x => Math.Pow(x, -0.9) * Math.Log(1 / x);
integrator.AbsoluteTolerance = 1e-8;
integrator.ConvergenceCriterion = ConvergenceCriterion.WithinAbsoluteTolerance;
integrator.UseExtrapolation = true;
result = integrator.Integrate(f2, 0, 1);
Console.WriteLine("Integrand with a singularity on a boundary:");
Console.WriteLine(" Value: {0}", integrator.Result);
Console.WriteLine(" Exact value: 100");
Console.WriteLine(" Status: {0}", integrator.Status);
Console.WriteLine(" Estimated error: {0}", integrator.EstimatedError);
Console.WriteLine(" Iterations: {0}", integrator.IterationsNeeded);
Console.WriteLine(" Function evaluations: {0}", integrator.EvaluationsNeeded);
Unlike all methods mentioned previously, this integrator can deal successfully with infinite integration
intervals. Specify Double.NegativeInfinity or Double.PositiveInfinity for one or both
bounds of the integration interval. Internally, the class transforms the infinite interval to a finite interval, and
performs the integration of the transformed target function over the transformed interval.
Func<double, double> f3 = x => Math.Exp(-x - x * x);
integrator.AbsoluteTolerance = 1e-8;
integrator.ConvergenceCriterion = ConvergenceCriterion.WithinAbsoluteTolerance;
result = integrator.Integrate(f3,
Double.NegativeInfinity, Double.PositiveInfinity);
Console.WriteLine("Exp(-x^2-x) on [-inf,inf]");
Console.WriteLine(" Value: {0}", integrator.Result);
Console.WriteLine(" Status: {0}", integrator.Status);
Console.WriteLine(" Estimated error: {0}", integrator.EstimatedError);
Console.WriteLine(" Iterations: {0}", integrator.IterationsNeeded);
Console.WriteLine(" Function evaluations: {0}", integrator.EvaluationsNeeded);
R. Piessens, E. de Doncker, C. Uberhuber and D. Kahaner,
QUADPACK, A Subroutine Package for Automatic
Integration
, Springer-Verlag, New York, 1983.
S.D. Conte and Carl de Boor, Elementary Numerical Analysis: An Algorithmic Approach, McGraw-Hill,
1972.