How quickly does a numerical integral converge — and when does it slow down?

For $\int_0^1x^2\,dx=1/3$, four trapezoids give $11/32$, an error of exactly $1/96$. Each doubling of $n$ divides the trapezoid error by exactly $4$ and the left-rule error by about $2$. Now choose √x: slower: the trapezoid error shrinks by only about $2^{1.5}$ per doubling, because $\sqrt x$ has an infinite slope at $0$.

Intuition: how much does each refinement buy?

Every rule replaces $f$ on each panel by something simpler — a constant for rectangle rules, a straight line for the trapezoid rule. Halving the panel width $h$ makes the replacement better. How much better depends on how closely the simple shape follows the curve, which depends on the rule and on how smooth $f$ is. Plotting $\log(\text{error})$ against the refinement level turns a rate $h^p$ into a straight line with slope $p$.

Definition

With $n$ panels of width $h=(b-a)/n$ and points $x_i=a+ih$:

and $R_n$ uses right endpoints. A rule has order $p$ for a class of functions if its error behaves like $Ch^p$. The observed order from two refinements is $p\approx\log_2(E_{n}/E_{2n})$.

Mechanics: what the explorer computes

For $x^2$, $x^3-x$ and $|x-1/3|$ every approximation and every error is an exact rational, computed for $n=1,2,4,\dots,1024$. For $\sin x$ and $\sqrt x$ the values are binary64, with references $1-\cos3$ and $2/3$ from their antiderivatives. Panels are drawn while they are wider than a pixel. The convergence plot and table come from the same computed errors.

Derivation: why the trapezoid rule has order 2

On one panel $[x_i,x_i+h]$, the trapezoid replaces $f$ by its linear interpolant. The interpolation error is $\tfrac12f''(\xi)(t-x_i)(t-x_i-h)$, and integrating over the panel gives $-\tfrac{h^3}{12}f''(\xi_i)$. Summing $n=(b-a)/h$ panels,

for some $\xi$, provided $f''$ is continuous. So halving $h$ quarters the error. A rectangle rule replaces $f$ by a constant and loses a factor $h$: order $1$. The midpoint rule is order $2$ with half the trapezoid’s constant and the opposite sign.

Worked example

$\int_0^1x^2\,dx=1/3$. With $n$ trapezoids the error is exactly $1/(6n^2)$: $n=1$ gives $1/2$ (error $1/6$), $n=2$ gives $3/8$ (error $1/24$), $n=4$ gives $11/32$ (error $1/96$). The midpoint error is exactly $-1/(12n^2)$, and the left-rule error is $-1/(2n)+1/(6n^2)$, which approaches order $1$.

When the rate changes

• Singular slope. $\sqrt x$ has $f'\to\infty$ at $0$; the $h^2$ bound does not apply and the trapezoid error behaves like $h^{1.5}$.
• A corner off the grid. $|x-1/3|$ has a corner that never falls on a grid point $k/2^m$; the first refinements behave irregularly before settling.
• Signed integrands. For $x^3-x$ on $[-1,2]$, positive and negative parts partly cancel in the exact value $9/4$; the rules must resolve both.
• No universal promise. Highly oscillatory or discontinuous integrands can make any fixed rule converge slowly or appear to converge to the wrong value for small $n$.

Executable reference code

The reference computes exact errors of all four rules with Fraction and recovers the observed orders $1$ and $2$ for $x^2$. The NumPy version measures observed orders for $\sin x$ and $\sqrt x$ and shows the drop from $2$ to about $1.5$.

"""Exact errors of four quadrature rules for x**2 on [0, 1]; observed orders 1 and 2 (III.15)."""
from fractions import Fraction as F
import math


def rule(name, f, a, b, n):
    h = (b - a) / n
    x = [a + i * h for i in range(n + 1)]
    if name == "left":
        return h * sum(f(t) for t in x[:-1])
    if name == "right":
        return h * sum(f(t) for t in x[1:])
    if name == "midpoint":
        return h * sum(f((s + t) / 2) for s, t in zip(x, x[1:]))
    return h * (f(x[0]) / 2 + sum(f(t) for t in x[1:-1]) + f(x[-1]) / 2)


f = lambda t: t * t
exact = F(1, 3)
for n in (1, 2, 4, 8, 16):
    assert rule("trapezoid", f, F(0), F(1), n) - exact == F(1, 6 * n * n)
    assert rule("midpoint", f, F(0), F(1), n) - exact == F(-1, 12 * n * n)
    assert rule("left", f, F(0), F(1), n) - exact == F(-1, 2 * n) + F(1, 6 * n * n)

order = lambda e1, e2: math.log2(abs(e1) / abs(e2))
e = [rule("left", f, F(0), F(1), 2**k) - exact for k in range(8, 10)]
assert abs(order(*e) - 1) < 0.01

"""Observed trapezoid orders: about 2 for sin on [0, 3], about 1.5 for sqrt on [0, 1] (III.15)."""
import numpy as np

def trapezoid(f, a, b, n):
    x = np.linspace(a, b, n + 1)
    return np.trapezoid(f(x), x)

def observed(f, a, b, exact, k):
    e1, e2 = abs(trapezoid(f, a, b, 2**k) - exact), abs(trapezoid(f, a, b, 2**(k + 1)) - exact)
    return np.log2(e1 / e2)

assert abs(observed(np.sin, 0, 3, 1 - np.cos(3), 8) - 2) < 0.05
assert abs(observed(np.sqrt, 0, 1, 2 / 3, 10) - 1.5) < 0.05

Explore

Compare the trapezoid and midpoint curves: parallel lines, the midpoint one lower by a factor $2$. On the signed integrand, switch to the left rule and watch the order drop to $1$. On sin x, push $k$ to $10$: the error approaches $10^{-6}$, still far above rounding level.

Common mistakes

Assuming more panels always gives the promised rate regardless of smoothness. Reading the error from a single $n$ instead of its behaviour under refinement. Confusing the number of panels with the number of function evaluations (the trapezoid rule uses $n+1$ points).

Numerical considerations

For smooth integrands, higher-order rules (Simpson’s rule, Gaussian quadrature) reach small errors with far fewer evaluations; adaptive methods refine only where the integrand needs it. At very large $n$ rounding in the sum eventually matters. Reliable software reports an error estimate, usually from comparing two refinements, exactly as the observed-order column does.

Engineering and ML connections

Numerical integration computes expectations, normalising constants, areas under curves and solutions of differential equations (one step of many ODE methods is a quadrature rule). Monte Carlo integration (Part IV) trades the deterministic rates here for a rate of $n^{-1/2}$ that does not depend on dimension.

Exercises

Summary and checkpoints

Rectangle rules converge with order 1, trapezoid and midpoint rules with order 2 — for smooth integrands. Singularities and corners lower the observed order; signed integrands require resolving cancellation. Measure convergence by refinement rather than trusting a single approximation. This completes Part III: the next Part turns to probability, where integrals of densities measure chance.