Visualizing Numerical Integrals

Trapezoidal Rule

Approximating a definite integral $\int_a^b f(x) \;{\rm d}x$ for $x \in [a, b]$ using the trapezoidal rule. Here we assume that $a \leq b$ . Let $n \in \mathbb{N}$ be a positive integer and set $\Delta_x = \frac{b-a}{n}$ . Then we have the approximation

$\int_a^b f(x) \; {\rm d}x \approx \sum_{i=0}^{n-1} \frac {f(x_i) + f(x_i + \Delta_x)} {2} \cdot \Delta_x$

where $x_i = a + i \cdot \Delta_x$ . When $i = 0$ we have $a = x_0$ , and similarly $b = x_n$ .

%auto
# the function we are integrating numerically
var("x")
f(x) = 1/2 * x^2 - 2*x + sin(2*x)

def next_hue(dx):
    """
    Generator object for repetitive calls in the rgbcolor argument.
    """
    x = 0
    while True:
        if x > 1:
            x = 0
        yield hue(x, 0.3)
        x += dx

%auto
@interact
def _(ab=range_slider(-4, 10, 0.1, (-2,6), "[a, b]"), dx=slider(0.1,3,0.1,1,"dx")):
    a, b = ab
    xvals = srange(a, b, dx, include_endpoint=True)
    fx = [f(x) for x in xvals]
    fg = plot(f, (x,a,b))
    pts = point2d(zip(xvals,fx), rgbcolor="#ff3030", pointsize=30)
    lines = [([i,0], [i,f(i)], [i+dx, f(i+dx)], [i+dx, 0]) for i in xvals[0:-1]]
    nh = next_hue(1.0 / len(xvals))
    areas = sum(map(lambda p: polygon(p, rgbcolor=nh.next()), lines))
    show(fg + areas + pts)
    areas = [(f(i) + f(i+dx)) / 2 for i in xvals[0:-1]]
    print "approximation: %s" % float(dx * sum(areas))
    print "numerical integral: %s"  % numerical_integral(f, a, b)[0]

static example for fixed values, the above one is interactive

approximation: 5.51923393088
numerical integral: 4.58458454354

Rectangle Rule

Approximating a definite integral $\int_a^b f(x) \; {\rm d} x$ with $x \in [a,b]$ using rectangles. Again, we assume that $a \leq b$ , $n \in \mathbb{N}$ is positive, and we set $\Delta_x = \frac {b - a} {n}$ . Then we have the approximation

$\int_a^b f(x) \; {\rm d} x \approx \frac{f(a) \cdot \Delta_x}{2} + \sum_{i=1}^{n-1} f(x_i) \cdot \Delta_x + \frac{f(b) \cdot \Delta_x}{2}$

where $x_i = a+ i \cdot \Delta_x$ for $1 \leq i \leq n-1$ .

%auto
@interact
def _(ab=range_slider(-4,10,0.1,(-2,6),"[a, b]"), dx=slider(0.1,3,0.1,0.5,"dx")):
    a, b = ab
    xvals = srange(a, b, dx, include_endpoint=True)
    fx = [f(x) for x in xvals]
    fg = plot(f, (x,a,b))
    pts = point2d(zip(xvals,fx), rgbcolor="#ff3030", pointsize=30)
    #lines = [([i,0], [i,f(i)], [i+dx, f(i+dx)], [i+dx, 0]) for i in xvals[0:-1]]
    lines = [([xvals[0], 0], [xvals[0], f(xvals[0])], [xvals[0] + dx/2, f(xvals[0])], [xvals[0] + dx/2, 0])]
    for i in xvals[1:-1]:
        x1, x2 = i-dx/2, i+dx/2
        lines.append(([x1,0], [x1, f(i)], [x2, f(i)], [x2, 0]))
    x1, x2 = xvals[-1]-dx/2, xvals[-1]
    lines.append(([x1,0], [x1, f(x2)], [x2, f(x2)], [x2, 0]))
    nh = next_hue(1.0 / len(xvals))
    areas = sum(map(lambda p: polygon(p, rgbcolor=nh.next()), lines))
    show(fg + areas + pts)
    areas = [f(i) for i in xvals[1:-1]]
    areas.append(f(xvals[0]) / 2)
    areas.append(f(xvals[-1]) / 2)
    print "approximation: %s " % float(dx * sum(areas))
    print "numerical integral: %s"  % numerical_integral(f, a, b)[0]

static example with fixed values, the above one is interactive

approximation: 4.81471226681 
numerical integral: 4.58458454354