Lesson 8 of 13 25 min

Pressure-Velocity Coupling

For incompressible flow, there's no explicit equation for pressure. The continuity equation $\nabla \cdot \mathbf{u} = 0$ constrains velocity but doesn't directly give pressure. This creates the fundamental challenge of incompressible CFD: pressure-velocity coupling.

The Problem

Why No Pressure Equation?

In compressible flow, we have:

  • Continuity → density
  • Momentum → velocity
  • Energy → temperature
  • Equation of state → pressure (from $\rho, T$)

In incompressible flow:

  • $\rho = \text{constant}$ → no density variation
  • Continuity becomes $\nabla \cdot \mathbf{u} = 0$ → a velocity constraint
  • No explicit equation for $p$!

Pressure appears in momentum but isn't the primary unknown — it acts as a Lagrange multiplier enforcing mass conservation.

The Coupled System

The incompressible Navier-Stokes equations:

$$\nabla \cdot \mathbf{u} = 0$$

$$\rho \frac{\partial \mathbf{u}}{\partial t} + \rho (\mathbf{u} \cdot \nabla) \mathbf{u} = -\nabla p + \mu \nabla^2 \mathbf{u}$$

The coupling:
  • Momentum needs $\nabla p$ to compute $\mathbf{u}$
  • Continuity needs $\mathbf{u}$ to be divergence-free
  • But we don't know either!

The Checkerboard Problem

Visualize why co-located grids can produce checkerboard pressure patterns.

Why Collocated Grids Are Tricky

On a collocated grid (pressure and velocity at same points), the pressure gradient uses every-other-point:

$$\left.\frac{\partial p}{\partial x}\right|_P \approx \frac{p_E - p_W}{2\Delta x}$$

This doesn't "see" pressure at point P! A checkerboard pattern satisfies:

$$p_E - p_W = 0 \quad \text{(if } p_E = p_W \text{)}$$

Even with wildly oscillating pressure: $p = [1, 0, 1, 0, 1, 0]$

The Solution: Rhie-Chow Interpolation

Add pressure gradient smoothing when computing face velocities:

$$u_f = \overline{u} - \left(\frac{\partial p}{\partial x}\right)_f^{\text{smooth}}$$

This couples adjacent pressure values, eliminating checkerboarding.

The SIMPLE Algorithm

Step through the SIMPLE algorithm interactively.

Semi-Implicit Method for Pressure-Linked Equations

SIMPLE (Patankar & Spalding, 1972) is the most widely used pressure-velocity coupling method.

Step 1: Guess Pressure Field

Start with initial pressure $p^\ast$ (from previous iteration or initial condition).

Step 2: Solve Momentum with Guessed Pressure

$$a_P u^\ast_P = \sum a_{nb} u^\ast_{nb} + b - \frac{\partial p^\ast}{\partial x} V_P$$

This gives velocity field $\mathbf{u}^\ast$ that generally doesn't satisfy continuity.

Step 3: Derive Pressure Correction Equation

If the true pressure is $p = p^\ast + p^\prime$, and true velocity is $\mathbf{u} = \mathbf{u}^\ast + \mathbf{u}^\prime$:

From momentum:

$$a_P u^\prime_P = -\frac{\partial p^\prime}{\partial x} V_P$$

(neglecting neighbor velocity corrections — the "SIMPLE approximation")

Substituting into continuity $\nabla \cdot \mathbf{u} = 0$:

$$\nabla \cdot \left( \frac{V_P}{a_P} \nabla p^\prime \right) = \nabla \cdot \mathbf{u}^\ast$$

This is a Poisson equation for $p^\prime$!

Step 4: Solve Pressure Correction

$$a_P p^\prime_P = a_E p^\prime_E + a_W p^\prime_W + a_N p^\prime_N + a_S p^\prime_S + b_P$$

Where $b_P$ contains the mass imbalance from $\mathbf{u}^\ast$.

Step 5: Correct Pressure and Velocity

$$p = p^\ast + \alpha_p p^\prime$$

$$u_P = u^\ast_P - \frac{V_P}{a_P} \frac{\partial p^\prime}{\partial x}$$

Under-relaxation $\alpha_p$ (typically 0.3) stabilizes the iteration.

Step 6: Solve Other Scalars

With updated $\mathbf{u}$, solve:

  • Energy equation
  • Turbulence equations ($k$, $\epsilon$)
  • Species transport

Step 7: Check Convergence

If residuals are small enough, stop. Otherwise, use updated $p$ as new $p^\ast$ and repeat.

SIMPLE Variants

SIMPLEC (Consistent)

Van Doormaal & Raithby (1984) improved SIMPLE by not neglecting neighbor velocity corrections:

$$a_P u^\prime_P = \sum a_{nb} u^\prime_{nb} - \frac{\partial p^\prime}{\partial x} V_P$$

Using $(a_P - \sum a_{nb}) u^\prime_P$ instead of $a_P u^\prime_P$:

Result: Faster convergence, can use $\alpha_p = 1$ in some cases.

PISO (Pressure-Implicit Split-Operator)

Issa (1986) designed for transient flows:

  • Predictor: One SIMPLE-like step
  • First corrector: Correct pressure and velocity
  • Second corrector: Additional correction step
Result: Second-order in time without sub-iterations.

Comparison

MethodIterationsUnder-relaxationBest For
SIMPLEManyRequired (0.3-0.7)Steady-state
SIMPLECFewerCan be higherSteady-state
PISOFew per timestepNot neededTransient

Under-Relaxation

Why Under-Relax?

The SIMPLE approximation (neglecting $\sum a_{nb} u^\prime$) overestimates corrections. Without under-relaxation, the solution oscillates or diverges.

How It Works

Explicit under-relaxation:

$$\phi^{new} = \phi^{old} + \alpha(\phi^{calc} - \phi^{old})$$

Implicit under-relaxation (modifies coefficients):

$$\frac{a_P}{\alpha} \phi_P = \sum a_{nb} \phi_{nb} + S + \frac{(1-\alpha)}{\alpha} a_P \phi^{old}_P$$

Typical Values

VariableUnder-relaxation $\alpha$
Pressure0.3
Momentum0.7
Turbulence0.8
Energy0.9-1.0
Too low: Slow convergence Too high: Oscillations or divergence

Implementation Details

Face Mass Flux

The key quantity is face mass flux $\dot{m}_f$:

$$\dot{m}_f = (\rho u)_f A_f$$

This must be consistent across the momentum-pressure correction cycle.

Boundary Conditions

Velocity inlet:
  • $\dot{m}_f$ specified
  • No pressure correction at inlet face
Pressure outlet:
  • $p$ specified
  • $p^\prime = 0$ at outlet

Monitoring Convergence

Track:

  • Scaled residuals (normalized imbalances)
  • Mass imbalance (should approach zero)
  • Key quantities (drag, lift, mass flow)

Alternative Methods

Projection Methods

  • Solve momentum without pressure → $\mathbf{u}^\ast$
  • Project onto divergence-free space:
$$\mathbf{u} = \mathbf{u}^\ast - \nabla \phi$$

Where $\nabla^2 \phi = \nabla \cdot \mathbf{u}^\ast$

  • Recover $p$ from $\phi$

Used in: DNS, LES, academic codes

Coupled Solvers

Solve momentum and continuity simultaneously:

$$\begin{bmatrix} A & G \\ D & 0 \end{bmatrix} \begin{bmatrix} u \\ p \end{bmatrix} = \begin{bmatrix} f \\ 0 \end{bmatrix}$$

Pros: Robust, no under-relaxation needed Cons: Larger system, more memory

Used in: Some commercial codes for difficult cases

Practical Tips

Convergence Issues

SymptomPossible CauseFix
Slow convergenceUnder-relaxation too lowIncrease carefully
Oscillating residualsUnder-relaxation too highDecrease $\alpha$
Mass imbalance stuckPoor mesh qualityFix mesh
Checkerboard pressureMissing Rhie-ChowEnable interpolation

Starting a Simulation

  • Initialize with potential flow or uniform field
  • Start with lower under-relaxation
  • Increase as solution develops
  • Monitor mass balance throughout

Key Takeaways

  • No pressure equation for incompressible flow — pressure couples velocity to continuity
  • SIMPLE derives a pressure correction equation from mass conservation
  • Checkerboard pressure requires Rhie-Chow interpolation on collocated grids
  • Under-relaxation essential because SIMPLE neglects neighbor corrections
  • SIMPLEC is faster; PISO is better for transient
  • Mass balance is the key convergence indicator

What's Next

With the discretized equations assembled, we need to solve them. The next lesson covers Iterative Solvers & Convergence — how to efficiently solve large sparse linear systems and monitor convergence.

Convection Schemes