Finite-Difference Pricing: Stability, Boundaries and the Illusion of Precision

The main implementation challenge is to distinguish a stable calculation from an accurate and economically consistent one. This reference combines grid construction, time-stepping schemes, exotic payoff treatment and error diagnostics.

Finite differences replace derivatives with local relationships between grid values. This turns a continuous valuation equation into a sequence of algebraic problems. The method is particularly useful when a contract has low-dimensional state variables and explicit exercise or boundary conditions.

Choose a coordinate system before choosing a scheme

Let $u_i^n$ denote the value at spatial node $x_i$ and time-to-maturity level $\tau_n$. With spacing $h$, central approximations are

$$u_x\approx\frac{u_{i+1}-u_{i-1}}{2h},\qquad u_{xx}\approx\frac{u_{i+1}-2u_i+u_{i-1}}{h^2}.$$

In words: Estimate slope from the two neighbouring values and curvature from the neighbours relative to the central value.

For sufficiently smooth functions, both have second-order spatial truncation error. Near payoff discontinuities or nonsmooth boundaries, that nominal order may not describe observed convergence.

A uniform underlying-price grid and a uniform log-price grid concentrate resolution differently. Log coordinates can simplify the diffusion coefficient in basic option models, while a price grid can make strikes and barriers easier to align. Nonuniform grids can concentrate nodes where value or Greeks change rapidly.

The grid should represent important contractual levels accurately. A barrier placed between nodes and rounded to the nearest node is a modified contract unless an appropriate interpolation or boundary treatment corrects the displacement.

The three standard time schemes

After spatial discretisation, write the evolution as $du/d\tau=Lu$, where $L$ is a matrix representing the differential operator and boundary contributions are handled consistently.

Explicit Euler uses

$$u^{n+1}=u^n+kLu^n,$$

where $k$ is the time step. It is inexpensive per step but conditionally stable.

Backward Euler uses

$$(I-kL)u^{n+1}=u^n.$$

It requires a linear solve and strongly damps high-frequency numerical components. For standard diffusion problems it has favourable stability properties but first-order time accuracy.

Crank–Nicolson averages the two operator evaluations:

$$(I-\tfrac12kL)u^{n+1}=(I+\tfrac12kL)u^n.$$

It is second-order in time for smooth problems. Its stability does not guarantee absence of oscillations or preservation of positivity for arbitrary grids, coefficients and nonsmooth inputs.

A concrete stability bound

For the heat equation $u_\tau=a u_{xx}$, with positive diffusion coefficient $a$, the explicit update is

$$u_i^{n+1}=\lambda u_{i-1}^n+(1-2\lambda)u_i^n+\lambda u_{i+1}^n, \qquad \lambda=\frac{ak}{h^2}.$$

In words: The next central value combines its current value and its two neighbours, weighted by how quickly values spread and the sizes of the time and space steps.

The familiar condition $0\leq\lambda\leq1/2$ makes these weights nonnegative and provides the basic one-dimensional explicit diffusion stability restriction. Halving the spatial spacing then requires roughly quartering the time step.

This bound cannot be transplanted unchanged to every financial PDE. State-dependent diffusion, drift, discounting, nonuniform grids and multiple dimensions alter the discrete operator and its restrictions. In an underlying-price Black–Scholes grid, diffusion grows with the square of the price coordinate, so the most restrictive node may be far from the strike.

Stability, consistency and convergence answer different questions

Consistency, stability and convergence test different properties

For suitable well-posed linear problems, consistency and stability underpin convergence. For nonlinear control equations, monotonicity and other conditions become important. It is unsafe to turn this into a universal statement that every implicit method automatically converges to the economically correct answer.

A stable coarse grid can produce a smooth but biased price. A price can look converged while gamma remains noisy because numerical differentiation amplifies grid error. Separate convergence evidence is needed for the quantities actually used in hedging and risk.

Nonsmooth payoffs need deliberate treatment

A vanilla option payoff has a kink at the strike. A digital payoff has a discontinuity. Starting a high-order scheme directly from these features can create oscillations, especially in sensitivities near maturity.

Short initial backward-Euler steps before switching to Crank–Nicolson provide damping. Payoff averaging over grid cells, suitable grid alignment and local refinement can also help. The appropriate combination depends on the contract and the required outputs.

Smoothing should approximate the original payoff under a controlled scheme rather than quietly redefine it. The same smoothing choice should be examined under grid refinement to confirm that it approaches the intended contract.

Boundary conditions are part of the computation

An infinite state domain must usually be truncated. For a call at a sufficiently large underlying price, the asymptotic value is approximately underlying price minus discounted strike, under the basic no-dividend assumptions. Moving the outer boundary farther away tests truncation error separately from mesh error.

A continuously monitored down-and-out option with no rebate has zero value at its lower knockout barrier. That boundary applies at every time level. A discretely monitored barrier should instead be enforced only at the monitoring dates.

For an American option, each time step must incorporate the exercise constraint. Simply solving the European system and checking exercise only at maturity misses the embedded right. Projection methods, penalty methods or complementarity solvers provide different implementations, with their own accuracy and convergence checks.

Path-dependent contracts enlarge the state

An Asian option depending on an arithmetic average requires a running average or accumulated-sum state in addition to current price. A lookback option requires a running extremum. These extra variables create multidimensional grids and can introduce degenerate directions because some state components have deterministic evolution conditional on the underlying.

Alternating-direction methods split multidimensional operations to reduce computational cost. Mixed derivative terms, such as those arising from correlated price and variance shocks, require careful treatment. Operator splitting adds another approximation whose stability and error depend on the scheme.

Finite elements can be useful for irregular domains or complex boundary structure. Multigrid methods accelerate certain large linear solves by addressing errors at several spatial scales. Neither technique removes the need for a correct financial formulation.

Grid Greeks make errors visible

Once values are computed, delta and gamma can be estimated by finite differences in the relevant state. The retained teaching example uses neighbouring option values 10, 12 and 15 at equally spaced underlying nodes with spacing 2. The central estimates are

$$\Delta\approx\frac{15-10}{2(2)}=1.25, \qquad \Gamma\approx\frac{15-2(12)+10}{2^2}=0.25.$$

In words: The neighbouring values give a central slope of one and a quarter and a central curvature of one quarter.

These are arithmetic grid illustrations, not a claim about a standard vanilla call: a delta above one would require contextual explanation and would violate the usual no-dividend vanilla-call bound in the basic model. Numerical examples should be checked for economic interpretation as well as arithmetic correctness.

A theta estimate requires a clear calendar-time versus time-to-maturity sign convention. Mixing the two changes the sign. All Greek outputs should state bump sizes, units, interpolation method and whether other calibration inputs were held fixed or rebuilt.

A refinement record is stronger than one attractive number

Run a sequence of spatial and temporal refinements while varying the outer boundary independently. Record price, key Greeks and computational cost. If doubling resolution barely changes price but materially changes gamma, the engine may be sufficient for valuation and insufficient for a gamma-sensitive hedge.

For an error expansion $V(h)=V^*+Ch^p+o(h^p)$, Richardson extrapolation can estimate a limiting value from two grids. It is useful only when the asymptotic regime and order $p$ are credible. Applying it before payoff and boundary errors settle can create spurious precision.

Independent checks for a finite-difference valuation engine

Work through one explicit grid update

For the heat-equation scheme, take illustrative neighbouring values 10, 12 and 15 with $\lambda=0.25$. The new central value is $0.25(10)+0.5(12)+0.25(15)=12.25$. The update is a weighted average, so it remains between the smallest and largest neighbouring values.

If $\lambda=0.75$, the central weight becomes minus 0.5. The same local values happen to give a finite result, but the scheme no longer has the nonnegative-weight property and can amplify high-frequency errors. A single apparently sensible step is therefore insufficient evidence of stability.

Financial operators also contain drift and discounting. Their weights need not be identical to this simple heat-equation example, and the admissible time step should be checked for the actual operator. Testing only one central node can miss instability near a domain edge where coefficients differ.

Interpolation is part of the reported answer

The market state at which a price is requested will not always coincide with a grid node. Interpolating the grid solution introduces another approximation. Linear interpolation can give a reasonable price while producing a discontinuous delta as the requested state crosses nodes.

Higher-order interpolation can smooth the output but may overshoot or violate payoff bounds near nonsmooth regions. The interpolation method should be selected with the required output in mind and tested alongside grid refinement.

The same issue arises in calibration. If an optimiser repeatedly requests values between nodes, interpolation noise can contaminate its objective and its derivatives. A small calibration residual may then partly reflect numerical artefacts rather than a stable fit to instrument prices.

Exercise projection needs a convergence check

For an American-style claim, a simple discrete-time procedure takes the larger of continuation and immediate exercise at each eligible date. With a coarse timeline, this represents exercise only at those dates. Refining the timeline improves the approximation to continuous exercise, subject to the spatial solver and obstacle treatment.

The exercise boundary can be more sensitive than the option price. A small value error near the point of indifference can shift the inferred boundary materially. If the application needs exercise decisions, boundary convergence should be reviewed directly rather than inferred from price convergence alone.

Dividend dates and other jumps in contractual state can change the optimal exercise region abruptly. The timeline should resolve those events explicitly. Smearing a discrete event across several ordinary time steps changes both the contract and the exercise decision.

Computational efficiency should preserve the evidence

For one-dimensional central-difference problems, the implicit matrix is often tridiagonal, allowing a linear-time solve in the number of spatial nodes. If coefficients and step sizes stay constant, some factorisations can be reused. These are algorithmic savings that preserve the underlying numerical method.

Adaptive grids, adaptive time steps and operator splitting change the approximation and require renewed error checks. Accelerating a solver is useful only if it preserves the accuracy of the outputs that matter. Price, Greeks, exercise boundary and barrier sensitivity can have different computational requirements, so the performance target should identify which outputs are being protected.

The direction of the initial damping steps also matters. In time-to-maturity coordinates they begin at the known terminal payoff and move toward today’s value. Applying the damping at the opposite end of the calculation may fail to suppress the nonsmooth high-frequency components when they enter the solution. A reproducible numerical specification should therefore identify the number and size of damping steps and the events after which damping is reapplied, rather than listing only the name of the time-stepping scheme.

Finite differences are most valuable when the grid, equation and contract remain inspectable. The objective is a controlled approximation with documented residual error, not merely a price produced by a stable solver.