About
In this article, we take advantage of this capability to solve an optimization problem that fully constrains three-dimensional displacements using the mortar method.
This approach offers a key advantage — it enables the treatment of various combinations of bodies without altering their mesh structures.
Algorithm
Nitche Method
Instead of enforcing the contact constraint through a Lagrange multiplier (which introduces extra unknowns), Nitsche’s method weakly imposes the constraint by modifying the weak form with consistent terms.
The key idea is to:
- – Add penalty-like terms to stabilize the interface,
- – Add consistent adjoint terms to maintain symmetry and consistency.
This results in an augmented bilinear form on the contact surface.
In the general form of force equilibrium equation, we will use below weak form
\int_{\Omega} \boldsymbol{\varepsilon}(\boldsymbol{v}) : \boldsymbol{\sigma}(\boldsymbol{u}) , d\Omega = \int_{\Omega} \boldsymbol{v} \cdot \boldsymbol{f} , d\Omega
+
\int_{\Gamma_t} \boldsymbol{v} \cdot \bar{\boldsymbol{t}} , d\GammaAnd we would obtain linear equation from the weak form
\begin{align}
\mathbf{K}u =F_{N}+F_{body}\\
\begin{bmatrix}
\mathbf{K}_1 & \mathbf{0} \\
\mathbf{0} & \mathbf{K}_2
\end{bmatrix}
\begin{bmatrix}
\mathbf{u}_1 \\[4pt]
\mathbf{u}_2
\end{bmatrix}
=
\begin{bmatrix}
\mathbf{f}_1 \\[4pt]
\mathbf{f}_2
\end{bmatrix}
\end{align}In the Nitsche method, a penalty-like matrix is added to the stiffness matrix obtained from the above weak formulation by using bilinear form on the contact surface as it is mentioned above,
\begin{aligned}
a_\Gamma(u_1, u_2; v_1, v_2)
&=\\
\int_{\Gamma}
\Big[
\frac{\alpha}{h} (u_1 - u_2) \cdot (v_1 - v_2)\\
- \frac{1}{2} \big( (\boldsymbol{t}_1(u_1) + \boldsymbol{t}_2(u_2)) \cdot (v_1 - v_2) \big)\\
- \frac{1}{2} \big( (\boldsymbol{t}_1(v_1) + \boldsymbol{t}_2(v_2)) \cdot (u_1 - u_2) \big)
\Big] ds
\end{aligned}
\boldsymbol{t}_i(u_i) = \boldsymbol{C} : \boldsymbol{\varepsilon}(u_i) \boldsymbol{n},
\qquad
\boldsymbol{\varepsilon}(u_i) = \tfrac{1}{2}(\nabla u_i + \nabla u_i^\mathsf{T})\\where
\boldsymbol{C} is the elasticity tensor,
\boldsymbol{n} is the unit normal vector on the interface ( \Gamma ),
h is the characteristic element size,
and \alpha is the penalty parameter.
From this form, we would get K_{Nitsche}.
(K+K_{Nitsche})u=F_{N}+F_{body}The meaning of 1st term is that it prevents jump of displacement on the interface. Additionally, the meaning of 2nd and 3rd terms are kind of energy (like compliance).

Lagrange Method
Formulation
\begin{aligned}
a_c(\mathbf{u}_1, \mathbf{u}_2; \mathbf{v}_\lambda)
&= \int_{\Gamma_c} \mathbf{v}_\lambda \cdot (\mathbf{u}_1 - \mathbf{u}_2) \, \mathrm{d}\Gamma \\[6pt]
&= \int_{\Gamma_c} \mathbf{v}_\lambda \cdot \mathbf{u}_1 \, \mathrm{d}\Gamma
- \int_{\Gamma_c} \mathbf{v}_\lambda \cdot \mathbf{u}_2 \, \mathrm{d}\Gamma
\end{aligned}\mathbf{B}_1 = \int_{\Gamma_c} N_\lambda^{\mathrm{T}} N_1 \, \mathrm{d}\Gamma,
\qquad
\mathbf{B}_2 = - \int_{\Gamma_c} N_\lambda^{\mathrm{T}} N_2 \, \mathrm{d}\Gamma\underbrace{
\begin{bmatrix}
\mathbf{K}_1 & \mathbf{0} & \mathbf{B}_1^{\mathrm{T}} \\
\mathbf{0} & \mathbf{K}_2 & \mathbf{B}_2^{\mathrm{T}} \\
\mathbf{B}_1 & \mathbf{B}_2 & \mathbf{0}
\end{bmatrix}
}_{\text{global stiffness matrix } \mathbf{A}_{\mathrm{aug}}}
\begin{bmatrix}
\mathbf{u}_1 \\[4pt]
\mathbf{u}_2 \\[4pt]
\boldsymbol{\lambda}
\end{bmatrix}
=
\begin{bmatrix}
\mathbf{f}_1 \\[4pt]
\mathbf{f}_2 \\[4pt]
\mathbf{0}
\end{bmatrix}\begin{bmatrix}
\mathbf{K}_1 & \mathbf{0} & \mathbf{B}_1^{\mathrm{T}} \\
\mathbf{0} & \mathbf{K}_2 & \mathbf{B}_2^{\mathrm{T}} \\
\mathbf{B}_1 & \mathbf{B}_2 & \mathbf{\epsilon }\mathbf{M}
\end{bmatrix}
\begin{bmatrix}
\mathbf{u}_1 \\[4pt]
\mathbf{u}_2 \\[4pt]
\boldsymbol{\lambda}
\end{bmatrix}
=
\begin{bmatrix}
\mathbf{f}_1 \\[4pt]
\mathbf{f}_2 \\[4pt]
\mathbf{0}
\end{bmatrix}