Midpoint Rule Integration for Improper Integrals

This documentation was generated with the assistance of AI. Please report any inaccuracies.

The extended midpoint rule is an open quadrature formula — it never evaluates the integrand at the endpoints of the integration interval — which makes it the right building block for improper integrals: endpoints that cannot be evaluated exactly, infinite or semi-infinite ranges, and several standard kinds of endpoint singularity. It is described in Press et al., 2007, Sections 4.1 and 4.4 (pages 161, 167–172), and implemented as MidPointQuadrature (and MidPointMatrixQuadrature) plus four derived quadratures that each apply a different change of variables.

Panel (a) contrasts closed formulas, which evaluate the integrand at both endpoints, with open formulas such as the midpoint rule, which use only interior points and so never need to evaluate the integrand exactly at a and b

Panel (a) above is the key distinction: because the midpoint rule never samples at or , it keeps working even when the integrand is undefined, infinite, or singular exactly at an endpoint — the situations this page’s five quadrature classes are built to handle.

What makes an integral "improper"

This package’s midpoint-based quadratures handle four situations:

  • The integrand approaches a finite limit at a finite endpoint but cannot actually be evaluated there (e.g. at ).

  • The upper limit is , or the lower limit is .

  • The integrand has an integrable singularity at a known endpoint (e.g. at ).

  • The integrand decays exponentially toward .

A fifth case — a singularity at an unknown location strictly between the limits — is out of scope for these quadratures; it needs an adaptive or variable-stepsize method instead. An integral that is actually infinite or has no limiting value at all (e.g. or ) is not improper but impossible — no substitution can fix it.

The extended midpoint rule

Like the extended trapezoidal rule, this formula’s error series (the second Euler-Maclaurin summation formula) contains only even powers of , which is what makes it a suitable base rule for Simpson’s Rule Integration and Romberg Integration as well. Unlike the trapezoidal rule, however, its evaluation points cannot be reused if the number of intervals is only doubled — tripling is required instead. This costs a modest, bounded amount of unnecessary work (a factor of about versus for doubling) in exchange for being an open formula:

  • Stage 1:

  • Stage : adds new interior points, alternating spacing between the new points, and combines them with the previous estimate.

Handling improper integrals by change of variables

Each situation above is handled by substituting a new variable that removes the problem, then applying the midpoint rule (still an open formula) to the transformed integrand. This package implements the substitution as a MidPointQuadrature subclass that overrides only the point at which the (identity) mapping is applied, reusing the same refinement loop:

flowchart TD Q["Quadrature (abstract): next()"] --> M["MidPointQuadrature: extended midpoint\nrefinement loop, calls func(x)"] M -->|"func(x) = f(x), identity"| M1["MidPointQuadrature itself\n(open formula, no singularity)"] M -->|"func(x) = f(1/x)/x^2"| M2["InfinityMidPointQuadrature\n(semi-infinite / infinite range)"] M -->|"func(x) = f(-log x)/x"| M3["ExponentialMidPointQuadrature\n(exponential decay to infinity)"] M -->|"func(x) = 2x f(a+x^2)"| M4["LowerSquareRootMidPointQuadrature\n(inverse sqrt singularity at a)"] M -->|"func(x) = 2x f(b-x^2)"| M5["UpperSquareRootMidPointQuadrature\n(inverse sqrt singularity at b)"]
Situation Quadrature class Substitution

Endpoint unreachable, integrand otherwise well-behaved

MidPointQuadrature

None — the open formula simply avoids the endpoints.

Upper limit (or lower limit ), any sign, integrand decays faster than

InfinityMidPointQuadrature

, mapping to .

Integrand at the lower bound

LowerSquareRootMidPointQuadrature

, mapping to .

Integrand at the upper bound

UpperSquareRootMidPointQuadrature

, mapping to .

Upper limit , integrand decays exponentially

ExponentialMidPointQuadrature

, i.e. , mapping to .

The square-root substitutions are the special, frequently occurring case of the more general power-law substitution used for an integrand diverging as , ; the (rather than a general power of ) factor that appears in func(x) is exactly what cancels the singularity, leaving a smooth transformed integrand for the midpoint rule to integrate.

If a lower bound is with an upper bound that is finite (or vice versa for a doubly-infinite range), split the integral at a convenient interior breakpoint and integrate each piece with the quadrature that matches it — e.g. MidPointQuadrature on the finite piece and InfinityMidPointQuadrature on the infinite one.

Library implementation

MidPointQuadrature implements the refinement loop and calls a protected double func(double x) hook — identity by default — that each subclass overrides to apply its substitution; the tripling refinement loop itself is inherited unchanged. MidPointMatrixQuadrature and its four substitution subclasses mirror this for matrix-valued integrands.

Every one of these quadrature rules can be paired with any of the package’s three integration methods — plain (IntegratorType.QUADRATURE), Simpson’s Rule Integration (IntegratorType.SIMPSON), or Romberg Integration (IntegratorType.ROMBERG) — through the corresponding QuadratureType value (MID_POINT, INFINITY_MID_POINT, LOWER_SQUARE_ROOT_MID_POINT, UPPER_SQUARE_ROOT_MID_POINT, EXPONENTIAL_MID_POINT) and the package’s Integrator.create(…​) factory:

final var integrator = Integrator.create(
        0.0, 1.0, listener, IntegratorType.ROMBERG, QuadratureType.LOWER_SQUARE_ROOT_MID_POINT);
final var result = integrator.integrate();
QuadratureType.EXPONENTIAL_MID_POINT always treats the upper limit as infinite — the corresponding constructors (e.g. ExponentialMidPointQuadratureIntegrator, RombergExponentialMidPointQuadratureIntegrator) take only a lower limit, and Integrator.create(…​) ignores whatever upper limit is passed when this quadrature type is selected.

Example

Integrating from to — an inverse square-root singularity at the lower bound:

final var listener = new SingleDimensionFunctionEvaluatorListener() {
    @Override
    public double evaluate(final double point) {
        return 1.0 / Math.sqrt(point);
    }
};

final var integrator = new LowerSquareRootMidPointQuadratureIntegrator(0.0, 1.0, listener);
final var result = integrator.integrate(); // approximately 2.0