CFD · Incompressible Flow · Numerical Analysis
📅 July 2026 ⏱ ≈ 10 min read 🎯 Advanced

Spurious modes in the incompressible Navier-Stokes equations

The incompressible Navier-Stokes equations couple velocity and pressure through a constraint — ∇·u = 0 (the divergence of the velocity field is zero everywhere, meaning fluid volume is neither created nor destroyed) — rather than an evolution equation for pressure. Discretise velocity and pressure on the same grid points without care, and the discrete system admits non-physical pressure patterns that satisfy the discrete divergence constraint exactly while being completely wrong: checkerboard oscillations that grow unchecked and never appear in the true continuous solution. Understanding why they arise — and how staggered grids and stabilised methods eliminate them — is essential before trusting any pressure field a solver produces.

TL;DR: When velocity and pressure are stored at the same grid points, incompressible-flow solvers can produce fake "checkerboard" pressure patterns that the equations can't detect or damp. The fix is to separate where pressure and velocity live (staggered MAC grids), reconstruct velocities more carefully (Rhie-Chow interpolation), or add stabilising terms in finite-element methods.

Written by MySimulator Team · Reviewed by MySimulator Editorial Review

1. Why pressure has no evolution equation

The incompressible Navier-Stokes system is:

∂u/∂t + (u·∇)u = −∇p/ρ + ν∇²u (momentum)
∇·u = 0 (continuity / incompressibility)

There is no ∂p/∂t term anywhere — pressure is not an independent physical field that evolves on its own; it is a Lagrange multiplier enforcing the divergence-free constraint at every instant. Taking the divergence of the momentum equation and using ∇·u = 0 yields the pressure Poisson equation ∇²p = −ρ∇·((u·∇)u), which must be solved consistently with the discrete velocity field at every time step — this is the coupling that a naive discretisation can get wrong.

2. The checkerboard pressure mode

On a collocated grid — velocity and pressure both stored at cell centres — the standard second-order central difference for the pressure gradient at node i uses only p_{i-1} and p_{i+1}:

(∂p/∂x)_i ≈ (p_{i+1} − p_{i-1}) / (2Δx)

This stencil is completely blind to p_i itself. Consider an alternating pressure field p_i = (−1)^i · A — a "checkerboard" pattern. Its central-difference gradient evaluates to zero at every node, identical to the gradient of a uniform (constant) pressure field. The discrete momentum equation cannot distinguish a real, physically meaningless checkerboard mode from a harmless constant offset — both are invisible to the discrete gradient operator, and nothing in the discrete system damps them out.

+A
−A
+A
−A
−A
+A
−A
+A
+A
−A
+A
−A

In 2D/3D the same failure occurs along both directions simultaneously, and the null space of the discrete gradient operator can be large. Round-off error or a slightly asymmetric boundary condition is enough to seed this mode, and because nothing in the equations penalises it, it can grow until it dominates the visible pressure field — a symptom every CFD practitioner eventually meets.

3. The LBB / inf-sup compatibility condition

The rigorous statement of "which velocity/pressure discretisation pairs are safe" is the Ladyzhenskaya-Babuška-Brezzi (LBB) condition, also called the inf-sup condition. For a saddle-point problem with velocity space Vh and pressure space Qh, stability requires the existence of β > 0, independent of the mesh size h, such that:

infq∈Q_h supv∈V_h (∇·v, q) / (‖v‖ ‖q‖) ≥ β > 0

Intuitively: for every non-trivial pressure field q (including a checkerboard candidate), the velocity space must contain a velocity field whose divergence "sees" q with a mesh-independent strength. Equal-order collocated velocity/pressure (Q1-Q1 or P1-P1 finite elements, or a naive finite-difference grid with both fields at the same points) generically fails LBB — the checkerboard mode is exactly the pressure field that the discrete divergence operator cannot "see" at any velocity.

4. Staggered (MAC) grids

Harlow & Welch's Marker-and-Cell (MAC) grid (1965) is the classical, elegant fix: store pressure at cell centres, but store the horizontal velocity component u on the vertical cell faces and the vertical component v on the horizontal cell faces:

(∂p/∂x)_{i+1/2,j} = (p_{i+1,j} − p_{i,j}) / Δx ← now uses adjacent p's, no skipped node

(∇·u)_{i,j} = (u_{i+1/2,j} − u_{i-1/2,j})/Δx + (v_{i,j+1/2} − v_{i,j-1/2})/Δy

On the staggered grid, the pressure gradient and divergence operators are exact discrete adjoints of one another and involve every neighbouring pressure value — the checkerboard mode now produces a large, correctly-signed gradient and is damped by the discrete Poisson solve. MAC grids satisfy LBB by construction and remain the standard choice for finite-difference/finite-volume incompressible solvers, including the projection-method solver used elsewhere on this site.

Trade-off: staggered storage complicates interpolation (velocity components live at different points than pressure and each other), boundary condition bookkeeping, and makes higher-order or unstructured-mesh extensions harder — which is why collocated grids remain popular despite requiring extra stabilisation.

5. Rhie-Chow interpolation for collocated grids

Rhie & Chow (1983) showed that a collocated grid can be made stable without physically staggering the storage, by momentum-interpolating the face velocity used in the mass-conservation (pressure) equation rather than simply averaging the neighbouring cell-centre velocities:

u_face = ū_face − D_face · [ (∂p/∂x)_face − (∂p/∂x)‾_face ]

ū_face: linear average of cell-centre velocities
D_face: momentum-equation coefficient (≈ Δt/ρ or 1/A_P)
(∂p/∂x)_face: true gradient across the face — sensitive to checkerboarding
(∂p/∂x)‾_face: interpolated average gradient — insensitive to it

The face-velocity correction term explicitly reintroduces sensitivity to the pressure difference across the face that the naive average would miss, effectively adding a controlled fourth-order pressure-smoothing term. This is the dominant technique in modern collocated finite-volume codes (OpenFOAM, most commercial CFD packages), because it keeps every field at the same location — simplifying unstructured meshes, multigrid, and boundary treatment — while still satisfying a discrete analogue of LBB.

6. Spurious modes in finite elements — stabilisation

The same phenomenon appears in finite-element Stokes/Navier-Stokes solvers. Equal-order P1-P1 (linear velocity, linear pressure) elements violate LBB. Two families of fixes are standard:

PSPG stabilisation term (added to weak continuity equation):
+ Σ_K τ_K ∫_K ∇q · [ −∇p + ρ((u·∇)u) − ν∇²u ] dΩ

τ_K: element-size-dependent stabilisation parameter
Vanishes when the strong-form residual is zero → consistent

7. Detecting spurious modes in practice

8. Practical implications for real-time solvers

Real-time WebGL fluid solvers (Jos Stam's "Stable Fluids" and its descendants, used elsewhere on this site) typically dodge the checkerboard problem by storing velocity on a staggered MAC-style grid inside the fragment shader ping-pong buffers, or by deliberately over-smoothing the pressure solve (a handful of Jacobi iterations rather than a converged solve) — a pragmatic trade-off that hides low-level checkerboard artefacts behind numerical dissipation, acceptable for visual simulation but not for engineering-accurate CFD.

🌊 Navier-Stokes Fluid (WebGL)

See a staggered-grid pressure-projection solver running live — no checkerboarding, real-time Jacobi pressure solve.

Open Fluid Sim →