Implementing numerical integration from scratch — trapezoidal, Simpson's: Worked examples

Before anything, we recall the two machines we are feeding inputs into. Both take a function $f$ , an interval $[a, b]$ , and a number of equal slices $n$ . The slice width is the node spacing $h = \frac{b - a}{n}$ , and the points where we sample the curve are the nodes $x_{i} = a + ih$ .

The error terms — where the signs come from

Before the examples we settle the two error formulas once, so every "Verify:" below can lean on them without hand-waving. We define the error as rule minus truth, $E = (numerical estimate) - \int_{a}^{b} f d x$ . A positive $E$ means the rule overshoots.

Where $+ \frac{1}{12} h^{3} f^{''}$ actually comes from (a Taylor sketch)

We do not want you to trust the local error $+ \frac{1}{12} h^{3} f^{''} (ξ)$ ; here is the argument in plain steps. Put one slice on $[- \frac{h}{2}, \frac{h}{2}]$ (width $h$ , centred at $0$ for symmetry), and let $m$ be the midpoint value.

Step A — expand the curve. A Taylor series says that near the midpoint, $f (x) = f (0) + f^{'} (0) x + \frac{1}{2} f^{''} (0) x^{2} + \dots$ Why? Taylor is the one tool that writes an unknown curve as "value + tilt + bend + …", exactly the pieces a straight chord can and cannot match.

Step B — integrate the true area. Over the symmetric slice the odd term $f^{'} (0) x$ integrates to $0$ , leaving $\int_{- h /2}^{h /2} f d x = f (0) h + \frac{1}{24} f^{''} (0) h^{3} + \dots$ Why? $\int_{- h /2}^{h /2} x^{2} d x = \frac{h ^{3}}{12}$ , and $\frac{1}{2} \cdot \frac{h ^{3}}{12} = \frac{h ^{3}}{24}$ .

Step C — integrate the chord (what trapezoid actually computes). The trapezoid uses the two endpoint heights $f (\pm h /2)$ , whose average is $\frac{1}{2} [f (- \frac{h}{2}) + f (\frac{h}{2})] = f (0) + \frac{1}{8} f^{''} (0) h^{2} + \dots$ Multiplying by the width $h$ gives the trapezoid area $f (0) h + \frac{1}{8} f^{''} (0) h^{3} + \dots$ Why? The chord matches value and tilt but over-accounts the bend, because it samples the curve only at the two far ends where a convex curve is highest.

Step D — subtract. Local error = chord − truth: $(\frac{1}{8} - \frac{1}{24}) f^{''} (0) h^{3} = \frac{2}{24} f^{''} (0) h^{3} = + \frac{1}{12} h^{3} f^{''} (ξ) .$ Why the plus? $\frac{1}{8} > \frac{1}{24}$ , so the chord's bend-term wins → positive local error, matching the convex-overshoot picture.

Now sum the $n = \frac{b - a}{h}$ slices: each contributes $\frac{1}{12} h^{3} f^{''}$ , and $n$ of them turns one factor of $h$ into $(b - a)$ :

The scenario matrix

Every cell below is a class of situation the topic can throw at you. The examples that follow are each tagged with the cell they cover.

#	Case class	What is special about it	Example
C1	Exact case — polynomial within the rule's degree	Error should be exactly zero (up to rounding)	Ex 1
C2	Convex curve, overshoot sign	Trapezoid over-estimates; which sign is the error?	Ex 2
C3	Concave curve, undershoot sign	Trapezoid under-estimates — the opposite sign	Ex 3
C4	Rule fails silently — degree above what it captures	Simpson is not exact; a residual error remains	Ex 4
C5	Degenerate / reversed input — $b = a$ , $f \equiv 0$ , $n = 1$ , or $a > b$	Zero-width, flat, single slice, reversed limits	Ex 5
C6	Odd $n$ trap for Simpson	Parity bug: what to do when $n$ is odd	Ex 6
C7	Convergence / limiting behaviour — halving $h$	Error must drop by $4 \times$ (trap) and $16 \times$ (Simpson)	Ex 7
C8	Real-world word problem — data only, no formula	You have measurements, not $f (x)$	Ex 8
C9	Exam twist — sign change / negative area under axis	Integral can be negative; "area" intuition misleads	Ex 9

Example 1 — the exact case (C1)

Worked example

$\int_{0}^{2} (3 x + 1) d x$ with the trapezoidal rule, $n = 1$ . Forecast: Guess before computing — will trapezoid be a little off, or dead exact?

Step 1. True value first. $\int_{0}^{2} (3 x + 1) d x = [\frac{3}{2} x^{2} + x]_{0}^{2} = 6 + 2 = 8$ . Why this step? We need the truth to compare against, and this integrand is easy.

Step 2. Set up trapezoid: $h = \frac{2 - 0}{1} = 2$ , nodes $x_{0} = 0, x_{1} = 2$ , so $f (0) = 1, f (2) = 7$ . Why this step? $n = 1$ means one single trapezoid — no interior nodes at all.

Step 3. Apply $T_{1} = \frac{h}{2} [f (x_{0}) + f (x_{1})] = \frac{2}{2} [1 + 7] = 8$ . Why this step? Straight-line integrand → the trapezoid's slanted roof is the curve, so no gap.

Verify: $T_{1} = 8 =$ true value. Trapezoid is exact for degree-1, and $3 x + 1$ is degree-1, so error must be zero — confirmed. Consistent with $E_{T} = \frac{( b - a ) h ^{2}}{12} f^{''}$ : here $f^{''} = 0$ for a line, so $E_{T} = 0$ . (Units: if $x$ were metres and $f$ metres, area is m², dimensionally sound.)

Recall Why did

$n = 1$ suffice? More slices would still give exactly 8 — refining a line changes nothing ::: because a straight line joining any two of its own points lies on the line, so every trapezoid is exact regardless of how many you use.

Example 2 — convex curve, trapezoid overshoots (C2)

Look at figure s01. The horizontal axis is $x$ on $[0, 1]$ , the vertical axis is $f (x)$ . The lavender curve is $x^{2}$ ; the green line is the two chords the trapezoids use. Takeaway: on a convex curve the chord lies above the curve, and the coral shaded slivers are the extra area the trapezoid wrongly counts — this is the overshoot made visible.

Figure — Implementing numerical integration from scratch — trapezoidal, Simpson's

Step 1. $h = \frac{1 - 0}{2} = 0.5$ ; nodes $0, 0.5, 1$ ; $f = 0, 0.25, 1$ . Why this step? Standard node table before plugging in.

Step 2. $T_{2} = \frac{0.5}{2} [0 + 1 + 2 (0.25)] = 0.25 \cdot 1.5 = 0.375$ . Why this step? Endpoints weight 1, the single interior node weights 2.

Step 3. Compare: $0.375 > 0.3333 \dots$ → overshoot by $+ 0.0417$ . Why this step? The sign of the error is the whole point of this case.

Verify: For any convex $f$ ( $f^{''} > 0$ ) the chord lies above the curve, so $T_{n} \geq \int$ . Using our fixed convention $E_{T} = + \frac{( b - a ) h ^{2}}{12} f^{''} (ξ)$ with $f^{''} = 2 > 0$ , the plus coefficient makes $E_{T} > 0$ — rule above truth. That matches the overshoot $+ 0.0417$ exactly in sign.

Example 3 — concave curve, trapezoid undershoots (C3)

Look at figure s02. Axes are the same ( $x$ horizontal on $[0, π]$ , $f (x)$ vertical). The lavender curve is $sin x$ ; the green line is the trapezoid chords. Takeaway: the chord now sags below the arch, and the coral slivers are the area the trapezoid misses — the exact mirror image of s01, which is why the sign of the error flips.

Step 1. $h = \frac{π}{2} \approx 1.5708$ ; nodes $0, \frac{π}{2}, π$ ; $f = 0, 1, 0$ . Why this step? Sample the arch at its base, peak, base.

Step 2. $T_{2} = \frac{π /2}{2} [0 + 0 + 2 (1)] = \frac{π}{4} \cdot 2 = \frac{π}{2} \approx 1.5708$ . Why this step? Only the peak node contributes; endpoints are zero.

Step 3. $1.5708 < 2$ → undershoot by $- 0.4292$ . Why this step? Concavity flips the sign compared to C2.

Verify: $f^{''} = - sin x < 0$ on $(0, π)$ , so with $E_{T} = + \frac{( b - a ) h ^{2}}{12} f^{''} (ξ)$ the plus coefficient times a negative $f^{''}$ gives $E_{T} < 0$ — rule below truth. Undershoot confirmed. Together C2 and C3 prove the sign of the trapezoid error equals the sign of $f^{''}$ : convex overshoots, concave undershoots.

Example 4 — Simpson is not exact (C4)

Look at figure s05. The lavender curve is $x^{4}$ ; the mint dashed curve is the parabola Simpson fits through the three nodes $(0, 0), (0.5, 0.0625), (1, 1)$ . Takeaway: the parabola cannot follow the sharp late upturn of $x^{4}$ , and the coral gap between them is the degree-4 error Simpson cannot absorb — this is why the answer is not exact.

Step 1. $h = 0.5$ ; nodes $0, 0.5, 1$ ; $f = 0, 0.0625, 1$ . Why this step? $x^{4}$ at $0.5$ is $0. 5^{4} = 0.0625$ — do not confuse with $x^{2}$ .

Step 2. $S_{2} = \frac{0.5}{3} [0 + 4 (0.0625) + 1] = \frac{0.5}{3} (1.25) = 0.208333 \dots$ Why this step? The 1-4-1 stencil, peak node weight 4.

Step 3. True value is $0.2$ , so the estimate overshoots: $E_{S} = 0.208333 - 0.2 = + 0.008333$ . Why this step? This is the case where the rule almost wins but leaks a fourth-degree error.

Verify — which error form, and why. Our whole interval $[0, 1]$ is a single Simpson pair (only $n = 2$ slices), so the natural tool is the per-pair error, not the composite one — using the composite $- \frac{( b - a ) h ^{4}}{180} f^{(4)}$ here would double-count the pair structure. With the per-pair form, width $2 h = 1$ and $f^{(4)} = 24$ : $E_{S}^{pair} = - \frac{( 2 h ) ^{5}}{2880} f^{(4)} (ξ) = - \frac{( 1 ) ^{5}}{2880} \cdot 24 = - 0.008333.$ This gives $E_{S} < 0$ under the truth-minus-rule grouping used to derive the classic $-$ sign. Our page convention is $E = rule - truth$ , which flips it to $+ 0.008333$ — matching the computed overshoot. The magnitude $0.008333$ is unambiguous; only the printed sign depends on the convention — exactly the sign trap flagged earlier. Either way, "Simpson is exact" holds only up to degree 3.

Example 5 — degenerate and reversed inputs (C5)

Step 1 — zero-width interval. $a = b = 2$ ⇒ $h = \frac{2 - 2}{n} = 0$ . Every term is multiplied by $h = 0$ , so $T_{n} = S_{n} = 0$ . Why this step? An interval of no width encloses no area — the formula must return $0$ , and it does because $h = 0$ zeroes the whole sum. No crash.

Step 2 — flat function. $f (x) \equiv 0$ on $[0, 5]$ . Every $f (x_{i}) = 0$ , so the weighted sum is $0$ ; $T = S = 0$ . True integral $\int_{0}^{5} 0 d x = 0$ . Exact. Why this step? Zero is a degree-0 polynomial — inside both rules' exactness, so no error possible.

Step 3 — single slice ( $n = 1$ ) and why Simpson can't. Trapezoid with $n = 1$ has no interior nodes; the sum $\sum_{i = 1}^{0}$ is empty and $T_{1} = \frac{h}{2} [f (a) + f (b)]$ — a valid single trapezoid. Simpson with $n = 1$ is mathematically impossible: Simpson's rule is the exact integral of the parabola through three nodes, but $n = 1$ gives only two nodes ( $x_{0}, x_{1}$ ). Two points determine a line, never a unique parabola, so there is no parabola to integrate — the rule is undefined, and the code's ValueError is the honest report of that. Why this step? Knowing it fails for a geometric reason (too few points), not a coding whim, stops off-by-one panics.

Step 4 — reversed limits. $\int_{1}^{0} x^{2} d x$ : here $a = 1, b = 0$ , so $h = \frac{0 - 1}{2} = - 0.5 < 0$ . Redo Example 2's trapezoid arithmetic with the negative $h$ : nodes $1, 0.5, 0$ , $f = 1, 0.25, 0$ , $T_{2} = \frac{- 0.5}{2} [1 + 0 + 2 (0.25)] = - 0.375$ . Why this step? The negative $h$ carries the sign for you, so the estimate is exactly $- 1 \times$ the forward answer.

Verify: (a) $0$ , (b) $0$ , (c) $T_{1}$ on $[0, 2]$ of $f = 3 x + 1$ gave $8$ in Example 1, (d) $- 0.375 = - T_{2}^{forward}$ , matching $\int_{1}^{0} = - \int_{0}^{1}$ . All consistent. The only genuine crash is $n = 0$ (then $h = (b - a) /0$ blows up), which is why production code must guard n >= 1.

Example 6 — the odd- $n$ Simpson trap (C6)

Step 1. $n = 3$ is odd. The 1-4-1 stencil eats slices two at a time; with 3 slices one is left with no partner parabola. Why this step? This is the structural reason parity matters, not a coding whim.

Step 2. The parent's code checks if n % 2 == 1: raise ValueError. So the call refuses loudly — the best outcome, better than a silent lie. Why this step? A raised error is a feature: it stops the #1 silent bug in its tracks.

Step 3 — the fix. Two clean options: (i) bump to the next even $n$ (here $n = 4$ ), or (ii) apply Simpson to the first even number of slices and handle the final dangling slice with a trapezoid. Simplest for an exam: choose $n$ even from the start. Why this step? You must be able to repair the case, not just detect it.

Verify: Redo $\int_{0}^{1} x^{2} d x$ with the fixed $n = 4$ : $h = 0.25$ , nodes $0, 0.25, 0.5, 0.75, 1$ , $f = 0, 0.0625, 0.25, 0.5625, 1$ . $S_{4} = \frac{0.25}{3} [0 + 1 + 4 (0.0625 + 0.5625) + 2 (0.25)] = \frac{0.25}{3} [1 + 4 (0.625) + 0.5] = \frac{0.25}{3} (4) = 0.33333 \dots$ — exact, as expected for a parabola. The even $n$ works flawlessly.

Example 7 — convergence / limiting behaviour (C7)

Look at figure s03. The horizontal axis is the number of subintervals $n$ (values $2, 4, 8$ ); the vertical axis is the trapezoid error on a logarithmic scale (so equal drops mean equal ratios). Takeaway: the three lavender dots fall in a straight line on the log axis, and each coral "÷4" label marks that doubling $n$ divides the error by exactly $4$ — the visual fingerprint of the $h^{2}$ factor.

Step 1. $n = 2$ : from Example 2, $T_{2} = 0.375$ , error $E_{T} = 0.375 - \frac{1}{3} = 0.041667$ . Why this step? Reuse the known value as the baseline.

Step 2. $n = 4$ : $h = 0.25$ , nodes $0, 0.25, 0.5, 0.75, 1$ , $f = 0, 0.0625, 0.25, 0.5625, 1$ . $T_{4} = \frac{0.25}{2} [0 + 1 + 2 (0.0625 + 0.25 + 0.5625)] = 0.125 [1 + 2 (0.875)] = 0.125 (2.75) = 0.34375$ . Error $= 0.010417$ . Why this step? One doubling of $n$ to test the ratio.

Step 3. $n = 8$ : computing the eight-node sum gives $T_{8} = 0.3359375$ , error $= 0.0026042$ . Why this step? A second doubling confirms the pattern is not a fluke.

Step 4 — conclude. Error ratios: $\frac{0.041667}{0.010417} = 4.00$ and $\frac{0.010417}{0.0026042} = 4.00$ . Exactly $4 \times$ per doubling. Why this step? This closes the demonstration: the halving of $h$ has been carried all the way to a numeric verdict.

Verify: The $4 \times$ drop is precisely the $h^{2}$ factor in $E_{T} = \frac{( b - a ) h ^{2}}{12} f^{''}$ : halving $h$ multiplies $h^{2}$ by $\frac{1}{4}$ . (This is exactly why Richardson extrapolation can combine $T_{n}$ and $T_{2 n}$ to cancel the $h^{2}$ term and manufacture Simpson.)

Example 8 — real-world word problem, data only (C8)

Step 1. Identify the integral. Distance is the area under the speed–time curve, $\int_{0}^{8} v d t$ , and our data spacing is uniform: $h = 2$ s. Why this step? Numerical integration is the only tool when you have samples, not a symbolic $f$ .

Step 2. Count slices: 5 points ⇒ $n = 4$ subintervals. $n = 4$ is even, so Simpson is allowed — use it for higher accuracy. Why this step? Points $= n + 1$ ; check parity before choosing Simpson.

Step 3 — Simpson. $S_{4} = \frac{h}{3} [v_{0} + v_{4} + 4 (v_{1} + v_{3}) + 2 v_{2}] = \frac{2}{3} [0 + 35 + 4 (10 + 30) + 2 (22)]$ $= \frac{2}{3} [35 + 160 + 44] = \frac{2}{3} (239) = 159.33$ m. Why this step? Odd-index nodes (10, 30) get weight 4; the even interior node (22) gets weight 2.

Step 4 — trapezoid cross-check. $T_{4} = \frac{2}{2} [0 + 35 + 2 (10 + 22 + 30)] = 1 \cdot [35 + 124] = 159$ m. Why this step? A second method within a hair of the first is strong evidence we did not blunder.

Verify: Both give $\approx 159$ m; they agree to $< 1 m$ , sensible for a car accelerating for 8 s. Units: $(m/s) \times s = m$ ✓ — an integral over time of speed yields distance, dimensionally correct.

Example 9 — exam twist: signed area (C9)

Look at figure s04. Horizontal axis $x$ on $[0, 2 π]$ , vertical axis $f (x)$ , with the $x$ -axis drawn in. The mint region is the positive hill on $[0, π]$ (area $+ 2$ ); the coral region is the negative valley on $[π, 2 π]$ (area $- 2$ ). Takeaway: numerical rules compute signed area — heights below the axis are negative and subtract — so the two regions cancel to $0$ .

Step 1. True value: $\int_{0}^{2 π} sin x d x = [- cos x]_{0}^{2 π} = - 1 - (- 1) = 0$ . Why this step? The truth exposes the trap: the answer is zero, not a large area.

Step 2. $h = \frac{2 π}{4} = \frac{π}{2}$ ; nodes $0, \frac{π}{2}, π, \frac{3 π}{2}, 2 π$ ; $f = 0, 1, 0, - 1, 0$ . Why this step? Note the negative sample $- 1$ at $\frac{3 π}{2}$ — that is where intuition fails newcomers.

Step 3. $S_{4} = \frac{π /2}{3} [0 + 0 + 4 (1 + (- 1)) + 2 (0)] = \frac{π}{6} \cdot 0 = 0$ . Why this step? The weight-4 peaks $+ 1$ and $- 1$ cancel; the answer is exactly $0$ .

Verify: $S_{4} = 0 =$ true value. If the exam wanted the geometric area (total shaded region regardless of sign), you would integrate $∣ sin x ∣$ or split at $x = π$ and add $\int_{0}^{π} = 2$ and $∣ \int_{π}^{2 π} ∣ = 2$ to get $4$ . Always ask: signed integral or geometric area? — they differ whenever the curve crosses the axis.

Connections

Parent topic — the derivations these examples exercise.
Riemann sums — the limit the exact cases (C1) reduce to.
Newton-Cotes formulas — why trapezoid caps at degree 1 and Simpson at degree 3 (C4).
Polynomial interpolation — Simpson integrates the interpolating parabola of Ex 4's figure.
Taylor series — origin of the $O (h^{2})$ and $O (h^{4})$ error terms derived and verified here.
Richardson extrapolation — the $4 \times$ error drop in Ex 7 is what makes it work.
scipy.integrate.quad — how production code handles the data/word-problem case (Ex 8) adaptively.

Worked examples — Implementing numerical integration from scratch — trapezoidal, Simpson's