Locking-Free Mixed Finite Element Discretization

Updated 6 March 2026

Locking-free mixed finite element discretization is a class of methods that mitigate volumetric locking by introducing auxiliary variables such as stress or pressure.
These methods employ tailored finite element spaces and inf-sup stable formulations, ensuring error estimates remain robust despite degenerating physical parameters.
They are widely applied to problems in elasticity, plate and shell theories, and poroelasticity, maintaining uniform convergence in high-contrast, nearly-incompressible scenarios.

Locking-free mixed finite element discretization denotes the class of mixed variational and finite element methods for elasticity (and related PDEs) whose stability, approximation, and convergence properties are provably uniform with respect to critical physical limits—most notably, the nearly-incompressible limit (volumetric locking, λ → ∞ in isotropic elasticity or t → 0 in plates). These methods guarantee that error constants, norm equivalences, and convergence rates are independent of degenerating parameters such as high contrast in coefficients or vanishing physical constants. This property is achieved through judicious choice of mixed variational formulations, stable saddle-point discretizations (fulfilling the Ladyzhenskaya–Babuška–Brezzi (LBB) or inf-sup conditions), and the use of specialized finite element spaces, often involving stress or pressure as independent variables in addition to displacement.

1. Mixed Variational Formulations and Locking Phenomena

Classical locking arises when single-field primal (displacement) finite element methods for linear elasticity, Reissner–Mindlin plates, or incompressible nonlinear elasticity fail to approximate the limiting constraint subspace as a material parameter degenerates (e.g., Poisson's ratio ν → ½ or plate thickness t → 0). This manifests in under-convergent or stagnant error rates for displacement and stress, spurious pressure oscillations, and stiff numerical response.

Mixed variational formulations, in which auxiliary variables—stress, pressure, rotation, or shear—are introduced as Lagrange multipliers or independent fields, render the problem into a constrained (saddle-point) structure. For linear elasticity in stress-displacement form, the continuous mixed system on a domain Ω ⊂ ℝⁿ (n = 2, 3) with isotropic elasticity tensor C (Lamé parameters λ, μ) is:

Equilibrium: ∇·σ = f in Ω
Constitutive: σ = C : ε(u)
Symmetry: σ symmetric
Kinematics: ε(u) = (∇u + ∇uᵗ)/2

with natural boundary conditions, and the mixed weak form seeks (σ, u) in appropriate H(div; S)-conforming stress and L²-displacement spaces, ensuring constraints are enforced at the discrete level (Chung et al., 25 Apr 2025, Carstensen et al., 14 Mar 2025).

For Reissner–Mindlin plate and shell models, mixed forms introduce the shear rotation and shear strain as auxiliary fields, restoring uniformity in the limit t → 0 (Ainsworth et al., 27 Jun 2025, Xia, 2015). In incompressible nonlinear elasticity, pressure appears as a Lagrange multiplier for the pointwise isochoric constraint, leading to a two-field or three-field saddle-point structure (Huang et al., 2017, Huang et al., 2018).

2. Inf-sup Stability and Discrete De Rham Complexes

Locking-free mixed discretizations are constructed to satisfy inf-sup stability (LBB condition) uniformly with respect to the critical parameters. For the mixed stress-displacement method, the finite element spaces Σ_h (for stress) and U_h (for displacement) are chosen such that

$\inf_{v_h\in U_h}\sup_{\tau_h\in\Sigma_h} \frac{( \operatorname{div}\,\tau_h, v_h )}{\|\tau_h\|_{H(\operatorname{div})}\|v_h\|_{L^2}} \geq \beta > 0$

with β independent of h and λ (Chung et al., 25 Apr 2025). This is achieved by using, for instance, the classical (stable) Arnold–Falk–Winther or Brezzi–Douglas–Marini families for symmetric stress, or related normal-normal continuous spaces (Carstensen et al., 14 Mar 2025). Typical pairings are Pk (symmetric tensors) for stress and P_{k–1} for displacement, ensuring exact or weak stress symmetry and div-conformity as needed.

The framework of the de Rham complex further guides the construction of locking-free elements for plates and shells. For Reissner–Mindlin plates, the stability depends on discrete complexes that exactly replicate the kernel structure of the continuous model, including harmonic forms arising in multiply-connected or mixed boundary condition domains (Ainsworth et al., 27 Jun 2025). A discrete reduction operator and careful tracking of commutative diagrams (nabla, rot, div) yield necessary and sufficient conditions for locking-free performance in all topologies.

For methods involving pressure or rotation as independent fields (Stokes, Biot, or nonlinear elasticity), the finite element pairs must be Stokes-stable, e.g., DP-Q2–P1, stabilized MINI, or iso-parametric Q2–P1 combinations, and support a uniform Fortin operator or an equivalent projection (Huang et al., 2017, Huang et al., 2018, Zhao et al., 2020).

3. Multiscale, Virtual, and Hybrid Mixed Schemes

Locking-free mixed finite element principles are compatible with advanced discretization concepts, including:

Multiscale Hybrid-Mixed (MHM): Combines coarse skeleton face unknowns (for tractions or fluxes) with local fine-scale solves. Stability and absence of locking are ensured by careful stabilization (e.g., least-squares terms), appropriate face polynomial degree choices, and uniform Fortin-like projections. The convergence and error constants remain λ-independent due to preservation of coarse-to-fine coupling under nearly-incompressible limits and high-contrast material coefficients (Gomes et al., 2024, Chung et al., 25 Apr 2025).
Virtual Element Methods (VEM): Generalize mixed conformity to arbitrary polygonal/polyhedral meshes, constructing virtual spaces via projections and stabilization terms. Locking-free estimates are achieved by ensuring k-consistency and uniform inf-sup via Fortin-type arguments, even when small or distorted edges are present (Lovadina et al., 5 Jan 2026, Leppe et al., 28 Jan 2026).
Primal Hybrid/Multi-field/Hybridizable Discontinuous Galerkin (HDG): Primal-hybrid-pressure methods couple broken displacement, pressure, and skeleton Lagrange multipliers using mesh-dependent norms and local equilibrium solvers. A postprocessing step recovers H(div)-conforming, locally equilibrated, weakly symmetric stress fields, with convergence constants independent of the incompressibility parameter (Taraschi et al., 29 Jan 2026).
Staggered Discontinuous Galerkin and MSMFE–MFMFE: In fully mixed poroelasticity (Biot) formulations, neither zero storage (c₀ → 0) nor incompressibility pose stability issues, thanks to mixed stress-rotation and velocity-pressure formulation, block-diagonalization, and local elimination procedures (Zhao et al., 2020, Ambartsumyan et al., 2020, Lee, 2014).
Mixed Petrov–Galerkin for Rod, Plate, and Shell Theories: In low-dimensional mechanical models subject to shear/membrane locking (Cosserat rods, Mindlin plates), mixed Hellinger–Reissner or saddle-point forms with independent stress/shear resultants, together with appropriate polynomial interpolation for kinematic and stress fields, guarantee locking-free performance at all aspect ratios and slenderness parameters (Herrmann et al., 2 Jul 2025, Xia, 2015).

4. Theoretical Error Analysis and Uniform Convergence

A characteristic feature of locking-free schemes is that stability, coercivity, and best-approximation error estimates remain robust as the critical parameter degenerates (e.g., λ, t, c₀). Typical error bounds for the mixed stress-displacement method (denoting A = C⁻¹):

$\|\sigma - \sigma_h\|_A + \|u - u_h\|_{L^2} \leq C h^s \|f\|_{L^2} \qquad (C\ \text{independent of}\ \lambda).$

For the primal-hybrid-pressure and VEM discretizations, similar estimates hold for the auxiliary variables (pressures, skeleton tractions), and the reconstructed H(div) or virtual stresses recover the sharp order with no volumetric locking (Carstensen et al., 14 Mar 2025, Leppe et al., 28 Jan 2026, Taraschi et al., 29 Jan 2026). In poroelasticity and Biot's consolidation, both semi-discrete and fully discrete error analysis demonstrates first-order convergence rates in all fields, uniformly in time and independent of the storage coefficient c₀→0 (Lee, 2014).

For eigenvalue problems, the convergence rate for frequencies or eigenvalues remains O(h^{2s}), with spectral pollution and all error constants unaffected by mesh geometry or Poisson ratio (Leppe et al., 28 Jan 2026). In nonlinear elasticity, the energy, displacement, determinant, and pressure errors exhibit the theoretical optimal rates—O(h²⁾ for energy, O(h) for displacement and pressure—regardless of the level of incompressibility or singularity due to cavitation (Huang et al., 2017, Huang et al., 2018).

Numerical evidence in all cases supports the absence of any locking effect, even for ν→0.49999 or in materials with coefficient contrast up to 10⁶ (Chung et al., 25 Apr 2025, Lovadina et al., 5 Jan 2026).

5. Implementation Structures, Reduction, and Solver Considerations

The algebraic structure of mixed discrete systems is typically block-saddle-point in stress, displacement, pressure, and potentially rotation or Lagrange multipliers. Efficient solver strategies are built on:

Block elimination or static condensation (eliminating mixed variables locally), yielding reduced symmetric positive-definite systems in primary variables (Ambartsumyan et al., 2020, Carstensen et al., 14 Mar 2025).
Preconditioned MINRES or block-LDL^T solvers (for symmetric indefinite systems) (Carstensen et al., 14 Mar 2025, Taraschi et al., 29 Jan 2026).
Exploitation of block-diagonal or multi-level features in vertex-based or face-based quadrature for hybrid or multiscale methods (Ambartsumyan et al., 2020, Gomes et al., 2024).
For nonlinear locking-free methods in large deformation elasticity, damped Newton solvers with Fortin-projected corrections and orientation-preserving constraints (Huang et al., 2017, Huang et al., 2018).

For VEM discretizations, all scalar products and consistent projections are computable from the degrees of freedom, and stabilization terms are chosen to preserve consistency and avoid spurious kernel modes (Lovadina et al., 5 Jan 2026, Leppe et al., 28 Jan 2026).

For mixed Cosserat rods, block-linearization of the coupled system in centerline, quaternion, and resultant stress/moment unknowns is performed at each Newton-Raphson iteration, ensuring favorable conditioning and singularity-free kinematics (Herrmann et al., 2 Jul 2025).

6. Applications: Heterogeneity, High Contrast, Polygonal and Holey Domains

Locking-free mixed finite element discretizations have demonstrated effectiveness for:

High-contrast, heterogeneous elastic media: Multiscale methods achieve uniform convergence and contrast-robustness by local oversampling and spectral basis truncation independent of k_max/k_min (Chung et al., 25 Apr 2025).
Irregular, polygonal, and non-simply connected meshes: Mixed VEM and MHM approaches naturally accommodate arbitrary geometric topologies, including domains with holes or mixed boundary conditions (Gomes et al., 2024, Ainsworth et al., 27 Jun 2025, Leppe et al., 28 Jan 2026, Lovadina et al., 5 Jan 2026).
Multiphysics and poroelasticity: Locking-free five-field formulations extend to coupled flow/solid systems, preventing pressure oscillations and ensuring correct incompressible and vanishing-storativity asymptotics (Zhao et al., 2020, Ambartsumyan et al., 2020, Lee, 2014).
Thin/Slender structural elements: Mixed methods are impervious to locking arising from thickness or slenderness parameters in Mindlin plates, Reissner-Mindlin shells, and Cosserat rods (Herrmann et al., 2 Jul 2025, Xia, 2015).

7. Summary Table of Representative Locking-free Mixed FE Schemes

Reference	Model/Equation Type	Key Mixed Spaces/Features
(Chung et al., 25 Apr 2025)	Linear elasticity (high contrast, multiscale)	Stress-displacement, Falk–Brezzi–Douglas–Marini pairs
(Carstensen et al., 14 Mar 2025)	Linear elasticity (2D, symmetric stresses)	Normal-normal continuous quadratic stress, split disp.
(Zhao et al., 2020)	Biot poroelasticity (polygonal)	5-fields: stress, disp., rotation, flux, pressure, DG
(Gomes et al., 2024)	Linear elasticity (multiscale, hybrid-mixed)	MHM trace-face unknowns, LS-stabilized Neumann solves
(Huang et al., 2017, Huang et al., 2018)	Incompressible nonlinear elasticity	DP-Q2–P1, Q2–P1, P2+–P1, Fortin operator, damped Newton
(Xia, 2015, Ainsworth et al., 27 Jun 2025)	Mindlin plates, Reissner-Mindlin	Mixed/three-field, reduction operator, MITC/RT/BDM pairs
(Lovadina et al., 5 Jan 2026, Leppe et al., 28 Jan 2026)	Elasticity, contact, eigenproblems	Mixed VEM, k-consistency, polygonal meshes
(Taraschi et al., 29 Jan 2026)	Primal-hybrid, pressure (linear elasticity)	Broken H¹ displacement, edge traction, local stress recov

All these methods exhibit uniform error estimates (no locking) and optimal convergence with respect to mesh size, independent of material degeneracy or geometric singularity.

References: