Non-Symmetric Nitsche Method

• Non-symmetric Nitsche method is a weak imposition technique that modifies or omits the adjoint symmetry term to enforce boundary and interface conditions without additional multiplier unknowns.
• It is applied in various settings—including finite element, unfitted interface, and domain decomposition—using penalty-based, penalty-free, or super-penalty formulations for enhanced stability.
• Advanced error analyses show that while the loss of adjoint consistency can affect L² convergence, tailored dual problems help restore optimal error estimates.

Searching arXiv for recent and foundational papers on the non-symmetric Nitsche method. The non-symmetric Nitsche method is a weak imposition technique for boundary or interface constraints in which the consistency term is retained but the adjoint or symmetry term is modified, sign-flipped, or omitted, so that the resulting discrete form is generally not symmetric. Across finite element, unfitted/interface, domain decomposition, contact, and boundary element settings, it appears in penalty-based, super-penalty, and penalty-free forms. Its recurring advantages are the elimination of Lagrange multiplier unknowns and, in several formulations, stability without a large penalty parameter; its classical drawback is loss of adjoint consistency and the resulting difficulty of proving optimal $L^2$ estimates by standard duality arguments (Benzaken et al., 2022, Chouly et al., 2011, Chen et al., 7 Oct 2025).

1. Algebraic structure and defining variants

A standard structural decomposition writes a Nitsche discretization as the sum of an interior bilinear form, a boundary or interface consistency term, a symmetry term, and a penalty term. In the abstract variational framework, the symmetric method is obtained by adding both the consistency contribution and its adjoint counterpart, whereas a non-symmetric method is obtained by modifying that last step. The tutorial literature makes this decomposition explicit and notes that the same construction steps can be used to build non-symmetric Nitsche methods (Benzaken et al., 2022).

For the Poisson problem with weak Dirichlet enforcement, a canonical non-symmetric form is

$$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$

with $c_0\ge 0$ and $\alpha\ge 1$. The non-symmetry is concentrated in the boundary pair

$$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$

which is not invariant under interchange of trial and test functions. In the penalty-free case $c_0=0$, this reduces to the classical continuous Galerkin non-symmetric Nitsche form studied without any boundary penalty term (Chen et al., 7 Oct 2025, Burman, 2011).

Different sign conventions coexist in the literature. In the hypersingular boundary element formulation for domain decomposition, the interface bilinear form contains the switch parameter $\sigma\in\{-1,1\}$; $\sigma=1$ gives the symmetric variant, while $\sigma=-1$ gives the non-symmetric or skew-symmetric variant. In that setting, the non-symmetric character is carried by the interface term

$$\frac{\sigma}{2}\langle [u_h],\, T_1 v_h - T_2 v_h\rangle_{\gamma},$$

and functions are allowed to be discontinuous across the interface $$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$0 (Chouly et al., 2011).

2. Stability mechanisms: coercivity, inf–sup, and penalty-free formulations

Non-symmetric Nitsche methods are associated with two distinct stability mechanisms. In some formulations the bilinear form is directly coercive in a mesh-dependent energy norm; in others, stability is recovered through an inf–sup argument rather than through coercivity.

A basic coercive example is the hypersingular boundary element method with $$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$1. There, discrete ellipticity holds for any $$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$2, whereas the symmetric case $$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$3 requires $$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$4. The same contrast appears in fitted and unfitted interface formulations: the non-symmetric choice trades symmetry for robustness with respect to the penalty parameter. In the unfitted Poisson interface method, the mesh-dependent norm

$$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$5

satisfies $$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$6, so the discrete problem is coercive on $$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$7. The coefficient-weighted averages

$$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$8

are used precisely because they yield estimates independent of the coefficient jump (Chouly et al., 2011, Chen et al., 14 Oct 2025).

The penalty-free branch follows a different route. For the continuous Galerkin method without boundary penalty, testing with $$a(u,v) := (\nabla u,\nabla v) - \langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma} + c_0\langle h^{-\alpha}u,v\rangle_{\Gamma},$$9 only gives control of the interior gradient. The key observation is therefore that boundary control can be recovered through an inf–sup condition rather than via a coercivity estimate. In the Poisson case this is achieved by constructing special patch functions $c_0\ge 0$0 with prescribed average normal flux on boundary patches $c_0\ge 0$1, and then testing with $c_0\ge 0$2. Closely related patch constructions appear in penalty-free domain decomposition and in compressible and incompressible elasticity, where the non-symmetric form is shown to satisfy a Banach–Nečas–Babuška type stability estimate even though no term of the form $c_0\ge 0$3 is present on the boundary or interface (Burman, 2011, Boiveau, 2015, Boiveau et al., 2014).

This suggests a genuine bifurcation inside the subject. One branch keeps a penalty and exploits direct coercivity; the other removes the penalty and replaces coercivity by a carefully designed inf–sup argument. The non-symmetric method occurs in both branches.

3. Consistency, adjoint consistency, and error estimation

Primal consistency is a defining feature. In the hypersingular interface problem, the exact solution satisfies $c_0\ge 0$4 and the interface condition $c_0\ge 0$5 on $c_0\ge 0$6, so the jump and interface terms are not required for consistency. The same pattern holds for the Poisson and elasticity formulations: the exact solution satisfies the discrete equations when inserted into the Nitsche bilinear form. What is typically lost in the non-symmetric variant is adjoint consistency (Chouly et al., 2011, Burman, 2011).

That loss dominated the classical $c_0\ge 0$7 theory. For the penalty-free non-symmetric Poisson method, the proved energy estimate is optimal,

$c_0\ge 0$8

but the $c_0\ge 0$9 estimate is only

$\alpha\ge 1$0

The same half-order loss appears in penalty-free domain decomposition for discontinuous material parameters and in penalty-free compressible elasticity. In the boundary element hypersingular setting, the main estimate is “almost quasi-optimal”: $\alpha\ge 1$1 with the logarithmic perturbation attributed to the Nitsche coupling (Burman, 2011, Boiveau, 2015, Boiveau et al., 2014, Chouly et al., 2011).

Recent analyses changed this picture by modifying the dual problem instead of the primal method. For weak Dirichlet imposition, the adjoint bilinear form

$\alpha\ge 1$2

satisfies $\alpha\ge 1$3. By defining the dual problem with $\alpha\ge 1$4 rather than with the standard adjoint PDE, one restores adjoint consistency at the variational level. This yields optimal $\alpha\ge 1$5 convergence,

$\alpha\ge 1$6

for super-penalty, standard-penalty, and penalty-free variants, including the practically important case $\alpha\ge 1$7 on general shape-regular meshes. An analogous tailored dual problem yields optimal $\alpha\ge 1$8 estimates in unfitted interface FEM, together with the optimal energy estimate

$\alpha\ge 1$9

The modern theory therefore no longer treats suboptimal $$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$0 convergence as intrinsic to non-symmetry; rather, it identifies the standard dual problem as the source of the gap (Chen et al., 7 Oct 2025, Chen et al., 14 Oct 2025).

4. Interfaces, non-matching meshes, and unfitted discretizations

The non-symmetric Nitsche method is especially prominent in interface and domain decomposition problems because it weakly enforces continuity without introducing multiplier spaces. In the hypersingular boundary element setting, the surface is decomposed into two polygonal subdomains $$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$1 and $$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$2 with interface $$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$3. The discrete space

$$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$4

contains functions that satisfy the homogeneous boundary condition on $$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$5 but are in general discontinuous across $$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$6. The Nitsche coupling acts directly on the jump $$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$7 and on interface operators $$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$8, so non-matching meshes are allowed and no Lagrange multiplier is required. The paper explicitly lists the resulting advantages: no additional unknown on the interface, no inf–sup condition among interface spaces, elliptic discrete problems, and the possibility of standard linear solvers (Chouly et al., 2011).

The same principle appears in fitted and unfitted Poisson interface problems with discontinuous coefficients. In the fitted penalty-free method, the interface bilinear form is

$$-\langle\nabla u\cdot n,v\rangle_{\Gamma} + \langle u,\nabla v\cdot n\rangle_{\Gamma},$$9

again with no interface penalty term. In the unfitted case the same interface form is combined with a ghost penalty acting on faces near the cut interface, not on $c_0=0$0 itself. This stabilization controls small cut elements while preserving the penalty-free character of the Nitsche coupling on the interface (Boiveau, 2015).

Unfitted interface FEM with piecewise polynomial spaces on $c_0=0$1 and $c_0=0$2 sharpens this picture. There the non-symmetric interface bilinear form is

$c_0=0$3

The special weights $c_0=0$4 are introduced to make the constants independent of the contrast, and a penalty-free variant adds a ghost penalty

$c_0=0$5

The interface literature thus uses non-symmetric Nitsche both with and without explicit penalty on the jump, and in both fitted and cut settings (Chen et al., 14 Oct 2025).

5. Extensions in mechanics and fluid dynamics

In contact mechanics, the non-symmetric Nitsche method is not limited to equality constraints. For the Signorini problem, the contact conditions are rewritten as

$c_0=0$6

and the resulting penalty-free non-symmetric formulation is

$c_0=0$7

The method uses nonconforming piecewise affine elements and achieves optimal error estimates in the energy norm (Burman et al., 2016).

Linear elasticity provides another major branch. For compressible elasticity, the nonsymmetric penalty-free boundary form is

$c_0=0$8

with $c_0=0$9 defined through the normal traction of the stress tensor. For incompressible elasticity, the same nonsymmetric Nitsche boundary treatment is coupled with Galerkin least-squares pressure stabilization, producing locking-free convergence in the mixed formulation (Boiveau et al., 2014).

In incompressible flow, a non-symmetric Nitsche method has been analyzed for the stationary Navier–Stokes equations with slip boundary conditions. There the normal component of the velocity on $\sigma\in\{-1,1\}$0 is enforced weakly through boundary consistency and penalty terms involving $\sigma\in\{-1,1\}$1, while the tangential Navier law remains as a Robin-type term. By contrast, the later Stokes–Poisson–Boltzmann work with Navier slip explicitly uses the symmetric variant of Nitsche’s method, which isolates the non-symmetric method as one choice inside a broader weak-imposition framework rather than as the universal formulation (Bansal et al., 2023, Agrawal et al., 14 Apr 2026).

A general mechanics formulation pushes the same idea further: every Nitsche method is viewed as a stabilized Lagrange multiplier method in which the multiplier has been eliminated. In that framework, equality constraints arise by dropping the positive-part operator in the reduced energy, while inequality constraints retain the term

$\sigma\in\{-1,1\}$2

This viewpoint does not privilege the symmetric method: it identifies symmetric and non-symmetric Nitsche variants as different eliminations of a stabilized saddle-point formulation (Gustafsson et al., 5 Mar 2026).

6. Numerical behavior, parameter choices, and solver implications

The numerical record shows a persistent tension between old theory and observed performance. In the penalty-free Poisson method, the $\sigma\in\{-1,1\}$3-error is essentially insensitive to the penalty parameter $\sigma\in\{-1,1\}$4, and the $\sigma\in\{-1,1\}$5-error is often slightly smaller for moderate positive $\sigma\in\{-1,1\}$6, even though the theory with $\sigma\in\{-1,1\}$7 gives only $\sigma\in\{-1,1\}$8. In fitted and unfitted domain decomposition, the proved rate is again $\sigma\in\{-1,1\}$9 in $\sigma=1$0, but the experiments with $\sigma=1$1 show approximately $\sigma=1$2. The penalty-free elasticity paper reports the same phenomenon: optimal $\sigma=1$3 convergence in theory and numerically optimal $\sigma=1$4 behavior in the tests (Burman, 2011, Boiveau, 2015, Boiveau et al., 2014).

Parameter sensitivity depends strongly on the chosen variant. In the hypersingular boundary element method, the skew-symmetric case $\sigma=1$5 converges for all tested $\sigma=1$6, whereas the symmetric case can fail to converge when $\sigma=1$7 is too small and requires $\sigma=1$8 scaling at least like $\sigma=1$9 in theory. In unfitted interface FEM, the standard non-symmetric penalty term $\sigma=-1$0 is not required to be “very large,” because stability already follows from the weighted flux structure and the choice of $\sigma=-1$1. In penalty-free variants, the practical stabilization task shifts from tuning a boundary penalty to designing the correct inf–sup or ghost-penalty mechanism (Chouly et al., 2011, Chen et al., 14 Oct 2025).

The algebraic consequences are equally consistent across the literature. Symmetric Nitsche preserves symmetry of the discrete operator for symmetric model problems; non-symmetric Nitsche yields non-symmetric matrices. In the hypersingular boundary element paper, the symmetric case admits conjugate gradients once the penalty is large enough, whereas the non-symmetric case calls for Krylov methods for non-symmetric systems such as GMRES or BiCGSTAB. The penalty-free Poisson formulation likewise produces a sparse non-symmetric matrix of the form $\sigma=-1$2, with no extra interface or multiplier unknowns (Chouly et al., 2011, Burman, 2011).

Recent numerical experiments finally align the $\sigma=-1$3 theory with longstanding computations. Three-dimensional tests for the weak Dirichlet problem confirm $\sigma=-1$4-error $\sigma=-1$5 and $\sigma=-1$6-error $\sigma=-1$7 for standard penalty, super-penalty, and penalty-free non-symmetric Nitsche variants. In the unfitted interface setting, the tailored-dual analysis is also reported to match the observed optimal $\sigma=-1$8 rates. This suggests that the modern interpretation of non-symmetric Nitsche is no longer dominated by its historical half-order $\sigma=-1$9 loss, but by a more precise distinction between primal consistency, adjoint consistency, and the choice of dual problem in the analysis (Chen et al., 7 Oct 2025, Chen et al., 14 Oct 2025).