A high-order nodally bound-preserving and mass-conservative method for linear fourth-order elliptic problems and its applications to nonlinear parabolic equations

Technical Summary: A High-Order Nodally Bound-Preserving and Mass-Conservative Method for Linear Fourth-Order Elliptic Problems and Its Applications

Problem Statement
The paper addresses the numerical simulation of fourth-order partial differential equations (PDEs), which arise in materials science and fluid dynamics (e.g., Cahn–Hilliard models, lubrication-type equations). These physical models possess intrinsic structural properties: pointwise bound-preservation (solutions must remain within specific physical limits, e.g., concentration $u \in [a, b]$), mass-conservation, and energy dissipation. Designing high-order numerical schemes that simultaneously preserve these properties remains a significant challenge. Existing high-order finite element methods often fail to ensure mass-conservation, while low-order methods or post-processing techniques may lack accuracy or rigorous convergence guarantees. Furthermore, for nonlinear parabolic problems, no existing high-order Finite Element Method (FEM) simultaneously preserves bounds, mass, and energy dissipation.

Methodology
The authors propose a unified framework based on variational inequalities to enforce nodal bound-preservation and mass-conservation.

1. Linear Fourth-Order Elliptic Problems:

   • The authors reformulate the linear fourth-order elliptic system as a variational inequality over a closed convex set of admissible functions.
   • They introduce a strictly convex minimization formulation in the $H^{-1}$-$H^1$ norms. This reformulation eliminates the auxiliary variable $w$ from the minimization problem, ensuring well-posedness.
   • A key technical innovation is the construction of a bound-preserving and mass-conservative quasi-interpolation operator. This operator is essential for deriving optimal error estimates.
   • The discrete scheme is solved efficiently using primal–dual active set methods.

2. Nonlinear Fourth-Order Parabolic Problems:

   • The elliptic framework is extended to time-dependent problems by combining the variational inequality with Backward Differentiation Formula (BDF) temporal discretizations (up to high order).
   • For problems with a potential term (e.g., Cahn–Hilliard), the authors integrate the Scalar Auxiliary Variable (SAV) technique. This enhances stability while maintaining the bound and mass constraints.
   • The resulting fully discrete schemes are shown to admit an equivalent convex minimization structure, ensuring existence and uniqueness of the solution.
   • For first-order temporal discretization ($k=1$), a modified energy stability result is established.

3. Nonlinear Second-Order Parabolic Problems:

   • To apply the fourth-order framework to second-order problems (e.g., porous medium equations), the authors introduce a consistent fourth-order regularization strategy.
   • The second-order equation is perturbed by a small fourth-order operator ($\epsilon \Delta^2 u$), where $\epsilon$ scales with the mesh size ($h^{p+1}$).
   • While a complete well-posedness theory for this regularized second-order case is not provided, the scheme is shown to preserve bounds and mass numerically.

Key Contributions

• High-Order Elliptic Schemes: The paper presents the first high-order FEM for linear fourth-order elliptic problems that is simultaneously nodally bound-preserving and mass-conservative. It provides a rigorous $H^1$-error estimate with optimal convergence rates under suitable regularity assumptions.
• Space-Time High-Order Parabolic Schemes: The authors develop space-time high-order schemes for nonlinear fourth-order parabolic problems that preserve bounds, mass, and (for first-order time stepping) modified energy. This fills a gap where no existing high-order FEM simultaneously satisfies all three structural properties.
• Regularization for Second-Order Problems: A practical strategy is proposed to treat nonlinear second-order parabolic problems within the fourth-order framework via consistent regularization, enabling bound and mass preservation where standard high-order FEMs fail.
• Efficient Solvers: The variational inequality formulations are shown to be equivalent to strictly convex minimization problems, allowing for efficient solution via primal–dual active set methods.

Numerical Results
Extensive numerical experiments validate the theoretical findings:

• Accuracy: The schemes achieve optimal convergence rates in $L^2$ and $H^1$ norms for smooth solutions, even on low-quality meshes with obtuse triangles.
• Robustness: The methods handle discontinuous data, singularities (e.g., lubrication-type equations with $\beta=0.5$), and low regularity (e.g., Barenblatt solutions for porous medium equations) without requiring ad-hoc regularization of the mobility or restrictive relaxation of bounds.
• Structure Preservation: Numerical tests confirm strict adherence to nodal bounds and mass conservation. For the Cahn–Hilliard equation with logarithmic potentials, the proposed method (with SAV) avoids non-physical oscillations observed in other bound-preserving methods (e.g., standard Lagrange multipliers) and maintains energy stability.
• Efficiency: The primal–dual active set solver converges rapidly, typically requiring fewer than 5 iterations, even near singularities.

Significance and Claims
The paper claims to address critical gaps in the literature regarding the simultaneous preservation of physical structural properties in high-order numerical methods.

• For fourth-order elliptic problems, it provides a rigorous high-order method with optimal error analysis that preserves both bounds and mass.
• For nonlinear fourth-order parabolic problems, it offers the first high-order FEM framework that simultaneously preserves bounds, mass, and energy dissipation (in the first-order case).
• For nonlinear second-order problems, it offers a practical, high-order strategy to enforce these properties, supported by strong numerical evidence despite the lack of a complete analytical well-posedness theory for the regularized case.

The authors emphasize that their framework is robust for challenging benchmarks involving singularities and low regularity, making it a viable tool for complex physical simulations where structural fidelity is paramount. Future work is identified as extending the error analysis to $L^2$ norms, providing rigorous convergence proofs for the parabolic schemes, and applying the framework to coupled multiphysics models.