OpenFOAM Incompressible Fluid Module (v14)

Overview of the Incompressible Fluid Module

The incompressibleFluid module in OpenFOAM v14 provides solvers for flows where the fluid density is treated as constant, valid when the Mach number remains below approximately 0.3. In this regime, pressure waves propagate effectively instantaneously relative to convective time scales, and the divergence of the velocity field vanishes:

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

This assumption is not merely a convenience; for most industrial applications involving water, oil, or low-speed air (HVAC, automotive external aerodynamics, turbomachinery casings), it yields orders-of-magnitude savings in computational cost with negligible accuracy loss. The incompressible solvers operate on a pressure-based formulation, solving a Poisson equation for pressure derived from the continuity and momentum equations. The core mathematical challenge is pressure-velocity coupling, since pressure does not appear explicitly in a standalone transport equation for incompressible flow.

Incompressible Flow :: Foundational theory, the derivation of the continuity equation, and the mathematical structure of the Navier-Stokes equations under the incompressibility constraint.

Steady-State: simpleFoam

The simpleFoam solver is the workhorse of industrial incompressible CFD. It solves the steady-state, incompressible, turbulent Navier-Stokes equations using the Simple algorithm (Semi-Implicit Method for Pressure-Linked Equations). The SIMPLE algorithm proceeds as follows:

1. Solve the momentum predictor equation to obtain a velocity field. 2. Construct and solve a pressure correction (Poisson) equation. 3. Correct the pressure and velocity fields to satisfy continuity. 4. Solve the turbulence transport equations. 5. Iterate until convergence.

The algorithm under-relaxes pressure and velocity corrections to maintain stability — typical under-relaxation factors are 0.3 for pressure and 0.7 for momentum. The momentum equation solved is:

$$ \frac{\partial}{\partial x_j}(\rho U_i U_j) = -\frac{\partial p}{\partial x_i} + \frac{\partial}{\partial x_j}\left(\mu_{eff}\frac{\partial U_i}{\partial x_j}\right) + \rho \mathbf{g} $$

where $\mu_{eff} = \mu + \mu_t$ includes both molecular and turbulent viscosity.

**Buoyancy treatment:** For natural and mixed convection, simpleFoam supports body force terms in the momentum equation. The Boussinesq approximation treats density as constant everywhere except in the buoyancy term:

$$ \rho = \rho_0(1 - \beta(T - T_0)) $$

where $\beta$ is the thermal expansion coefficient. This is accurate for temperature differences where $\beta\Delta T \ll 1$ — roughly $\Delta T < 30\text{-}40\text{K}$ in air.

**Steady Reference Frame (SRF):** Rotating machinery can be modeled by specifying a rotating reference frame. The momentum equation gains Coriolis and centrifugal source terms:

$$ \mathbf{S}_{Cor} = -2\rho\,\boldsymbol{\omega}\times\mathbf{U} $$

$$ \mathbf{S}_{cent} = \rho\,\boldsymbol{\omega}\times(\boldsymbol{\omega}\times\mathbf{r}) $$

These are added as fvVectorMatrix source terms. The SRF approach assumes the flow relative to the rotating frame is steady — a reasonable approximation for fully filled passages far from blade tips.

Transient Solvers: pimpleFoam and pisoFoam

**pimpleFoam** combines the PIMPLE algorithm (a blend of SIMPLE and PISO), which allows the SIMPLE outer correction loop to operate with transient time advancement. The PIMPLE loop structure:

The outer corrector loop permits larger Courant numbers ($Co > 1$, typically 5-50), at the cost of an iterative inner correction sequence per time step. Each outer iteration modifies pressure and velocity to approach divergence-free continuity. This is ideal for RANS and LES simulations where time step size is dictated by physics (streak shedding, vortex formation) instead of numerical stability.

**pisoFoam** follows the pure PISO algorithm without outer corrector loops. Each time step consists of:

1. Momentum predictor. 2. Pressure equation. 3. Velocity correction. 4. HSE (Helmholtz separation equation) corrector for transient terms: $\frac{\partial \phi}{\partial t} + \nabla\cdot(\mathbf{U}\phi)$.

Because there are no outer correctors, pisoFoam requires small time steps with CFL < 1 (typically 0.3-0.5) to maintain accuracy. It excels in laminar transient flows and as a benchmark tool. The canonical test case for pisoFoam is the driven cavity flow (Lid-driven cavity) at $Re = 100$, where exact or highly-resolved numerical solutions exist for comparison.

**Pressure-Velocity Coupling :: Detailed analysis of the SIMPLE, PISO, and PIMPLE algorithms, their stability regions, and when to select each.

Turbulence Modeling in Incompressible Flow

The incompressible solvers employ the Boussinesq eddy-viscosity hypothesis, relating the Reynolds stresses to the mean strain rate:

$$ -\overline{u_i'u_j'} = \nu_t\left(\frac{\partial U_i}{\partial x_j} + \frac{\partial U_j}{\partial x_i}\right) - \frac{2}{3}k\delta_{ij} $$

where $\nu_t = k/\omega$ (Spalart-Allmaras, k-omega) or $\nu_t = C_\mu k^2/\epsilon$ (k-epsilon).

Initialization of turbulence quantities is critical:

Here $I$ is the turbulence intensity (typically 1%-10% for engineering), and $l_{turb}$ is the turbulence length scale (typically 5%-7% of the characteristic length). Poor initialization leads to initial transients that can delay convergence by orders of magnitude.

**Turbulence Models :: Comparative survey of RANS models,LES, and hybrid approaches, with guidance on model selection.

**Numerical Schemes :: Discretization schemes for convection, diffusion, and time derivatives in incompressible solvers.

Low-Mach Number Limit

The incompressible assumption breaks down when the Mach number exceeds approximately 0.3. At this point, density variations become significant enough to affect the pressure-velocity coupling structure itself. Pressure waves no longer propagate instantaneously, and the acoustic time scale becomes relevant to the flow physics. The transition criterion is:

$$ Ma = \frac{U}{a} > 0.3 $$

where $a = \sqrt{\gamma R T}$ is the speed of sound. When this threshold is crossed, the solver should switch to a compressible formulation (see Compressible Fluid Module).

Comparison with ANSYS Fluent

Fluent's segregated solver offers more options for turbulence model selection, near-wall treatment (Enhanced Wall Treatment, All y+), and multiphase coupling, while OpenFOAM's incompressible solvers benefit from complete transparency, extensibility, and the ability to modify the discretization schemes at compile time.

References

• Versteeg, H. & Malalasekera, W. (2007). An Introduction to Computational Fluid Dynamics: The Finite Volume Method (2nd ed.). Pearson. Chapters 5-7 cover SIMPLE, PISO, and pressure-velocity coupling.
• Jasak, H. (1996). Error Analysis and Estimation for the Finite Volume Method with Applications to Fluid Flow. PhD Thesis, Imperial College London.
• Weller, H.G., Tabor, G., Jasak, H. & Fureby, C. (1998). A tensorial approach to computational continuum mechanics using object-oriented techniques. Computers in Physics, 12(6), 620-631.
• OpenFOAM Foundation. (2023). OpenFOAM v2212/v14 User Guide.
• Ferziger, J.H. & Peric, M. (2002). Computational Methods for Fluid Dynamics (3rd ed.). Springer.