Mixed Finite Element Method Overview

• Mixed Finite Element Method is a numerical technique that simultaneously approximates key field variables (e.g., flux and potential) of PDEs.
• It employs H(div) and H(curl) conforming spaces to directly impose conservation laws and maintain stability through inf-sup conditions.
• MFEMs are applied in elasticity, porous media flow, and multiphysics simulations, offering robust, locally conservative approximations confirmed by theoretical and numerical studies.

The mixed finite element method (MFEM) is a class of finite element discretizations for partial differential equations (PDEs) where the solution is sought in terms of multiple, physically distinct field variables (e.g., flux and potential, stress and displacement, pressure and velocity), frequently formulated as a saddle-point system. MFEMs are characterized by their simultaneous approximation of primal and dual fields using finite element spaces that conform to the natural functional context of the PDE—most notably, H(div) and H(curl) conforming spaces—in contrast to standard finite element methods that only approximate scalar or vector fields in L² or H¹ spaces. The approach allows direct imposition of conservation laws and often delivers locally conservative approximations of physically relevant quantities such as fluxes or stresses.

1. Mixed Variational Formulation and Theoretical Foundations

At the heart of MFEMs is the mixed (or saddle-point) variational formulation. For instance, for second-order elliptic problems such as Darcy flow or Poisson's equation, the system is recast so that both the scalar potential $u$ and the flux $q$ are primary unknowns. Formally, for a domain $\Omega$ with appropriate boundary conditions and a symmetric, uniformly positive definite tensor $a(x)$, the continuous mixed weak form reads: find $q \in H(\text{div}; \Omega)$, $u \in L^2(\Omega)$ such that

$$(a q, v) - (\nabla \cdot v, u) = \langle g, v \cdot n \rangle_{\partial\Omega}\quad \forall v \in H(\text{div}; \Omega)$$

$$(\nabla \cdot q, w) = (f, w)\quad \forall w \in L^2(\Omega)$$

where $g$ models boundary data and $f$ is a source term. The choice of Sobolev spaces governing flux and potential is dictated by the structure of the PDE and conservation principles. The inf-sup (Babuška-Brezzi) stability condition is essential for well-posedness and convergence; MFEMs traditionally inherit stability through careful selection of function spaces and discrete versions of the inf-sup condition.

2. Discrete Spaces, Element Choices, and Quadrature

Discrete MFEMs require finite element spaces compatible with the variational formulation. Common choices:

• H(div)-conforming spaces: Raviart-Thomas (RT), Brezzi-Douglas-Marini (BDM), and variants for velocity/flux or stress.
• Scalar potential spaces: $P_k$ discontinuous or continuous Lagrange spaces.
• Hybrid and multipoint variants: Weak Galerkin (WG-MFEM) methods use discontinuous piecewise polynomials for both interior and “boundary” (facet/edge-based) traces, allowing highly flexible meshes, including arbitrary polygons and polyhedra (Wang et al., 2012).
• Rotation/skew-symmetric spaces: For elasticity, auxiliary skew-symmetric fields are discretized on P₀ or P₁ spaces (Ambartsumyan et al., 2018).

Quadrature techniques can be used to localize interactions:

• Vertex/element mass lumping: E.g., vertex-based quadrature in MSMFE methods localizes stress and rotation DOFs to vertices, enabling cell-centered schemes and local elimination of auxiliary variables (Ambartsumyan et al., 2018).

Discretization can be carried out on general simplicial, quadrilateral, or polygonal meshes. On quadrilaterals, the lowest-order BDM₁ elements are adapted along with suitable quadrature rules (vertex or trapezoidal) to facilitate localized algebraic systems and guarantee optimal approximation properties (Ambartsumyan et al., 2018).

3. Algebraic Structure, Local Elimination, and Hybridization

MFEMs yield block-structured algebraic systems typical of saddle-point problems. Schematically, with field vectors $u_h$, $q_h$, the system can be written as: $$\begin{pmatrix} A & B^T \ B & 0 \end{pmatrix} \begin{pmatrix} q_h \ u_h \end{pmatrix} = \begin{pmatrix} g \ f \end{pmatrix}$$ Here, $A$ is the mass/stiffness operator for the flux variable and $B$ the coupling/divergence operator. For problems involving more than two fields (e.g. elasticity with weak symmetry), additional blocks enter for auxiliary (e.g., rotation) variables (Ambartsumyan et al., 2018).

A salient feature is the possibility of local elimination of flux or auxiliary variables, yielding cell-centered, symmetric positive-definite systems for scalar unknowns (Ambartsumyan et al., 2018, Ambartsumyan et al., 2018):

• MSMFE-0 uses localized quadrature for stress, eliminating stress unknowns around each vertex.
• MSMFE-1 extends quadrature to asymmetry bilinear forms, enabling further elimination of rotation variables; often the final system involves only cell-centered displacements.

Hybridization and static condensation approaches, common particularly in weak Galerkin mixed FEMs, move the continuity constraints (normal flux across element boundaries) from being strongly enforced to weakly enforced via face-based Lagrange multipliers. The global system is then assembled for the multipliers, typically yielding a symmetric positive-definite Schur complement of significantly reduced size (Mu et al., 2015). After solution, local post-processing recovers flux and potential variables per element.

4. Stability, Convergence, and Superconvergence Properties

MFEMs are predicated on the satisfaction of discrete inf-sup conditions, crucial for avoiding spurious modes and guaranteeing stable convergence. For example:

• MSMFE methods use macroelement inf-sup arguments and mass-lumped quadrature to demonstrate coercivity and discrete stability (Ambartsumyan et al., 2018, Ambartsumyan et al., 2018).
• WG-MFEM and its hybridizations achieve optimal order convergence and inf-sup stability on general meshes via stabilization terms and weak divergence (Wang et al., 2012, Mu et al., 2015).

Convergence rates:

• First order in natural norms: $$||q - q_h||_{H(\text{div})} + ||u - u_h||_{L^2} \leq C h$$ for flux and potential in MSMFE schemes on simplicial/quadrilateral meshes, and in WG-MFEM (Ambartsumyan et al., 2018, Ambartsumyan et al., 2018, Wang et al., 2012).
• Higher order and superconvergence: Second-order superconvergence $||u - R u_h||_{L^2} \leq C h^2$ for cell-centered displacement in MSMFE, and $O(h^{k+2})$ for WG-MFEM multiplier approximations (Ambartsumyan et al., 2018, Mu et al., 2015).
• Hybridization often yields superconvergence of the face-based multipliers and facilitates optimal error control.

5. Numerical Applications and Practical Aspects

MFEMs provide robust discretization strategies for diverse problems:

• Elasticity: MSMFE methods apply to linear elasticity equations with weakly enforced stress symmetry, delivering locking-free performance in the nearly incompressible regime (Poisson ratios close to 0.5) and resolving discontinuities in heterogeneous materials (Ambartsumyan et al., 2018, Ambartsumyan et al., 2018).
• Darcy and Porous Media Flows: MFEMs are a preferred framework for locally conservative approximations of flow in heterogeneous porous materials. The vertex-quadrature MSMFE approach is inspirational for Darcy flow multipoint methods (Ambartsumyan et al., 2018).
• General Elliptic Problems: WG-MFEM and its hybridized variants offer high flexibility and local conservation for second-order elliptic PDEs on arbitrary meshes, including nonconvex domains and curved boundaries, albeit requiring careful boundary treatments for optimal accuracy on curved domains (Wang et al., 2012, Liu et al., 2022).

Numerical experiments routinely confirm theoretical rates:

• First-order convergence in all variables and second-order superconvergence at cell centers (MSMFE methods).
• Robust performance in near-incompressible and high-contrast media.
• Efficiency gains via local elimination and cell-centered formulations, yielding reduced algebraic complexity.

6. Extensions and Modern Developments

Recent advances in MFEM research extend the methodology to:

• Multiphysics coupling: MFEMs form the basis for simulation frameworks integrating elasticity, poroelasticity, and fluid-solid interactions.
• High-order, adaptive, and multilevel methods: Hierarchical enrichment, stabilization via bubble functions, and adaptive algorithms enhance the flexibility and accuracy even on unstructured or distorted meshes.
• Hybridization frameworks: The hybridized WG-MFEM formulation reduces global system size and enhances parallelizability, notably in the context of complex geometries and nonconforming meshes.
• Fast solvers and preconditioners: The cell-centered, SPD structure emerging from multipoint/hybridized MFEMs is highly amenable to scalable block-diagonal or Schur-complement-based iterative solvers.

7. Summary Table: MSMFE and WG-MFEM Key Attributes

Scheme	Discrete Spaces	Quadrature/Elimination	Convergence
MSMFE-0	BDM₁ (stress), P₀ (rot), P₀ (disp)	Vertex quadrature, σ elimination	O(h) (all vars), O(h²) u
MSMFE-1	BDM₁ (stress), P₁ (rot), P₀ (disp)	Quadrature in asymmetry bilinear form	O(h) (all vars), O(h²) u
WG-MFEM	Pₖ(T) flux, Pₖ(e) flux-bdry, Pₖ₊₁(T)	Element-local weak divergence, hybrid	O(h^{k+1}), O(h^{k+2}) u

• "A multipoint stress mixed finite element method for elasticity on simplicial grids" (Ambartsumyan et al., 2018)
• "A multipoint stress mixed finite element method for elasticity on quadrilateral grids" (Ambartsumyan et al., 2018)
• "A Weak Galerkin Mixed Finite Element Method for Second-Order Elliptic Problems" (Wang et al., 2012)
• "A Hybridized Formulation for the Weak Galerkin Mixed Finite Element Method" (Mu et al., 2015)
• "A Weak Galerkin Mixed Finite Element Method for second order elliptic equations on 2D Curved Domains" (Liu et al., 2022)

The mixed finite element paradigm continues to be developed for advanced PDE models across mechanical, hydrological, and multiphysics contexts, with strong theoretical backing and extensive computational validation.