CFD Solver Methods — Pressure-Velocity Coupling and Discretization

For incompressible flow, the density is constant, so there's no thermodynamic equation of state linking density, pressure, and temperature. Pressure appears in the momentum equation as a gradient term, but there is no dedicated pressure transport equation — instead, pressure is implicitly constrained by the divergence-free condition on velocity ($\nabla \cdot \mathbf{u} = 0$). You have to find a pressure field that, when you solve momentum, gives you a velocity field that satisfies continuity. It's a constraint rather than a direct equation. The SIMPLE algorithm is the most widely used approach to handle this: it iterates between a momentum "predictor" step and a pressure "corrector" step until both are satisfied simultaneously.

Sure. SIMPLE stands for Semi-Implicit Method for Pressure-Linked Equations, developed by Patankar and Spalding in 1972. Here are the three stages repeated until convergence. Step 1: solve the momentum equation using the current (guessed) pressure field to get a predicted velocity u* that does not yet satisfy continuity. Step 2: derive a pressure correction equation by requiring the corrected velocity to be divergence-free — this gives a Poisson equation for the pressure correction p'. Step 3: update both velocity and pressure using the corrections. Then repeat. The under-relaxation factors control how aggressively you update between iterations — too aggressive and you diverge; too conservative and convergence is slow.

Great practical question. SIMPLE is for steady-state only — it relies on under-relaxation to damp out the time-evolution and find a steady solution. PISO (Pressure Implicit with Splitting of Operators) is for transient flows: it does multiple pressure corrections within each time step to satisfy continuity tightly, with no under-relaxation. It requires a small Courant number (Co < 1) to be stable. PIMPLE is OpenFOAM's combination of PISO outer loops with SIMPLE-like under-relaxation — it can handle Courant numbers up to 5–10, making it much more practical for LES and transient RANS where the time step is set by physics rather than stability. SIMPLEC is an improved SIMPLE with a tighter pressure-velocity coupling that allows larger under-relaxation for pressure (up to 1.0), improving convergence speed.

Under-relaxation limits how much any field changes between SIMPLE iterations: $\phi^{new} = \phi^{old} + \alpha(\phi^{computed} - \phi^{old})$. With $\alpha = 1$ you apply the full correction — fast convergence but risk of divergence. With $\alpha = 0.1$ you take tiny steps — very stable but painfully slow. The pressure and velocity factors are the most sensitive. The rule of thumb taught in every CFD course: $\alpha_p + \alpha_U \approx 1$, so the default $\alpha_U = 0.7$, $\alpha_p = 0.3$ is a safe starting point. For highly non-linear or buoyancy-driven flows, try $\alpha_U = 0.5$, $\alpha_p = 0.2$. For simple attached flow, you can push to $\alpha_U = 0.9$ to speed things up. Turbulence quantities (k, ε, ω) typically use 0.5–0.7.

The residuals in OpenFOAM are normalized initial residuals — the ratio of the current equation imbalance to the initial imbalance. There is no universal threshold, but a good rule of thumb: continuity residual below 1e-3, momentum and turbulence below 1e-4 for engineering accuracy. For aeroacoustics or precise heat transfer, push to 1e-5 or below. More importantly, also monitor your key output variable — drag coefficient, pressure drop, maximum temperature — and check that it has stopped changing. It's common for all residuals to look "converged" but the actual quantity of interest is still drifting if the convergence is asymptotic.

GAMG is Geometric Algebraic MultiGrid — it's dramatically faster for the pressure Poisson equation because pressure information propagates across the whole domain (it's an elliptic equation). A simple iterative solver like Gauss-Seidel handles nearby cells fast but takes ages to propagate information across the domain — that's the "multigrid problem." GAMG builds successively coarser representations of the mesh, solves cheaply on the coarse levels to spread information globally, then refines back. For a million-cell pressure field, GAMG can be 50–100x faster than a non-multigrid solver. Velocity and turbulence are advection-dominated (hyperbolic-like) — they don't need multigrid; smoothSolver or PBiCGStab with an ILU preconditioner works well and is simpler to configure.