2. Riemann sums and the area problem

📖 Reading · 15 min
💡 Every code box below is live — edit it and hit Run.

Forget antiderivatives for two lessons. Here is a completely separate question:

What is the area under a curve?

Nothing so far answers it. Rectangles have area base × height, but a curve's height won't hold still. So do what §0.0 did — approximate with something you can compute, then take the approximation to its limit.

The construction

To find the area under f from a to b:

1. Partition. Cut [a,b] into n subintervals. Equal widths are simplest:

\Delta x = \frac{b-a}{n}, \qquad x_i = a + i\Delta x

2. Sample. Pick a point x_i^* in each subinterval and use f(x_i^*) as that rectangle's height.

3. Sum.

S_n = \sum_{i=1}^{n}f(x_i^*)\,\Delta x

4. Take the limit.

\int_a^b f(x)\,dx = \lim_{n\to\infty}\sum_{i=1}^n f(x_i^*)\Delta x

That limit, when it exists, is the definite integral. The elongated S is a stylised sum — Leibniz's notation, and the dx is what \Delta x became.

Which sample point?

Three standard choices:

  • Left endpoint: x_i^* = x_{i-1}
  • Right endpoint: x_i^* = x_i
  • Midpoint: x_i^* = \frac{x_{i-1}+x_i}{2}

For an increasing function, left endpoints undershoot every rectangle and right endpoints overshoot, so

L_n \le \int_a^b f \le R_n

— a genuine bracket, and a squeeze in the sense of §1.3. For a decreasing function the roles swap.

The midpoint rule is much better than either, and not for the obvious reason. It isn't merely "between" them — its errors cancel. On each subinterval the midpoint rectangle cuts the curve, losing area on one side and gaining it on the other, and to first order those exactly balance. Midpoint error is O(1/n^2) against the endpoint rules' O(1/n), which §4.11 quantifies.

For a continuous function on a closed interval, the choice of sample point doesn't affect the limit — that's part of what integrability means, and §15.4 makes it precise.

Doing one exactly

\int_0^1 x^2\,dx

with right endpoints. Here \Delta x = \frac1n and x_i = \frac in:

R_n = \sum_{i=1}^n\left(\frac in\right)^2\cdot\frac1n = \frac{1}{n^3}\sum_{i=1}^n i^2

Use the closed form \sum_{i=1}^n i^2 = \frac{n(n+1)(2n+1)}{6}:

R_n = \frac{n(n+1)(2n+1)}{6n^3} = \frac{2n^3+3n^2+n}{6n^3} = \frac13+\frac{1}{2n}+\frac{1}{6n^2}

\lim_{n\to\infty}R_n = \frac13

So \int_0^1x^2dx = \frac13, established with no antiderivative anywhere.

Notice the error term. R_n - \frac13 = \frac{1}{2n} + O(1/n^2) — the error is proportional to \frac1n, so doubling n halves it. That's the O(1/n) claim, visible in the algebra.

Also notice how much work that was, for the easiest possible integrand. Doing \int_0^1 e^{x^2}dx this way is hopeless — there's no closed form for the sum. That's the problem the Fundamental Theorem solves, two lessons from now, and it's why this lesson feels laborious. It's supposed to.

Signed area

If f<0 the rectangles have negative height, so their contribution is negative. The definite integral computes signed area: regions below the axis count against you.

\int_0^{2\pi}\sin x\,dx = 0

not because there's no area, but because the hump above [0,\pi] exactly cancels the trough below [\pi,2\pi]. For actual geometric area you'd integrate |f| — a distinction §5.0 depends on, and the same one as displacement versus distance in §2.10.

The general definition

Equal widths and endpoint sampling are conveniences. The real definition allows any partition and any sample points:

\int_a^b f\,dx = \lim_{\|P\|\to0}\sum_{i=1}^n f(x_i^*)\Delta x_i

where \|P\| = \max\Delta x_i is the mesh, the widest subinterval. Requiring the mesh to shrink (rather than merely n\to\infty) prevents a partition that refines in one place while leaving a fat interval elsewhere.

f is Riemann integrable on [a,b] if this limit exists and is the same for every such choice. Continuous functions are integrable; so are bounded functions with finitely many discontinuities. The Dirichlet function from §1.2 is not — sampling rationals gives 1, irrationals gives 0, and the limit depends on choices it shouldn't. §15.4 gives the exact criterion.

Doing it in Python

The three rules, converging:

def riemann(f, a, b, n, kind="left"):
    w = (b - a) / n
    if kind == "left":
        pts = [a + i*w for i in range(n)]
    elif kind == "right":
        pts = [a + (i+1)*w for i in range(n)]
    else:
        pts = [a + (i+0.5)*w for i in range(n)]
    return sum(f(p) for p in pts) * w

f = lambda x: x*x

print(f"{'n':>7} {'left':>12} {'right':>12} {'midpoint':>12}")
for n in (4, 10, 100, 1000, 10000):
    print(f"{n:>7} {riemann(f,0,1,n,'left'):>12.8f} "
          f"{riemann(f,0,1,n,'right'):>12.8f} {riemann(f,0,1,n,'mid'):>12.8f}")

print(f"\ntrue value 1/3 = {1/3:.8f}")
print("left brackets from below, right from above, midpoint sits far closer")

The error orders, made explicit:

def riemann(f, a, b, n, kind):
    w = (b - a) / n
    off = {"left": 0.0, "right": 1.0, "mid": 0.5}[kind]
    return sum(f(a + (i + off) * w) for i in range(n)) * w

f = lambda x: x*x
true = 1/3

print(f"{'n':>7} {'left err':>13} {'ratio':>7} {'mid err':>13} {'ratio':>7}")
pl = pm = None
for n in (10, 20, 40, 80, 160, 320):
    el = abs(riemann(f,0,1,n,"left") - true)
    em = abs(riemann(f,0,1,n,"mid") - true)
    print(f"{n:>7} {el:>13.3e} {(pl/el if pl else 0):>7.2f} "
          f"{em:>13.3e} {(pm/em if pm else 0):>7.2f}")
    pl, pm = el, em

print("\ndoubling n halves the endpoint error (ratio 2, first order)")
print("and quarters the midpoint error (ratio 4, second order)")

The exact algebra, confirmed:

import sympy as sp

n, i = sp.symbols('n i', positive=True, integer=True)

R = sp.summation((i/n)**2 * (1/n), (i, 1, n))
print(f"R_n = {sp.simplify(R)}")
print(f"    = {sp.expand(sp.simplify(R))}")
print(f"limit as n -> oo: {sp.limit(R, n, sp.oo)}")

Signed area:

from math import sin, pi

def riemann(f, a, b, n=200000):
    w = (b - a) / n
    return sum(f(a + (i + 0.5) * w) for i in range(n)) * w

print(f"int of sin from 0 to pi     : {riemann(sin, 0, pi):>10.6f}")
print(f"int of sin from pi to 2pi   : {riemann(sin, pi, 2*pi):>10.6f}")
print(f"int of sin from 0 to 2pi    : {riemann(sin, 0, 2*pi):>10.6f}")
print(f"int of |sin| from 0 to 2pi  : "
      f"{riemann(lambda x: abs(sin(x)), 0, 2*pi):>10.6f}")

print("\nthe signed integral is 0; the geometric area is 4.")
print("below-axis regions subtract, exactly like displacement vs distance.")

Sample-point choice not mattering — until it does:

import random

random.seed(3)

def riemann_random(f, a, b, n):
    """Sample an arbitrary point in each subinterval."""
    w = (b - a) / n
    return sum(f(a + (i + random.random()) * w) for i in range(n)) * w

f = lambda x: x*x
print("random sample points inside each subinterval:")
for n in (100, 10_000, 1_000_000):
    print(f"  n={n:<9} {riemann_random(f, 0, 1, n):.8f}")

print(f"\nstill converging to {1/3:.8f}: for a continuous f the choice is irrelevant.")
print("for the Dirichlet function it decides the answer, and no limit exists.")

Worked example

Estimate \int_1^3\frac1x\,dx with n=4 using all three rules, and bracket the true value.

\Delta x = \frac{3-1}{4} = 0.5, with partition points 1, 1.5, 2, 2.5, 3.

Left endpoints (1, 1.5, 2, 2.5):

L_4 = 0.5\left(1 + \tfrac{1}{1.5} + \tfrac12 + \tfrac{1}{2.5}\right) = 0.5(1+0.6667+0.5+0.4) = 1.2833

Right endpoints (1.5, 2, 2.5, 3):

R_4 = 0.5\left(0.6667+0.5+0.4+0.3333\right) = 0.9500

Midpoints (1.25, 1.75, 2.25, 2.75):

M_4 = 0.5\left(0.8+0.5714+0.4444+0.3636\right) = 1.0897

The bracket. \frac1x is decreasing, so left endpoints overestimate and right endpoints underestimate — the reverse of the increasing case:

0.9500 \le \int_1^3\frac{dx}{x} \le 1.2833

True value: \ln 3 - \ln 1 = \ln 3 = 1.0986.

Inside the bracket ✓, and the midpoint estimate 1.0897 is off by 0.81%, while the endpoint rules are off by 13.5% (R_4) and 16.8% (L_4). Four rectangles, and the midpoint rule is already more than sixteen times more accurate than either endpoint rule.

Getting the bracket's direction right requires noticing the function decreases. Always check monotonicity before claiming which side you're on.

Your turn

1. Estimate \int_0^2x^2dx with n=4 right endpoints.

2. For \int_0^1x\,dx with n right endpoints, find R_n in closed form and its limit.

3. Why does \int_{-1}^{1}x^3dx = 0 without any computation?

4. For a decreasing f, which of L_n, R_n overestimates?

Solutions

1. \Delta x = 0.5, right endpoints 0.5, 1, 1.5, 2:

R_4 = 0.5\left(0.25 + 1 + 2.25 + 4\right) = 0.5(7.5) = \boxed{3.75}

True value is \frac83 \approx 2.667, so this overestimates badly — x^2 is increasing, and n=4 is coarse. The midpoint rule with the same 4 rectangles gives 2.625, off by 1.6%.

2. \Delta x = \frac1n, x_i = \frac in:

R_n = \sum_{i=1}^n\frac in\cdot\frac1n = \frac{1}{n^2}\sum_{i=1}^n i = \frac{1}{n^2}\cdot\frac{n(n+1)}{2} = \frac{n+1}{2n} = \frac12+\frac{1}{2n}

\lim_{n\to\infty}R_n = \boxed{\frac12}

Which must be right: the region is a triangle with base 1 and height 1. Note the error \frac{1}{2n} is positive and O(1/n) — first order, as promised.

3. x^3 is odd (§0.2), and the interval [-1,1] is symmetric about the origin. The area above the axis on [0,1] exactly mirrors the area below on [-1,0], so the signed contributions cancel.

The general rule, worth having:

\int_{-a}^{a}f = 0 \text{ if } f \text{ odd}, \qquad \int_{-a}^af = 2\int_0^af \text{ if } f \text{ even}

This is why §0.2 bothered with parity — it turns hard integrals into trivial ones, and it's how §8.5's Fourier series decides which coefficients vanish.

4. For a decreasing function, L_n overestimates and R_n underestimates.

On each subinterval a decreasing function is highest at the left edge, so a left-endpoint rectangle is taller than the region it stands on. The reverse of the increasing case — which is exactly the worked example's \frac1x.

The mnemonic to avoid: don't memorise "left is smaller". Memorise "the sample point at the higher end of the function overestimates", and read the monotonicity off the problem.

Check yourself in code

Compare the three Riemann rules on \int_0^1x^2dx = \frac13.

For n = 4, 10, 100, 1000, print the left, right and midpoint sums to 6 decimals.

Print exactly this:

n=4      left=0.218750  right=0.468750  mid=0.328125
n=10     left=0.285000  right=0.385000  mid=0.332500
n=100    left=0.328350  right=0.338350  mid=0.333325
n=1000   left=0.332834  right=0.333834  mid=0.333333
def riemann(f, a, b, n, offset):
    """offset 0.0 = left, 1.0 = right, 0.5 = midpoint."""
    w = (b - a) / n
    return sum(f(a + (i + offset) * w) for i in range(n)) * w

f = lambda x: x * x

for n in (4, 10, 100, 1000):
    # print the three sums for this n
    print(f"n={n:<6} ...")
def riemann(f, a, b, n, offset):
    """offset 0.0 = left, 1.0 = right, 0.5 = midpoint."""
    w = (b - a) / n
    return sum(f(a + (i + offset) * w) for i in range(n)) * w

f = lambda x: x * x

for n in (4, 10, 100, 1000):
    print(f"n={n:<6} left={riemann(f,0,1,n,0.0):.6f}  "
          f"right={riemann(f,0,1,n,1.0):.6f}  mid={riemann(f,0,1,n,0.5):.6f}")

Chop the region into rectangles, sum their areas, and let the widths shrink: the limit is the definite integral. Endpoint rules bracket the answer for a monotone function and converge at O(1/n); the midpoint rule's errors cancel and it converges at O(1/n^2). The integral is signed, so regions below the axis subtract. And computing one directly requires a closed form for the sum — which exists for x^2 and almost nothing else.

Next: what that limit actually asserts, and which functions have one.