Skip to content

Summary of Governing Equations and Numerical Implementation

The Theory of Porous Media provides a coupled description of saturated soils, accounting for both the solid skeleton and the pore fluid. Based on the previous derivations, the mechanical and hydraulic behaviour is described by two primary balance equations:


Overview of Final Balance Equations

The governing equations of the system are:

  • Balance of linear momentum (mixture):

    \[ \rho^{tot} \mathbf{b} + \nabla \cdot \left( \boldsymbol{\sigma} - \bar{p}^w \boldsymbol{\delta} \right) - \rho^s \bar{\mathbf{a}}^s - \rho^w \bar{\mathbf{a}}^w = \mathbf{0} \]
  • Mass balance of pore water:

    \[ n (\nabla \cdot \bar{\mathbf{v}}^s) + \frac{n}{\bar{K}^w} \frac{\mathrm{d}^s \bar{p}^w}{\mathrm{d}t} + \frac{1}{\bar{K}^w} \nabla \bar{p}^w \cdot \mathbf{w}^w + \nabla \cdot \mathbf{w}^w = 0 \]
  • Generalised Darcy law (from balance of linear momentum of the pore water):

    \[ \mathbf{w}^w = \frac{k^w}{\bar{\mu}^w} \mathbf{K}^s \left( -\nabla \bar{p}^w + \bar{\rho}^w (\mathbf{b} - \mathbf{a}^w) \right) \]
  • Evolution equation of porosity (from mass balance of the solid):

    \[ \frac{\mathrm{d} n}{\mathrm{d}t} = (1 - n) \frac{\mathrm{d} \varepsilon^v}{\mathrm{d}t} \]

These four equations are strongly coupled through the porosity \(n\), solid deformation \(\bar{\mathbf{v}}^s\), Darcy velocity \(\mathbf{w}^w\), and pore pressure \(\bar{p}^w\).


Finite element formulation

For finite element implementation, the problem is reformulated using the following primary field variables:

  • Displacement of the solid skeleton: \(\mathbf{u}\),
  • Pore water pressure: \(p^w\).

These fields are discretised at the finite element nodes as shown in Figure 1.


Figure 1: Stable and unstable two-phase finite elements for saturated soils discretising the solid displacement and the pore water pressure at the nodes.

Their simultaneous solution leads to a mixed formulation, as the displacement and pressure fields are governed by different physical laws (momentum and mass conservation, respectively) and generally require different interpolation spaces.

Mixed formulation

A mixed formulation arises when solving a system where the primary variables belong to different function spaces. This is typical in incompressible or nearly incompressible problems (as in saturated soils), where:

  • The displacement field \(\mathbf{u}\) requires \(C^0\)-continuous shape functions,
  • The pressure field \(p^w\) is governed by a constraint (e.g. incompressibility or mass conservation), which introduces a saddle-point structure into the discrete system.

This structure must satisfy the Ladyzhenskaya–Babuška–Brezzi (LBB) condition (also known as the inf-sup condition) to ensure numerical stability. Violating this condition may lead to non-physical oscillations in the pressure field.

Finite element choices and stabilisation

To obtain a stable and convergent solution, the following strategies are commonly used:

  • Taylor–Hood elements: Use quadratic interpolation for displacement and linear interpolation for pressure (e.g. \(P_2\)\(P_1\) elements), which satisfy the inf-sup condition.
  • Stabilised mixed methods: Additional stabilisation terms (e.g. pressure projection, bubble enrichment, residual-based methods) allow equal-order interpolation (e.g. \(P_1\)\(P_1\)) while maintaining stability.

A detailed derivation of mixed formulations, stability criteria, and stabilisation techniques goes beyond the scope of this lecture. Interested students are referred to the lecture Stabilized Finite Element Methods for Computational Fluid Dynamics offered by the Institute of Mechanics at TU Darmstadt.


Boundary Conditions

To complete the problem formulation, appropriate boundary conditions must be specified for both mechanical and hydraulic fields.

Boundary conditions fall into two categories:

  • Essential (Dirichlet) boundary conditions: where the primary field variable is prescribed directly,
  • Natural (Neumann) boundary conditions: where the associated flux or force is specified. These emerge naturally from the weak form of the governing equations.

Mechanical boundary conditions (from momentum balance):

  • Prescribed displacements: \(\mathbf{u} = \bar{\mathbf{u}}\) on the boundary \(\Gamma_u\) (Dirichlet),
  • Prescribed tractions (total stress): \(\boldsymbol{\sigma}^{tot} \cdot \mathbf{n} = \bar{\mathbf{t}}\) on the boundary \(\Gamma_t\) (Neumann),

with:

\[ \bar{\mathbf{t}} = \left( \boldsymbol{\sigma} - \bar{p}^w \boldsymbol{\delta} \right) \cdot \mathbf{n}. \]

This expression is crucial: the total external traction includes both effective stress and pore water pressure. The portion of the load carried by the solid skeleton depends directly on the local pore pressure. This formulation enables the simulation of phenomena such as consolidation and undrained response.

Hydraulic boundary conditions (from water mass balance):

  • Prescribed pore pressure: \(p^w = \bar{p}^w\) on \(\Gamma_p\) (Dirichlet),
  • Prescribed fluid flux: \(\mathbf{w}^w \cdot \mathbf{n} = \bar{q}\) on \(\Gamma_q\) (Neumann),

These conditions control fluid exchange across the boundary, including:

  • Impermeable surfaces (\(\bar{q} = 0\)),
  • Infiltration boundaries (\(\bar{q} > 0\)),
  • Drained boundaries (\(p^w = 0\)).

Constitutive Models and Material Parameters

To close the system of balance equations, appropriate constitutive relationships must be introduced for both the mechanical and hydraulic behaviour.

Effective stress–strain relationship (solid skeleton):

  • Elastic or elasto-plastic laws, e.g.:
  • Parameters for such models are typically:
    • Elastic modulus, Poisson's ratio,
    • Friction angle, cohesion, dilation angle,
    • Hardening or softening parameters.

Hydraulic permeability or conductivity:

  • May be isotropic or anisotropic,
  • May evolve with deformation or porosity (e.g. Kozeny–Carman-type relationships).

Compressibility of pore water:

  • Represented by the bulk modulus \(\bar{K}^w\),
  • Important for pore pressure dissipation and storage,
  • Often assumed constant for saturated water (\(\bar{K}^w \approx 2.2\) GPa), but may decrease significantly if air is entrapped (even in small amounts).

These models govern how the porous medium responds to loading and fluid flow. Their selection and calibration are essential for accurate prediction of coupled hydro-mechanical processes.