Numerical PDEs (Finite Differences)

Solve heat, wave and Laplace equations on a grid. Learn grid discretisation, finite-difference stencils, explicit (FTCS) and implicit schemes, Crank–Nicolson, stability conditions (r = αΔt/Δx² ≤ ½, CFL), iterative solution of Laplace's equation, and verification against exact solutions — with a demonstration of instability.

When geometry or boundary conditions defeat separation of variables, engineers discretise the domain into a grid and replace derivatives with finite differences — the foundation of CFD and thermal simulation codes.

An analogy

Finite differences turn a smooth field into a chessboard: each square's temperature is updated from its neighbours' values, like people adjusting to the average opinion of those around them.

Discretisation

Replace the domain by nodes; derivatives by differences (central in space). Boundary conditions set edge nodes (Dirichlet) or ghost nodes (Neumann).

Explicit vs implicit

Explicit schemes update each node directly (cheap per step) but have stability limits. Implicit schemes solve a linear system each step (tridiagonal in 1D — Thomas algorithm) but allow large steps.

Stability conditions

Heat FTCS: r ≤ ½. Wave: CFL C ≤ 1 (information cannot travel more than one cell per step). Von Neumann analysis derives these limits.

Verification and beyond

Compare with exact solutions and refine grids (error ∝ Δx², Δt for FTCS). Finite volumes (conservation, CFD) and finite elements (complex geometry) extend the ideas.

Formulas

• FTCS (explicit heat): uᵢⁿ⁺¹ = uᵢⁿ + r(uᵢ₊₁ⁿ − 2uᵢⁿ + uᵢ₋₁ⁿ), stable if r ≤ ½
• Crank–Nicolson: Average explicit and implicit — unconditionally stable, O(Δt², Δx²)
• Wave equation (CFL): C = cΔt/Δx ≤ 1
• Laplace 5-point stencil: u_{i,j} = (u_{i+1,j} + u_{i−1,j} + u_{i,j+1} + u_{i,j−1})/4

Worked example: Stable time step

Given: α = 1.11 × 10⁻⁴ m²/s (copper), Δx = 1 cm

Find: Max explicit Δt.

• Δt ≤ 0.5 × 10⁻⁴/1.11 × 10⁻⁴ = 0.45 s

Answer: Δt ≤ 0.45 s

Worked example: Laplace iteration

Given: Square plate: top at 100 °C, others 0 °C; one interior node

Find: Centre temperature.

• u = (100 + 0 + 0 + 0)/4

Answer: 25 °C (finer grids converge to the same value by symmetry)

Explicit heat equation: stable vs unstable time step

import numpy as np

# Explicit FTCS for u_t = a u_xx, u(0)=u(1)=0, u(x,0)=sin(pi x)
a, nx = 1.0, 21
dx = 1.0 / (nx - 1)
for r in (0.4, 0.6):                       # r = a dt / dx^2, stable if r <= 0.5
    dt = r * dx**2 / a
    x = np.linspace(0, 1, nx)
    u = np.sin(np.pi * x)
    u[nx // 2] += 1e-6                    # tiny disturbance (e.g. round-off)
    steps = int(round(0.2 / dt))
    for _ in range(steps):
        u[1:-1] = u[1:-1] + r * (u[2:] - 2 * u[1:-1] + u[:-2])
    exact = np.exp(-np.pi**2 * 0.2) * np.sin(np.pi * 0.5)
    print(f"r = {r}: u(0.5, 0.2) = {u[nx // 2]: .4e}   exact {exact:.4e}")
# r = 0.4: u(0.5, 0.2) =  1.3812e-01   exact 1.3891e-01
# r = 0.6: u(0.5, 0.2) = -6.6412e+11   exact 1.3891e-01

With r = 0.4 the scheme matches the exact solution to 0.6 %. With r = 0.6 a one-in-a-million disturbance grows to 10¹¹ — violating the stability limit destroys the solution no matter how accurate the formula looks.

Common mistakes

• Wrong: Ignoring stability limits. — Right: Check r ≤ ½ or CFL ≤ 1, or use implicit schemes.
• Wrong: Trusting a single grid. — Right: Do grid-refinement studies.

Applications

• Thermal simulation: Transient heating of components.
• CFD: Finite-volume flow solvers.
• Acoustics & seismics: Wave propagation with CFL limits.

Key takeaways

• Grid + difference stencils.
• Explicit: cheap but r ≤ ½ / CFL ≤ 1.
• Implicit/CN: stable, solve systems.
• Verify with exact solutions and refinement.

Going further

• Von Neumann analysis: Fourier-mode stability proofs.
• Finite volume method: Conservative CFD discretisation.
• ADI methods: Efficient implicit 2D/3D schemes.