Finite difference method

The finite difference method replaces continuous derivatives with algebraic quotients evaluated on a grid of points. For a function f(x), the forward difference approximates f′(x) as [f(x+h)−f(x)]/h, the backward difference uses [f(x)−f(x−h)]/h, and the central difference gives [f(x+h)−f(x−h)]/(2h). The central difference has truncation error of order h², making it more accurate for smooth functions.¹ Higher-order derivatives are obtained by combining these stencils, such as the second derivative [f(x+h)−2f(x)+f(x−h)]/h². These approximations are derived from Taylor series expansions, and the method's accuracy depends on the step size h and the order of the stencil.

To solve an ordinary or partial differential equation, the domain is discretized into a grid, and each derivative is replaced by a finite difference approximation. For example, the heat equation ∂u/∂t = α ∂²u/∂x² becomes a system of algebraic equations linking values at neighboring time and space points. Explicit schemes, such as the forward-time central-space (FTCS) method, compute new values directly from old ones but are conditionally stable, requiring a small time step relative to the spatial step (the Courant–Friedrichs–Lewy condition).² Implicit schemes, like the backward-time central-space (BTCS) method, solve a linear system at each step and are unconditionally stable, though more computationally expensive. The method also handles boundary conditions by incorporating them into the stencil equations, and it extends naturally to nonlinear problems through iterative techniques.

The finite difference method traces its roots to Leonhard Euler in the 18th century, who used difference quotients to approximate solutions of differential equations. In the early 20th century, Lewis Fry Richardson applied finite differences to weather prediction, and the method gained prominence with the advent of digital computers in the 1950s. Variants include the Crank–Nicolson method, which averages explicit and implicit schemes for second-order accuracy in time, and the leapfrog method, used in wave equations. For complex geometries, the finite difference method is often less flexible than the finite element method, but it remains efficient for structured grids. Modern extensions include compact finite differences, which achieve high accuracy with small stencils, and summation-by-parts operators that ensure stability for hyperbolic problems.

Beyond standard applications, the finite difference method has niche uses in option pricing, where the Black–Scholes equation is discretized on a grid of asset prices and time. In seismology, it models wave propagation through heterogeneous media, and in astrophysics, it simulates stellar structure and accretion disks. A subtle edge case is the treatment of discontinuities, where naive finite differences produce oscillations; techniques like flux limiters or shock-capturing schemes are required. Another overlooked detail is the use of non-uniform grids, where the step size varies to concentrate points near sharp gradients, improving efficiency. The method also appears in image processing for edge detection, where discrete Laplacians approximate second derivatives. Historically, Richardson's 1922 work on weather forecasting was famously inaccurate due to insufficient data, but it laid the groundwork for numerical weather prediction.³

Stencil

The set of grid points used in a finite difference approximation.

Truncation error

The error introduced by approximating a derivative with a finite difference, proportional to a power of the step size.

Explicit scheme

A time-stepping method where new values are computed directly from known values.

Implicit scheme

A time-stepping method that requires solving a system of equations at each step.