PaperlayerView exact revision on arXiv

Source-linked AI summary

A hybridization-derived interior penalty Nitsche method for heterogeneous and anisotropic diffusion

Vincent Fontaine

arXiv:2610.29214v1math.NAcs.NA

TL;DR

Heterogeneous, anisotropic diffusion can involve arbitrarily large coefficient contrast, motivating a method and analysis for that setting. The paper derives an H1-conforming IPN method by closed-form elimination of a stabilized skeletal multiplier and combines boundary and interior penalties. Its strongly imposed boundary comparison finds that transmissibility weighting retains near-optimal energy-convergence rates over tested meshes and contrasts, while removing diffusion scaling delays convergence.

  • Problem

    The paper addresses heterogeneous, anisotropic diffusion with potentially arbitrarily large contrast, for which robust analysis is a motivation.

    Source: §2.1 · p. 3

  • Method

    The method condenses a stabilized skeletal multiplier in closed form to obtain an H1-conforming formulation with Nitsche boundary and transmissibility-weighted interior flux-jump penalties.

    Source: §1 · pp. 2–3

  • Results

    Over tested meshes and contrasts, transmissibility weighting retains energy-convergence rates close to optimal order, while removing diffusion scaling delays observed energy convergence.

    Source: §6 · p. 21

  • Takeaways & Limitations

    The closed-form numerical normal flux is single-valued across the mesh skeleton.

    Source: §1 · p. 2

    For general discontinuous κ, piecewise elliptic regularity may fail and its constant is generally not contrast-uniform, so Theorem 4.2 does not apply in that regime.

    Source: §4.2 · p. 16

Abstract

from arXiv · show

We develop an interior penalty Nitsche finite element formulation for heterogeneous and anisotropic diffusion problems, which is derived from an interior-penalty stabilized Lagrange multiplier method. A distinctive feature of the derived method is that it simultaneously combines two popular stabilization strategies: (i) a boundary penalty mechanism similar to Nitsche's methods and (ii) an interior penalty mechanism as for Continuous Interior Penalty (CIP) methods. Furthermore, the normal component of the numerical flux is single-valued across the mesh skeleton and available in closed form. We establish the equivalence between the hybridized Lagrange multiplier formulation and its condensed Interior Penalty Nitsche (IPN) form, and we deliver the stability and convergence analysis of the latter. For sufficiently regular solutions, we establish energy-norm error estimates for the three symmetry variants. Under a piecewise elliptic regularity assumption on the adjoint problem, the $L^2$ estimate is optimal for the symmetric variant and exhibits a half-order loss for the incomplete and nonsymmetric variants; its constant may depend on the diffusion coefficients through the adjoint regularity constant. This single-valued numerical normal flux, available in closed form, provides a starting point for locally conservative $H(\mathrm{div})$-conforming flux reconstruction, left for future work. Numerical experiments validate the theoretical findings.

1. Introduction

The paper derives an H1-conforming IPN method by eliminating a stabilized skeletal Lagrange multiplier in closed form. The resulting method combines boundary and interior penalties, with a single-valued numerical flux and stability, error, and numerical analyses.

1. Introduction

  • The construction extends stabilized multiplier and hybridization ideas, while avoiding the delicate inf-sup-compatible pairing required by unstabilized multiplier methods.The article proceeds through the formulation, analysis, numerical experiments, and concluding perspectives, including locally conservative H(div)-flux correction.Source: §1 · pp. 1–3
  • Unlike EDG, the method uses a continuous H1-conforming primal field and a discontinuous face-wise broken skeletal trace; unlike HHO, it eliminates the trace and retains a global primal problem.The boundary and interior penalty mechanisms are jointly inherited from the reduced hybrid structure.Source: §1 · p. 2
  • The IPN method combines weak Dirichlet enforcement through a Nitsche boundary penalty with interior penalization of normal-flux jumps, while retaining the unknowns of continuous Galerkin.It is obtained by algebraically removing skeletal multipliers using numerical traces expressed solely in terms of the primal variable.Source: §1 · p. 2
  • Closed-form, face-local elimination yields a single-valued skeletal normal flux and a transmissibility-weighted flux-jump family in an H1-conforming formulation.The work analyzes stability and a priori energy- and L2-norm errors for three symmetry variants and studies contrast robustness and transmissibility weighting numerically.Source: §1 · pp. 2–3

2. Preliminaries

The preliminaries specify a heterogeneous, anisotropic diffusion problem and its conforming mesh-based finite element setting. They also introduce trace operators that support the hybridized and primal formulations.

2.1. Model problem

  • The model is diffusion with homogeneous Dirichlet data, a symmetric piecewise-constant tensor that may jump across material interfaces, and potentially arbitrarily large contrast.The weak solution is unique in H1_0(Ω).Source: §2.1 · p. 3

2.2. Discrete setting

  • The discrete setting uses shape-regular conforming simplicial meshes, a continuous H1-conforming primal polynomial space, and a single-valued face-wise broken trace space.The mesh skeleton distinguishes interior from boundary interfaces, and the trace space enforces homogeneous Dirichlet data skeletally.Source: §2.2 · pp. 3–4

2.3. Discrete trace operators

  • The trace operators define DG jumps and weighted means, plus an HDG local jump measuring the discrepancy between the primal field and skeletal trace.H1 conformity makes the primal trace single-valued across interior faces, so admissible weighted averages there coincide with that trace.Source: §2.3 · pp. 4–5

3. The interior penalty Nitsche method

The method is built by hybridizing a Nitsche formulation, recovering its traces in closed form, and condensing the skeletal unknown to obtain the primal IPN scheme. The reduction yields a single-valued numerical normal flux and couples boundary and interior penalty mechanisms.

3.1. The hybridizable Nitsche formulation

  • The auxiliary hybridizable Nitsche formulation uses the primal field and a skeletal trace, with a symmetrization parameter controlling the adjoint-consistency term.The bilinear form is symmetric when ϵ = 1 and nonsymmetric otherwise.Source: §3.1 · pp. 5–6
  • Positive stabilization makes each face-local skeletal problem uniquely solvable; facewise constant stabilization additionally gives a pointwise closed-form elimination.The explicit elimination avoids assembling or inverting a skeletal system.Source: §3.1 · p. 6

3.2. Inherited numerical traces

  • The recovered numerical trace and flux depend only on the primal solution, and the numerical normal flux is single-valued across the mesh skeleton.This follows from the discrete flux-transmission requirement.Source: §3.2 · pp. 6–7

3.3. Closed-form elimination: the IPN method

  • Substituting the closed-form trace reduces the hybridized problem to the primal IPN method, analytically realizing the skeletal Schur complement.The reduction requires strictly positive, facewise constant stabilization and avoids explicit skeletal-block assembly and inversion.Source: §3.3 · pp. 7–8
  • The IPN family has symmetric, nonsymmetric, and incomplete variants; its interior term penalizes normal diffusive-flux jumps, extending the local CG stencil to face-neighbours.The nonsymmetric variant supplies the positively stabilized flux-jump term, while the symmetric variant has a sign-reversed term; the incomplete variant has no interior mechanism.Source: §3.3 · pp. 8–9, Figure 1 · p. 9
  • For the strongly imposed boundary variant, the family comprises CG, positively stabilized, and sign-reversed flux-jump schemes; the three are not identical to the weak-boundary solutions.For weakly imposed Dirichlet data, single-valuedness gives global conservation: total outward numerical boundary flux equals the integrated source; that argument does not apply to the strongly imposed variant.Source: §3.3 · pp. 9–10

4. Stability and a priori error analysis

The section analyzes stability and a priori energy- and L2-error estimates for the IPN method, identifying penalty conditions and the roles of coefficient regularity and symmetry. It establishes an optimal L2 estimate for all three symmetry variants under piecewise elliptic regularity and explains how trace elimination couples interior and boundary stabilization.

4.1. Stability analysis

  • The three IPN variants are coercive and yield a unique discrete solution under sufficient penalty conditions; N-IPN needs no positive stability threshold, while S-IPN on simplices requires α0 > 2(d + 1).The stability constant is independent of mesh size, polynomial degree, and diffusion tensor; the incomplete variant has the sufficient threshold α0 > η0.Source: §4.1 · p. 12

4.2. Convergence analysis of the native discrete variable

  • Normal continuity of the exact diffusive flux makes the interior stabilization residual vanish, yielding strong consistency and Galerkin orthogonality.The interior mechanism contributes no residual to the Galerkin orthogonality relation.Source: §4.2 · p. 13
  • Under the stated piecewise regularity of the exact solution, the energy-norm error estimate is independent of coefficient jumps and local anisotropy ratios, with no quasi-monotonicity assumption.The estimate’s displayed weights and broken regularity seminorm retain dependence on the diffusion problem.Source: §4.2 · p. 14
  • With piecewise elliptic regularity of the adjoint, the L2-error rate is optimal for ϵ = 1 and loses at most half an order for other variants.The hidden constant depends on the piecewise elliptic-regularity constant, which may depend on diffusion coefficients through contrast.Source: §4.2 · p. 15
  • For general discontinuous diffusion coefficients, the adjoint regularity assumption may fail and its constant may not be contrast-uniform, so the L2 estimate need not apply.In such cases, piecewise adjoint regularity may be only H1+γ, with γ potentially arbitrarily small as contrast increases; when γ ≤ 1/2, the required L2 gradient traces need different treatment.Source: §4.2 · p. 16
  • The term h in (50) arises from weak, nonconforming Dirichlet enforcement through ηh, and vanishes when uh is sought in the conforming space Vh,0 ⊂ H1_0(Ω).In the conforming variant, the corresponding optimal estimate holds for all three ϵ variants under assumption (49).Source: §4.2 · p. 16
  • Under piecewise elliptic regularity assumption (49), all three variants ϵ ∈ {0, ±1} satisfy the optimal estimate ∥u − uh∥0,Ω ≲ κ0 h^µ+1 |u|µ+1,Th.The estimate applies to the conforming variant as well as the symmetric, incomplete, and nonsymmetric variants.Source: §4.2 · p. 16
  • Trace elimination propagates ϵ to both interior flux-jump and boundary Nitsche terms: ϵ = −1 gives unconditional coercivity but loses adjoint consistency, while ϵ = 1 ensures adjoint consistency and optimal L2 error under the stated regularity, with a stronger penalty threshold.This coupling is inherited from the auxiliary hybrid formulation.Source: §4.2 · p. 16

5. Numerical experiments

The experiments test convergence and contrast sensitivity on a smooth anisotropic checkerboard and the singular Kellogg problem, then isolate the effect of transmissibility weighting in the conforming variant.

5.1. Test 1: heterogeneous anisotropic diffusion with a smooth solution

  • For k = 1–4, the three variants have energy errors within 0.25%, while their L2 curves are compatible with order k + 1 and S-IPN has the smallest errors.The L2 observations are empirical because the adjoint regularity assumption is not established for this checkerboard.Source: §5.1 · p. 17, Figure 3 · p. 18
  • Across contrasts up to 10^12, relative energy errors stay within [1.5, 2.4]×10^-4 for all three variants.At h = 1/64 and k = 2, the errors are essentially flat across the tested range.Source: §5.1 · p. 17, Figure 4 · p. 18

5.2. Test 2: heterogeneous diffusion with a singular solution (Kellogg)

  • In the Kellogg problem, the energy rate is close to γ ≈ 0.5354 across variants, while the empirical L2 rate is about 1.1 and condition (49) is violated.Higher polynomial degrees lower the error constant without changing the limiting energy slope.Source: §5.2 · pp. 18–19, Table 1 · p. 19

5.3. The conforming variant: convergence and the role of the transmissibility weighting

  • At λ = 10^3, transmissibility weighting gives observed orders 1.96 in energy and 3.39 in L2; the L2 observation is empirical.The checkerboard adjoint regularity assumption required for the theoretical L2 estimate is not established.Source: §5.3 · p. 20
  • At contrast 10^8, contrast-blind weighting yields 2.56 times the finest-mesh energy error; transmissibility weighting maintains observed energy orders of 1.93–2.00.The table compares the ratio of contrast-blind to transmissibility-weighted energy errors across contrasts and mesh sizes.Source: §5.3 · p. 20, Table 2 · p. 21

6. Conclusion

The paper derives and analyzes a condensed IPN method for heterogeneous, anisotropic diffusion, combining boundary and transmissibility-weighted flux-jump penalties. Transmissibility weighting retains near-optimal observed energy rates, while flux reconstruction and L2 conclusions remain bounded by stated analysis requirements.

6. Conclusion

  • The H1-conforming IPN method is obtained by condensing a stabilized whole-skeleton Lagrange multiplier formulation, and its analysis establishes equivalence, stability, and energy- and L2-norm error estimates for three symmetry variants.Its energy estimate uses elementwise interpolation and requires no quasi-monotonicity of the diffusion coefficients.Source: §6 · p. 21
  • The L2 benchmark results are empirical because the required piecewise adjoint regularity is not established for the smooth checkerboard and fails for the Kellogg configuration.The stated L2 estimate depends on that regularity assumption.Source: §6 · p. 21
  • Across tested meshes and contrasts, transmissibility weighting retained energy convergence rates close to optimal, while removing diffusion scaling delayed observed convergence.The comparison concerns the strongly imposed boundary variant.Source: §6 · p. 21
  • The closed-form single-valued normal flux is a starting point for locally conservative H(div)-conforming reconstruction, whose local equilibration and estimator robustness require further analysis.A guaranteed Prager–Synge a posteriori estimator is also left for future work.Source: §6 · p. 21

Reproducibility

The article’s Zenodo archive provides the code and reference data accompanying the paper, along with materials for reproducing the numerical results.

Reproducibility

  • The Zenodo archive includes Python/NGSolve implementations, benchmark drivers, figure-generation scripts, and reproduction instructions.It also contains the article’s code and reference data.Source: Reproducibility · p. 21
Loading 2610.29214v1…