Finite difference method

The formula

The whole method rests on Taylor's theorem. Slide a function a small step h and difference the results:

Forward:   f'(x) ≈ [ f(x+h) − f(x) ] / h                 error O(h)
Backward:  f'(x) ≈ [ f(x) − f(x−h) ] / h                 error O(h)
Central:   f'(x) ≈ [ f(x+h) − f(x−h) ] / (2·h)           error O(h²)

Reading: the first derivative is a slope, so take a rise over a run using nearby samples. The central form straddles the point and cancels the leading error term, so it is one order more accurate for the same spacing.

Second derivative (the one that shows up in every diffusion and wave problem):

f''(x) ≈ [ f(x+h) − 2·f(x) + f(x−h) ] / h²               error O(h²)

Reading: curvature is the difference of two slopes. This three-point stencil is the workhorse of heat and wave codes.

The 1D heat equation ∂u/∂t = α·∂²u/∂x² becomes, with forward time and central space (FTCS):

u[j, n+1] = u[j, n] + r·( u[j−1, n] − 2·u[j, n] + u[j+1, n] )
where   r = α·Δt / Δx²

Reading: tomorrow's temperature at a node is today's value plus a nudge toward the average of its two neighbors. That scheme is stable only when r ≤ 1/2. That 1/2 is the 1D limit; in 2D and 3D on equal spacing it tightens to r ≤ 1/4 and r ≤ 1/6 (the general condition is α·Δt·Σ(1/Δxᵢ²) ≤ 1/2).

Where you meet it

Thermal analysis at a review board. Somebody drops a transient heat map of an avionics box or a rocket nozzle liner on the screen. Under the hood of most in-house thermal solvers is a finite-difference (or finite-volume) grid marching temperature forward in time. When the reviewer asks "did your time step respect the stability limit," they are asking about that r condition — and since a box or a liner is a 3D model, the answer had better be r ≤ 1/6, not the 1D 1/2.
FDTD on the RF bench. Antenna and radome engineers around Redstone run FDTD codes (Ansys, CST, or open-source Meep) to predict radiation patterns and shielding. Every one of them stakes a Yee grid over the geometry and leapfrogs E and H fields in time.
Structural transients and CFD. Explicit dynamics for a blast or impact case, and many early CFD marching schemes, are finite differences dressed up. The Courant number you set in the input deck is the same stability constraint under a different name.
Quick-and-dirty field checks. Column of thermocouple readings across a wall, and you want the heat flux — a central difference of the temperature profile gives dT/dx on the spot without a solver.

How it works

The method converts calculus into arithmetic on a grid. That is its whole appeal and its whole trap.

The accuracy trap is order. A forward difference of sin(x) at x = 1 with h = 0.1 gives 0.4974 against the true 0.5403 — off by 0.043. The central difference gives 0.5394, off by only 0.0009, roughly fifty times better for the same two samples. Halve h on a central scheme and the error drops by four (O(h²)), while a forward scheme only halves it. Engineers reach for central differences on interior points and accept lopsided one-sided stencils at boundaries.

The stability trap is worse, because it does not warn you gently — it detonates. Run the FTCS heat scheme at r = 0.4 for 200 steps and the solution decays smoothly, peak amplitude around 0.14. Bump to r = 0.6 and the same 200 steps produce a peak near 5·10¹¹. Nothing about the code changed except one number crept past 1/2. This is the mistake people make: they refine the spatial grid to get a prettier answer, Δx shrinks, r = α·Δt/Δx² explodes past 1/2, and the run blows up. Refine space and you must shrink time as the square — cut Δx in half and Δt must drop by four. Von Neumann stability analysis, which tracks how each Fourier mode grows per step, is the tool that predicts exactly where the cliff sits.

FDTD carries the same constraint under the name Courant number. On a cubic Yee grid, c·Δt ≤ Δx/√3 in three dimensions (Δx/√2 in 2D, Δx in 1D). Exceed it and the fields grow without bound. The Yee grid's own cleverness is staggering: E-field components sit on cell edges, H-field components on the faces, and the two leapfrog a half step apart in time. That interleaving makes each curl a natural central difference and keeps the scheme second-order without any extra bookkeeping.

Two more limits worth knowing. Explicit schemes (the ones above) are cheap per step but stability-capped on Δt; implicit schemes (Crank–Nicolson, backward Euler) stay stable at any step but make you solve a linear system every step — you trade a matrix solve for freedom from the CFL cliff. And finite differences want structured grids. Curved boundaries and complex geometry are where the finite element and finite volume methods took the field away from pure FDM, because they conform to messy shapes that a rectangular stencil cannot.

History

The bones are old. Approximating derivatives by differences goes back to Newton and Euler, and Carl Runge was numerically attacking differential equations by 1900. But the moment finite differences became a real theory for partial differential equations is dated precisely: 1928, when Richard Courant, Kurt Friedrichs, and Hans Lewy published their paper on solving problems of mathematical physics by finite differences ^[1]^[2]. Buried in it was the condition that now carries all three names — the CFL condition — a necessary limit on the time step relative to the grid spacing for a scheme to converge. They proved it before there was a single computer to blow up.

The tools sharpened during and after World War II, when machines finally made the arithmetic practical. John von Neumann developed the Fourier-mode stability analysis that lets you check a scheme cheaply, and the 1950s work of Friedrichs, Lax, and Wendroff built out the theory for hyperbolic problems like shocks and waves ^[1].

The branch most Huntsville RF engineers touch daily came in 1966. Kane Yee, publishing in IEEE Transactions on Antennas and Propagation, took Maxwell's curl equations and discretized them on the staggered space-and-time grid that now bears his name ^[3]^[4]. For over a century after Maxwell, electromagnetics had been solved almost entirely in the frequency domain; Yee marched the actual fields forward in time instead. The method sat quietly until Allen Taflove, in 1980, validated it on wave-penetration problems and coined the acronym FDTD ^[3]. That name, and Yee's little staggered cell, are now baked into nearly every commercial electromagnetics solver.

Related tools

/tools/heat-conduction
/tools/thermal-expansion
/tools/heat-exchanger-duty
/tools/antenna-gain-beamwidth
/tools/frequency-wavelength

Sources

Written by HE in our own words from the cited sources — engineering judgment included, your stamp still required. All entries →