Chapter 18
20 min read
Section 155 of 353

Double Integrals over Rectangles

Multiple Integrals

Learning Objectives

By the end of this section, you will be able to:

  1. Understand how the definite integral extends from one variable to two variables, and what double integrals represent geometrically
  2. Construct double Riemann sums to approximate the volume under a surface
  3. Apply Fubini's theorem to evaluate double integrals as iterated integrals
  4. Compute double integrals over rectangular regions using both integration orders
  5. Calculate volumes, average values, and physical quantities using double integrals
  6. Connect double integrals to machine learning applications including probability distributions and expected values

The Big Picture: From Lines to Surfaces

"Double integrals let us sum infinitely many infinitesimally small contributions over a two-dimensional region—finding volumes, masses, and probabilities that single integrals cannot reach."

In single-variable calculus, we learned that the definite integral abf(x)dx\int_a^b f(x)\,dx represents the signed area under a curve. Now we take a momentous step: extending integration to functions of two variables.

For a function f(x,y)f(x, y) defined over a region RR in the xyxy-plane, the double integral Rf(x,y)dA\iint_R f(x, y)\,dA represents the signed volume between the surface z=f(x,y)z = f(x, y) and the region RR.

Why Double Integrals Matter

Double integrals appear throughout science, engineering, and data science:

  • Physics: Mass, center of mass, moments of inertia of laminae
  • Engineering: Stress distributions, heat flow, fluid dynamics
  • Probability: Joint probability distributions and expected values
  • Machine Learning: Bayesian inference, kernel methods, multivariate integration
  • Computer Graphics: Rendering, texture mapping, illumination models

Historical Context

The theory of multiple integration developed gradually in the 17th and 18th centuries. Isaac Newton and Gottfried Leibniz established single-variable calculus, but extending their ideas to multiple dimensions required new insights.

Leonhard Euler (1707-1783) was among the first to systematically work with double integrals, applying them to problems in mechanics and astronomy. Joseph-Louis Lagrange (1736-1813) developed the theory further in his work on celestial mechanics.

The modern rigorous treatment came from Guido Fubini (1879-1943), who proved the fundamental theorem bearing his name: under suitable conditions, a double integral can be evaluated as two successive single integrals.

Fubini's Legacy

Fubini's theorem is one of the most practically important results in multivariable calculus. It transforms the daunting task of integrating over a region into manageable iterated single integrals—a reduction that makes computation possible.


The Volume Problem

Consider a surface z=f(x,y)z = f(x, y) hovering above a rectangular region R=[a,b]×[c,d]R = [a, b] \times [c, d] in the xyxy-plane. If f(x,y)0f(x, y) \geq 0, we want to find the volume of the solid that lies:

  • Above the rectangle RR
  • Below the surface z=f(x,y)z = f(x, y)

This is the two-dimensional analogue of finding the area under a curve. Just as we used Riemann sums to approximate area, we will use double Riemann sums to approximate volume.

The Fundamental Question

Given a surface z=f(x,y)z = f(x, y) over a rectangle R=[a,b]×[c,d]R = [a, b] \times [c, d]

What is the volume between the surface and the rectangle?


The Double Integral: Definition

We define the double integral using a limiting process, just as we did for single integrals.

Double Riemann Sums

Step 1: Partition the rectangle. Divide [a,b][a, b] into nn subintervals and [c,d][c, d] into mm subintervals, creating n×mn \times m subrectangles RijR_{ij}.

Step 2: Choose sample points. In each subrectangle RijR_{ij}, pick a sample point (xij,yij)(x_{ij}^*, y_{ij}^*).

Step 3: Form the Riemann sum. Each term f(xij,yij)ΔAijf(x_{ij}^*, y_{ij}^*) \cdot \Delta A_{ij} represents the volume of a rectangular box with base RijR_{ij} and height f(xij,yij)f(x_{ij}^*, y_{ij}^*).

Double Riemann Sum

Sn,m=i=1nj=1mf(xij,yij)ΔAijS_{n,m} = \sum_{i=1}^{n} \sum_{j=1}^{m} f(x_{ij}^*, y_{ij}^*) \Delta A_{ij}

where ΔAij=ΔxiΔyj\Delta A_{ij} = \Delta x_i \cdot \Delta y_j is the area of each subrectangle

Step 4: Take the limit. As the partition becomes finer (more and smaller subrectangles), the Riemann sum approaches the double integral:

Definition: The Double Integral

Rf(x,y)dA=limP0i=1nj=1mf(xij,yij)ΔAij\iint_R f(x, y)\,dA = \lim_{\|P\| \to 0} \sum_{i=1}^{n} \sum_{j=1}^{m} f(x_{ij}^*, y_{ij}^*) \Delta A_{ij}

The limit exists if f is continuous on R (and more generally, if f is bounded with discontinuities of measure zero)

In the notation dAdA, we write dA=dxdydA = dx\,dy or dA=dydxdA = dy\,dx, representing an infinitesimal area element.

Interactive: 3D Riemann Sums

Explore how rectangular boxes approximate the volume under a surface. Increase the number of partitions to see the approximation improve:

Riemann Sum

Approximation:9.0000
Exact value:5.3333
Error:3.6667
% Error:68.75%
Total boxes:
4 × 4 = 16
Key insight: As we increase the number of partitions, the sum of the box volumes approaches the true volume under the surface. This limit is the double integral.

Interactive: Volume Under a Surface

Visualize the double integral as the volume between a surface and a rectangular region. Adjust the function and bounds to explore different scenarios:

Function:
Drag to rotate

Understanding the Visualization

  • Blue rectangle: Region R on the xy-plane
  • Colored surface: Graph of z = f(x, y)
  • Yellow lines: Connect region to surface
  • The double integral gives the volume between the surface and the region R
R f(x, y) dA = Volume

Iterated Integrals

While the definition using Riemann sums is conceptually important, we need a practical way to compute double integrals. The key insight: we can evaluate a double integral by performing two single integrals in succession.

Fubini's Theorem

Fubini's Theorem (for Rectangles)

If f(x,y)f(x, y) is continuous on R=[a,b]×[c,d]R = [a, b] \times [c, d], then:

Rf(x,y)dA=cd[abf(x,y)dx]dy=ab[cdf(x,y)dy]dx\iint_R f(x, y)\,dA = \int_c^d \left[\int_a^b f(x, y)\,dx\right]dy = \int_a^b \left[\int_c^d f(x, y)\,dy\right]dx

In words: We can integrate first with respect to one variable (treating the other as a constant), then integrate the result with respect to the remaining variable. The order doesn't matter—both give the same answer.

OrderNotationProcess
dx dy∫_c^d ∫_a^b f(x,y) dx dy1. Integrate over x (y constant) 2. Integrate result over y
dy dx∫_a^b ∫_c^d f(x,y) dy dx1. Integrate over y (x constant) 2. Integrate result over x

Choosing the Order

While both orders give the same result for rectangles, sometimes one order leads to simpler antiderivatives. If the inner integral looks difficult, try the other order!

Interactive: Iterated Integrals

See how double integrals are computed step by step using iterated integrals. Watch the cross-sections as we integrate in each direction:

Integrate x first, then y

R f(x, y) dA = ∫0102 (x² + y) dx dy

Step 1: Inner integral (fix y, integrate over x)

02 (x² + y) dx = [x³/3 + xy]02 = 8/3 + 2y

Step 2: Outer integral (integrate over y)

01 (8/3 + 2y) dy = [8y/3 + y²]01 = 11/3

Result: 11/3 ≈ 3.667
Fubini's Theorem:

For continuous f on a rectangle, the order of integration can be switched. Both orders give the same result.


Properties of Double Integrals

Double integrals share many properties with single integrals. Let ff and gg be integrable over a region RR:

Fundamental Properties

1.
Linearity: R[f+g]dA=RfdA+RgdA\iint_R [f + g]\,dA = \iint_R f\,dA + \iint_R g\,dA
2.
Scalar multiplication: RcfdA=cRfdA\iint_R cf\,dA = c\iint_R f\,dA
3.
Additivity over regions: If R=R1R2R = R_1 \cup R_2 with no overlap, then RfdA=R1fdA+R2fdA\iint_R f\,dA = \iint_{R_1} f\,dA + \iint_{R_2} f\,dA
4.
Comparison: If fgf \leq g on R, then RfdARgdA\iint_R f\,dA \leq \iint_R g\,dA
5.
Area formula: R1dA=Area(R)\iint_R 1\,dA = \text{Area}(R)

Average Value

The average value of f(x,y)f(x, y) over a region RR is defined as:

favg=1Area(R)Rf(x,y)dAf_{avg} = \frac{1}{\text{Area}(R)} \iint_R f(x, y)\,dA

This generalizes the single-variable formula favg=1baabf(x)dxf_{avg} = \frac{1}{b-a}\int_a^b f(x)\,dx. It represents the "average height" of the surface over the region.

Interpretation

If you were to spread the volume under the surface uniformly over the region R, the average value would be the height of the resulting flat surface.


Worked Examples

Example 1: Basic Double Integral

Evaluate R(x2+y)dA\iint_R (x^2 + y)\,dA where R=[0,2]×[0,1]R = [0, 2] \times [0, 1].

Solution: Using Fubini's theorem with dx dy:

0102(x2+y)dxdy\int_0^1 \int_0^2 (x^2 + y)\,dx\,dy

Inner integral (integrate over x, y is constant):

02(x2+y)dx=[x33+xy]02=83+2y\int_0^2 (x^2 + y)\,dx = \left[\frac{x^3}{3} + xy\right]_0^2 = \frac{8}{3} + 2y

Outer integral (integrate over y):

01(83+2y)dy=[8y3+y2]01=83+1=113\int_0^1 \left(\frac{8}{3} + 2y\right)dy = \left[\frac{8y}{3} + y^2\right]_0^1 = \frac{8}{3} + 1 = \frac{11}{3}

Example 2: Volume Under a Paraboloid

Find the volume under z=4x2y2z = 4 - x^2 - y^2 over the square R=[0,1]×[0,1]R = [0, 1] \times [0, 1].

Solution:

V=0101(4x2y2)dxdyV = \int_0^1 \int_0^1 (4 - x^2 - y^2)\,dx\,dy

Inner integral:

01(4x2y2)dx=[4xx33xy2]01=413y2=113y2\int_0^1 (4 - x^2 - y^2)\,dx = \left[4x - \frac{x^3}{3} - xy^2\right]_0^1 = 4 - \frac{1}{3} - y^2 = \frac{11}{3} - y^2

Outer integral:

01(113y2)dy=[11y3y33]01=11313=103\int_0^1 \left(\frac{11}{3} - y^2\right)dy = \left[\frac{11y}{3} - \frac{y^3}{3}\right]_0^1 = \frac{11}{3} - \frac{1}{3} = \frac{10}{3}

Applications

Physics: Mass and Center of Mass

For a thin plate (lamina) occupying region R with density function ρ(x,y)\rho(x, y) (mass per unit area):

QuantityFormulaPhysical Meaning
MassM = ∬_R ρ(x,y) dATotal mass of the lamina
First moment (x)M_y = ∬_R x·ρ(x,y) dAMoment about the y-axis
First moment (y)M_x = ∬_R y·ρ(x,y) dAMoment about the x-axis
Center of mass(x̄, ȳ) = (M_y/M, M_x/M)Balance point of the lamina

Probability: Joint Distributions

For continuous random variables X and Y with joint probability density function f(x,y)f(x, y):

Probability of a region: P((X,Y)R)=Rf(x,y)dAP((X, Y) \in R) = \iint_R f(x, y)\,dA

Normalization: f(x,y)dxdy=1\iint_{-\infty}^{\infty} f(x, y)\,dx\,dy = 1

Expected value: E[g(X,Y)]=g(x,y)f(x,y)dxdyE[g(X, Y)] = \iint_{-\infty}^{\infty} g(x, y)f(x, y)\,dx\,dy


Machine Learning Connection

Double integrals are essential in machine learning, appearing whenever we work with multivariate distributions or compute expectations:

Applications in ML

1. Bayesian Inference

Computing posterior probabilities requires integrating over parameter spaces:

P(θD)=P(Dθ)P(θ)P(Dθ)P(θ)dθP(\theta | D) = \frac{P(D|\theta)P(\theta)}{\iint P(D|\theta')P(\theta')\,d\theta'}

2. Gaussian Mixture Models

Marginalizing over latent variables involves double integrals over continuous spaces.

3. Kernel Density Estimation

Normalizing kernel estimates requires integrating the estimated density over the domain.

4. Feature Engineering

Computing aggregate statistics over 2D feature spaces uses double integrals or their discrete approximations.

Numerical Integration

In practice, most ML applications use numerical methods (Monte Carlo, quadrature) to approximate double integrals. Understanding the theory helps you choose appropriate methods and estimate errors.


Python Implementation

Let's implement double integrals computationally:

Numerical Approximation

Computing Double Integrals Numerically
🐍double_integral_numerical.py
3Riemann Sum Function

This function computes a double Riemann sum using the midpoint rule for sample points. We partition the rectangle into n×m subrectangles.

16Area Element

dA = dx × dy is the area of each subrectangle. The sum Σ f(x*ᵢⱼ, y*ᵢⱼ) dA approximates the double integral.

21Midpoint Rule

We use the center of each subrectangle as the sample point. This typically gives better accuracy than corner points.

33SciPy Integration

scipy.integrate.dblquad uses adaptive quadrature to compute double integrals to high precision.

36 lines without explanation
1import numpy as np
2from scipy import integrate
3
4def double_integral_riemann(f, x_range, y_range, n=100, m=100):
5    """
6    Approximate double integral using Riemann sums.
7
8    f: function of (x, y)
9    x_range: tuple (x_min, x_max)
10    y_range: tuple (y_min, y_max)
11    n, m: number of partitions in x and y directions
12    """
13    x_min, x_max = x_range
14    y_min, y_max = y_range
15
16    dx = (x_max - x_min) / n
17    dy = (y_max - y_min) / m
18    dA = dx * dy
19
20    total = 0.0
21    for i in range(n):
22        for j in range(m):
23            # Midpoint sample
24            x = x_min + (i + 0.5) * dx
25            y = y_min + (j + 0.5) * dy
26            total += f(x, y) * dA
27
28    return total
29
30# Example: Volume under f(x,y) = x² + y² over [0,2]×[0,2]
31f = lambda x, y: x**2 + y**2
32
33# Using our Riemann sum
34approx = double_integral_riemann(f, (0, 2), (0, 2), n=50, m=50)
35print(f"Riemann sum approximation: {approx:.6f}")
36
37# Using scipy for comparison
38result, error = integrate.dblquad(f, 0, 2, 0, 2)
39print(f"SciPy result: {result:.6f}")
40print(f"Exact value (32/3): {32/3:.6f}")

Symbolic Iterated Integrals

Symbolic Integration with SymPy
🐍iterated_integral.py
3Symbolic Double Integral

SymPy allows us to compute double integrals exactly and symbolically, showing each step of the iterated integral.

12Inner Integral

First, integrate f(x, y) with respect to x, treating y as a constant. This gives a function of y alone.

16Outer Integral

Then integrate the result with respect to y to get the final answer. This is Fubini's theorem in action.

31Verification

By Fubini's theorem, integrating in the opposite order (dy then dx) should give the same result, confirming our calculation.

33 lines without explanation
1import sympy as sp
2
3def symbolic_double_integral(f_expr, x, y, x_bounds, y_bounds):
4    """
5    Compute double integral symbolically using SymPy.
6
7    f_expr: SymPy expression for f(x, y)
8    x, y: SymPy symbols
9    x_bounds: tuple (x_min, x_max)
10    y_bounds: tuple (y_min, y_max)
11    """
12    # Inner integral: integrate over x first
13    inner = sp.integrate(f_expr, (x, x_bounds[0], x_bounds[1]))
14    print(f"Inner integral (∫ f dx): {inner}")
15
16    # Outer integral: integrate over y
17    result = sp.integrate(inner, (y, y_bounds[0], y_bounds[1]))
18    print(f"Outer integral (∫∫ f dx dy): {result}")
19
20    return result
21
22# Define symbols and function
23x, y = sp.symbols('x y', real=True)
24f = x**2 + y**2
25
26print("Computing ∬ (x² + y²) dA over [0,2]×[0,2]")
27print("=" * 50)
28
29# Integrate in order dx dy
30result1 = symbolic_double_integral(f, x, y, (0, 2), (0, 2))
31print(f"\nResult: {result1} = {float(result1):.4f}")
32
33# Verify by integrating in opposite order dy dx
34print("\nVerifying with opposite order (dy dx)...")
35inner2 = sp.integrate(f, (y, 0, 2))
36result2 = sp.integrate(inner2, (x, 0, 2))
37print(f"Result: {result2} = {float(result2):.4f}")

Applications: Mass and Center of Mass

Physical Applications
🐍physics_applications.py
3Mass Calculation

For a lamina with density ρ(x, y), the total mass is the double integral of density over the region.

17Center of Mass

The center of mass coordinates involve moments: x̄ = M_y/M and ȳ = M_x/M, where M_x = ∬ y·ρ dA and M_y = ∬ x·ρ dA.

24First Moment

The moment about the y-axis (M_y) is the integral of x weighted by density. This determines the x-coordinate of the center of mass.

45 lines without explanation
1import numpy as np
2from scipy import integrate
3
4def compute_mass(density, x_range, y_range):
5    """
6    Compute mass of a lamina with variable density.
7    Mass = ∬_R ρ(x, y) dA
8    """
9    result, _ = integrate.dblquad(
10        density,
11        x_range[0], x_range[1],
12        lambda x: y_range[0],
13        lambda x: y_range[1]
14    )
15    return result
16
17def compute_center_of_mass(density, x_range, y_range):
18    """
19    Compute center of mass (x̄, ȳ) of a lamina.
20    x̄ = (1/M) ∬_R x·ρ(x, y) dA
21    ȳ = (1/M) ∬_R y·ρ(x, y) dA
22    """
23    M = compute_mass(density, x_range, y_range)
24
25    # Moment about y-axis (for x̄)
26    Mx = lambda x, y: x * density(x, y)
27    moment_x, _ = integrate.dblquad(
28        Mx, x_range[0], x_range[1],
29        lambda x: y_range[0], lambda x: y_range[1]
30    )
31
32    # Moment about x-axis (for ȳ)
33    My = lambda x, y: y * density(x, y)
34    moment_y, _ = integrate.dblquad(
35        My, x_range[0], x_range[1],
36        lambda x: y_range[0], lambda x: y_range[1]
37    )
38
39    return moment_x / M, moment_y / M
40
41# Example: Density ρ(x,y) = 1 + x + y on [0,1]×[0,1]
42density = lambda x, y: 1 + x + y
43
44mass = compute_mass(density, (0, 1), (0, 1))
45x_bar, y_bar = compute_center_of_mass(density, (0, 1), (0, 1))
46
47print(f"Mass: {mass:.4f}")
48print(f"Center of mass: ({x_bar:.4f}, {y_bar:.4f})")

Test Your Understanding

Test Your Understanding

Score: 0/8

Question 1 of 8

What does the double integral ∬_R f(x, y) dA represent geometrically when f(x, y) ≥ 0?


Summary

The Core Idea

The double integral Rf(x,y)dA\iint_R f(x, y)\,dA extends single integration to two dimensions. When f0f \geq 0, it represents the volume under the surface z=f(x,y)z = f(x, y) and above the region R.

Key Formulas

ConceptFormula
Double integral∬_R f(x,y) dA
Riemann sumΣᵢ Σⱼ f(x*ᵢⱼ, y*ᵢⱼ) ΔAᵢⱼ
Iterated integral (dx dy)∫_c^d ∫_a^b f(x,y) dx dy
Iterated integral (dy dx)∫_a^b ∫_c^d f(x,y) dy dx
Average valuef_avg = (1/Area(R)) ∬_R f dA
Area of R∬_R 1 dA = Area(R)

Key Takeaways

  1. Definition: The double integral is the limit of double Riemann sums as the partition becomes infinitely fine.
  2. Fubini's theorem: For continuous f on a rectangle, we can evaluate the double integral as two successive single integrals in either order.
  3. Interpretation: When f ≥ 0, the double integral gives volume. When f represents density, it gives mass.
  4. Properties: Double integrals are linear and additive over regions—they inherit all the properties of single integrals.
  5. Applications: From physics (mass, moments) to probability (expected values) to machine learning (Bayesian inference).
From Area to Volume:
"Double integrals let us sum over two-dimensional regions—the natural extension that opens the door to volumes, densities, and multivariate probability."
Looking Ahead: This section covered double integrals over rectangles. In the next section, we'll extend to general regions—where the limits of integration depend on the other variable, allowing us to integrate over triangles, circles, and any bounded region.
Loading comments...