A high order stabilization-free virtual element method for general second-order elliptic eigenvalue problem

High-Order Stabilization-Free Virtual Element Methods for Second-Order Elliptic Eigenvalue Problems

Introduction

This work introduces a high-order stabilization-free virtual element method (VEM) for general second-order elliptic eigenvalue problems. Classical VEMs are characterized by the presence of polynomially-consistent and stabilization terms, with the latter ensuring stability and convergence but at the cost of analytical inconvenience and potentially increased error in some settings. The stabilization-free approach aims to formulate VEMs omitting this stabilization component, achieving a more transparent convergence analysis while capturing accuracy and reducing computational artifacts associated with the stabilization term.

The proposed method targets eigenproblems of the form:

$$\nabla \cdot ( -\mathbf{\mathcal{K}} \nabla u + \boldsymbol{\beta} u ) + \gamma u = \lambda u \quad \text{in}~\Omega, \quad u=0~\text{on}~\partial \Omega,$$

with general, possibly non-symmetric coefficients, and encompasses non-self-adjoint cases. The authors develop both the theoretical foundation (well-posedness, optimal a priori error bounds) and demonstrate practical performance through a diverse set of numerical experiments.

Methodology

The main innovation is the construction of a high-order, stabilization-free VEM discretization, based on the polynomial projection operator framework introduced in prior works for scalar elliptic equations. For a given polygonal decomposition of the domain $\Omega$ into elements $E \in \mathcal{T}_h$, local virtual element spaces of arbitrary order $k \geq 2$ are built via:

• Vertex-based values and moments on element boundaries and interiors define the degrees of freedom.
• The elliptic and $L^2$ projections onto $\mathbb{P}_k(E)$ are computable from these degrees of freedom.
• The key polynomial projection $\Pi_\mathcal{P}^{0,E} \nabla$ projects the gradient onto a tailored polynomial-plus-curl space, replacing the classical need for a stabilization form.

The global bilinear form assembles such projected local terms for the diffusion, convection, and reaction components:

$$\mathcal{B}_h(u_h, v_h) = \sum_E \left( a_h^E(u_h,v_h) + b_h^E(u_h,v_h) + c_h^E(u_h,v_h) \right),$$

where all terms are computable solely from the prescribed degrees of freedom, with no explicit non-polynomial stabilization.

Error Analysis

A detailed a priori error estimation reveals that the method is optimally convergent in both the $H^1$-seminorm and the $L^2$-norm, provided appropriate regularity of the eigenspace. The discrete solution operator is proven to converge to its continuous analogue in norm, and the convergence analysis employs spectral approximation theory for non-self-adjoint compact operators.

Let $\Omega$0 denote an eigenvalue of multiplicity $\Omega$1, and $\Omega$2, $\Omega$3 the associated exact and discrete eigenspaces. The critical estimates are:

• Gap in energy norm: $\Omega$4
• Gap in $\Omega$5 norm: $\Omega$6
• Eigenvalue error: $\Omega$7

These rates are shown to be optimal.

Numerical Experiments

Extensive computations validate the theoretical results. The method is implemented for multiple high-order cases ($\Omega$8) on polygonal meshes of varying complexity, including convex quadrilaterals, pentagons, and concave octagons.

The utilized mesh types are shown below.

Figure 1: Different polygonal discretizations of the unit square.

Case 1: Isotropic Diffusion, Mild Convection

For $\Omega$9, convergence plots for the first five eigenvalues confirm that the stabilization-free VEM achieves optimal rates across all mesh types and polynomial orders.

Figure 2: Error curves for $E \in \mathcal{T}_h$0 over three mesh types, showing optimal convergence for all eigenvalues.

Figure 3: Error curves for $E \in \mathcal{T}_h$1 highlight the persistence of optimal convergence, even for higher polynomial order approximations.

Figure 4: Error curves for $E \in \mathcal{T}_h$2 illustrate the method's efficacy for high-order approximation, with the rare exception of minor stagnation for very high-order eigenvalue approximations on some meshes.

Case 2: Strong Convection

For increased convection ($E \in \mathcal{T}_h$3), results remain robust, further affirming the method's stability in advection-dominated regimes. While initial mesh coarseness can yield sub-optimal convergence for certain eigenpairs, refinement eliminates this discrepancy.

Figure 5: For $E \in \mathcal{T}_h$4, the third eigenvalue on coarse meshes converges sub-optimally but achieves optimal order upon refinement.

Figure 6: For $E \in \mathcal{T}_h$5, optimal convergence recovered across all cases.

Figure 7: $E \in \mathcal{T}_h$6 yields uniformly optimal rates, affirming robustness for higher polynomial order.

Case 3: Strong Anisotropy — SFVEM vs. SVEM

A direct comparison between the stabilization-free VEM (SFVEM) and the standard VEM (SVEM, which employs an explicit stabilization term) is performed for an anisotropic case ($E \in \mathcal{T}_h$7). Both methodologies converge at the optimal rate, but the stabilization-free form consistently outperforms standard VEM in terms of absolute error for a fixed number of degrees of freedom.

Figure 8: SFVEM yields significantly lower errors than SVEM, especially notable for highly anisotropic diffusion.

Implications and Future Directions

The stabilization-free VEM consolidates the advantages of VEMs on polygonal, possibly non-convex meshes with increased accuracy and simplified error analysis in the absence of ad hoc stabilization artifacts. The robust convergence in both isotropic and anisotropic, convection-dominated, and high-order polynomial settings marks its suitability for demanding applications in computational mechanics, wave propagation, and eigenvalue problems involving complex geometries and operator coefficients.

On the theoretical side, the method paves the way for simplified a posteriori error analysis (as the troublesome stabilization term is eliminated), and aids in the development of adaptive strategies, multi-physics couplings, and extension to three-dimensional or mixed formulations. On the practical side, this approach enables more efficient solution processes, with demonstrated competitive or superior accuracy compared to stabilization-based methods for a fixed computational budget.

Potential avenues for future study include rigorous assessment of the method in non-smooth domains, extension to non-linear or parameter-robust eigenproblems, and integration with domain decomposition or multigrid solvers tailored for the stabilization-free VEM framework.

Conclusion

This work provides a systematic study and validation of a high-order stabilization-free virtual element method for general second-order elliptic eigenvalue problems (2604.03727). The approach yields optimal convergence for both eigenvalues and eigenspaces on highly general polygonal meshes, without any need for artificial stabilization. Comparative studies establish the method's superiority in accuracy and efficiency relative to traditional VEMs with added stabilization, reinforcing the practical and theoretical significance of the stabilization-free paradigm for modern PDE eigenvalue computation.