A Mixed Extended Virtual Element Method for Elliptic Interface Problems on Polygonal Meshes

Abstract: We propose a lowest-order (H(\operatorname{div}))-conforming mixed extended virtual element method for elliptic interface problems on interface-unfitted polygonal meshes. The flux and pressure are approximated by subdomain-wise extended (H(\operatorname{div}))-VEM spaces and by piecewise constants, respectively. On cut elements, the computable polynomial projection is defined on the whole background element and then restricted to the two subdomains. Compared with BDM-type polynomial spaces, the mixed VEM space contains a non-polynomial component, which gives rise to additional consistency terms on cut elements. To control these terms, we use an enhanced kernel stabilization on cut elements and an interface normal-flux average in the mixed coupling. A corrected interface-flux penalty and a local divergence ghost penalty are added to obtain cut-position-independent stability without using a volume div-div augmentation. We prove continuity, a discrete inf-sup condition, and an optimal first-order error estimate in a mesh-dependent norm. The constants are independent of the mesh size and of the position of the interface relative to the background mesh, but may depend on the coefficient contrast.

• The paper develops a lowest-order H(div)-conforming mixed extended VEM that alleviates meshing constraints for elliptic interface problems on unfitted polygonal meshes.
• It employs polynomial projections and innovative stabilization terms, including ghost penalties, to robustly enforce transmission conditions.
• Numerical tests confirm optimal convergence rates and cut-position robustness even under high coefficient contrasts and complex geometries.

Mixed Extended Virtual Element Method for Elliptic Interface Problems on Polygonal Meshes

Introduction

This work develops a lowest-order $H(\text{div})$-conforming mixed extended Virtual Element Method (VEM) tailored for elliptic interface problems on interface-unfitted polygonal meshes. Elliptic interface problems pervade multiphysics applications, including composite heat conduction, heterogeneous porous media flow, and multiphase electromagnetic transmission, where solutions exhibit reduced regularity and discontinuous coefficients. Conventional interface-fitted FEMs can deliver optimal error estimates but are hampered by mesh-generation complexities for evolving or curved interfaces. Unfitted strategies alleviate meshing constraints but entail challenges surrounding the robust imposition of transmission conditions.

The proposed approach employs fully polynomial-based projections and stabilization mechanisms that ensure robust performance regardless of interface position or coefficient contrast. This method integrates ideas from mixed finite element theory, extended/immersed discretizations, and the virtual element paradigm, leading to a scheme that is both practical for complex geometries and possesses theoretical guarantees for stability and optimal convergence.

Methodological Framework

The key methodological innovation lies in the design of an $H(\text{div})$-conforming mixed extended VEM that leverages background polygonal meshes, irrespective of interface alignment. The unknowns are the pressure, discretized by piecewise constants, and the flux/velocity, discretized via a subdomain-wise extended $H(\text{div})$-VEM space. Notably, the VEM velocity space inherently contains a non-polynomial kernel component, necessitating new consistency and stabilization arguments absent in classical polynomial spaces.

Polynomial Projection and Flux Coupling

On background mesh elements intersected by the interface ("cut cells"), the computable $L^2$-projection onto polynomials is performed on the entire element and then restricted to subdomains demarcated by the interface. The coupling to the pressure exploits the fact that the normal component of the flux on the interface is computable via degrees of freedom, even in the absence of explicit basis functions inside elements.

The mixed formulation utilizes interface normal-flux averages (denoted $M_r V_h$), ensuring computability of flux jumps relevant for transmission conditions. An additional interface penalty term, based on a corrected interface-flux ($F_r$), is introduced to ensure stability independently of the spatial position of the interface within cut cells.

Kernel Stabilization and Ghost Penalties

The non-polynomial kernel component $(I - Q_k) V_h$ generically pollutes consistency on cut elements. To control these artifacts, the authors employ an enhanced kernel stabilization term ($G_{\text{ker}}$) specifically on cut elements. Additional stabilization terms include:

• Local divergence ghost penalty ($G_{\text{div}}$): Augments stability by controlling the divergence norm on extended subdomains, circumventing the need for volume-based div-div augmentations.
• Pressure-jump penalties ($J_1$, $H(\text{div})$0): These terms enforce the weak continuity of pressures across transmission and interface edges, vital for maintaining a discrete inf-sup condition in the mixed setting.

Mesh-Dependent Norms and Scalings

Norms encompassing contributions from physical residuals, kernel penalties, flux penalties, and interface-jump penalties are adopted for the convergence analysis. All stabilization and penalty terms are proved to be robust with respect to mesh size ($H(\text{div})$1) and do not deteriorate as the interface approaches degenerate configurations (small cut fractions).

Analytical Results

Well-Posedness and Stability

The discrete bilinear form is shown to be continuous and satisfies a discrete inf-sup (LBB) condition, with constants independent of $H(\text{div})$2 and the interface position relative to the mesh. The pressure and velocity coupling—through reconstructed interface flux averages and ghost-penalty terms—guarantee stability with respect to both physical and mesh-induced approximation errors.

Error Estimates

Optimal first-order a priori error estimates for the mesh-dependent norm are established, contingent on standard regularity of the exact solution ($H(\text{div})$3, $H(\text{div})$4, $H(\text{div})$5 and continuity of $H(\text{div})$6 across the interface). The analysis accommodates arbitrary polynomial-degree VEM spaces but is herein instantiated for the lowest-order case. Crucially, the error constants do not depend on $H(\text{div})$7 or the relative position of the interface, but may—by necessity—depend on the contrast of the diffusion coefficients across the interface.

Numerical Validation

Comprehensive experiments validate the theoretical findings:

• Convergence: For canonical straight and circular interfaces, and for high-contrast coefficients (e.g., $H(\text{div})$8), first-order convergence in the mixed error and second-order convergence in projected velocity are observed.
• Cut-position Robustness: The method shows negligible sensitivity to interface position, even as the interface approaches element vertices (i.e., as cut fractions become arbitrarily small), with mixed errors remaining stable and unaffected.
• Complex Interfaces and Coefficient Contrasts: On curved and nontrivial interfaces (e.g., heart-shaped), the scheme preserves optimal rates. Tests with a range of coefficient contrasts, including moderate and extreme cases, confirm that the method remains quantitatively robust, corroborating analytical predictions.

Theoretical and Practical Implications

This method fills a critical gap in the discretization of elliptic interface problems on general polygonal meshes without mesh-interface alignment. The employment of virtual elements enables high geometric flexibility, and the stabilization framework ensures that accuracy and stability are not compromised by interface geometry or mesh irregularity. Theoretical implications include:

• Generalizability: The construction serves as a blueprint for higher-order and multi-physics couplings within the VEM paradigm when interfaces are unfitted.
• Divergence and Transmission Robustness: The ghost penalty approach and interface coupling strategy are extensible to other unfitted-mesh formulations requiring strict mass conservation (e.g., mixed Darcy flow, incompressible Stokes interface problems).
• Contrast Dependence: While explicit independence from coefficient contrast is unattainable, the method manages contrast dependencies through explicit interface weighting and stabilization, reflected in the robust numerical performance.

Practically, the approach is amenable to application in computational mechanics, geosciences, and materials science, especially in dynamically evolving interface scenarios or when employing agglomerated grid techniques.

Future Directions

Potential extensions include:

• Adaptive Mesh Refinement: Embedding the method within an adaptive framework exploiting local error indicators, potentially leveraging interface-aware refinement strategies.
• Non-Linear and Time-Dependent Problems: Generalizations to nonlinear material laws or time-dependent interface evolution (e.g., Stefan or phase-field models) are natural extensions.
• Higher-Order and Nonconforming Variants: Systematic study of higher-order mixed extended VEMs and their implications for transmission accuracy.

Conclusion

The mixed extended VEM developed in this work establishes an effective, rigorously analyzed, and practically validated approach for elliptic interface problems on unfitted polygonal meshes. Its ability to deliver optimal convergence, cut-position robustness, and geometric flexibility distinguishes it from existing methodologies. The analytical and numerical results emphasize its applicability as a foundational tool for complex interface problems in mixed variational formulations.