Constraints¶

Once a model has assembled its matrix \(A\) and load \(\mathbf{b}\), the constraints are imposed before the system is handed to a solver. They come in two kinds: single-point constraints, which fix one degree of freedom to a value (Dirichlet or essential conditions), and multi-point constraints, which tie several degrees of freedom together (hanging nodes from non-conforming refinement, periodic boundaries).

Single-point (Dirichlet) constraints¶

A Dirichlet condition prescribes degrees of freedom directly. Group the unknowns into the free set \(f\) and the constrained set \(c\), so the assembled system reads

\[ \begin{pmatrix} A_{ff} & A_{fc} \\ A_{cf} & A_{cc} \end{pmatrix} \begin{pmatrix} \mathbf{x}_f \\ \mathbf{x}_c \end{pmatrix} = \begin{pmatrix} \mathbf{b}_f \\ \mathbf{b}_c \end{pmatrix} , \qquad \mathbf{x}_c = \mathbf{g} , \]

where \(\mathbf{g}\) is the vector of prescribed values imposed on the constrained degrees of freedom, for example a fixed potential or temperature on a boundary.

μfem (through MFEM) does not drop the constrained block to solve a smaller system. It keeps the full size and eliminates the constrained degrees of freedom symmetrically in place,

\[ \begin{pmatrix} A_{ff} & 0 \\ 0 & I \end{pmatrix} \begin{pmatrix} \mathbf{x}_f \\ \mathbf{x}_c \end{pmatrix} = \begin{pmatrix} \mathbf{b}_f - A_{fc}\,\mathbf{g} \\ \mathbf{g} \end{pmatrix} . \]

Reading this block by block:

The free equations keep \(A_{ff}\) and pick up the correction \(-A_{fc}\,\mathbf{g}\) on the right-hand side. The coupling block \(A_{fc}\) is then zeroed, because the constrained values it multiplied are now known and have moved to the right-hand side.
The constrained rows are overwritten by the trivial equation \(\mathbf{x}_c = \mathbf{g}\). Their original entries \(A_{cf}\), \(A_{cc}\) and \(\mathbf{b}_c\) are discarded, so \(A_{cf}\) is not used further. For a symmetric matrix \(A_{cf} = A_{fc}^{\mathsf{T}}\), so nothing is lost: that same coupling already reaches the free equations through the \(-A_{fc}\,\mathbf{g}\) term.

Zeroing the coupling column \(A_{fc}\), rather than leaving it in place, is what keeps a symmetric \(A\) symmetric, which matters for solvers such as conjugate gradient.

Multi-point constraints¶

Hanging nodes and periodic boundaries instead tie degrees of freedom to one another, for example a slave DOF to a master, \(x_{\text{slave}} = x_{\text{master}}\) (or a weighted combination). These are written as a projection \(P\) that reconstructs the full degree-of-freedom vector from its independent (master) part, filling in each slave as its master combination. Applying it symmetrically transforms both sides of the system:

\[ P^{\mathsf{T}} A\, P\;\mathbf{x} = P^{\mathsf{T}} \mathbf{b} . \]

The operator becomes \(P^{\mathsf{T}} A\, P\) and the right-hand side becomes \(P^{\mathsf{T}} \mathbf{b}\). Each eliminated slave row is given a unit diagonal so the operator stays non-singular.

The constrained system¶

Imposing the constraints does not just modify \(A\); it also produces the two vectors the solver actually receives: a constrained solution vector \(\mathbf{X}\) and a constrained right-hand side \(\mathbf{B}\).

\(\mathbf{X}\) already carries the prescribed values on the constrained degrees of freedom, and holds the current guess (or zero) on the free ones. It doubles as the starting vector for an iterative solve.
\(\mathbf{B}\) is the load after the essential contributions have been folded into it, exactly the \(\mathbf{b}_f - A_{fc}\,\mathbf{g}\) and \(P^{\mathsf{T}}\mathbf{b}\) terms above.

The solver then computes the free entries of \(\mathbf{X}\) from the constrained operator acting on \(\mathbf{X} = \mathbf{B}\). Afterwards the result is mapped back onto the full finite-element field: the eliminated slave and essential degrees of freedom are reinserted from their masters and prescribed values, so the recovered solution satisfies every constraint exactly. How the constrained system is solved, directly or iteratively, is covered in Linear Solvers.