A method for bounding high-order finite element functions: Applications to mesh validity and bounds-preserving limiters

TL;DR

Bounding extrema of high-order finite element polynomials is challenging due to loose global bounds and ill-conditioned basis transformations. The authors propose a robust bounding-box framework: for each basis function $\phi_i$, compute optimal bounding functions $L^i(\mathbf{x})$ and $U^i(\mathbf{x})$ so that $L^i(\mathbf{x})\leq \phi_i(\mathbf{x})\leq U^i(\mathbf{x})$; then bound any $u_h(\mathbf{x})=\sum_i u_i\phi_i(\mathbf{x})$ by combining these envelopes and refining via a low-order projection to form $u_{LO}$ and high-order fluctuations $u_h'$. The method yields local bounds with $O(M^{-2})$ convergence in 1D and extends naturally to tensor-product and simplex elements, enabling on-the-fly, order-agnostic bounds. It is demonstrated in mesh validity checks (ensuring $|\mathbf{J}|>0$ across high-order curved elements) and continuously bounds-preserving limiters for hyperbolic systems, where it provides tighter, subcell-aware bounds than Bernstein-based approaches and reduces unnecessary refinements. Compared with traditional convex-hull methods, the proposed bounding boxes offer tighter, more localized bounds, improving robustness and efficiency in high-order simulations and related applications in graphics and multi-physics simulations.

Abstract

We introduce a novel method for bounding high-order multi-dimensional polynomials in finite element approximations. The method involves precomputing optimal piecewise-linear bounding boxes for polynomial basis functions, which can then be used to locally bound any combination of these basis functions. This approach can be applied to any element/basis type at any approximation order, can provide local (i.e., subcell) extremum bounds to a desired level of accuracy, and can be evaluated efficiently on-the-fly in simulations. Furthermore, we show that this approach generally yields more accurate bounds in comparison to traditional methods based on convex hull properties (e.g., Bernstein polynomials). The efficacy of this technique is shown in applications such as mesh validity checks and optimization for high-order curved meshes, where positivity of the element Jacobian determinant can be ensured throughout the entire element, and continuously bounds-preserving limiters for hyperbolic systems, which can enforce maximum principle bounds across the entire solution polynomial.

Figures (27)

• Figure 1: Example of $C^0$ upper/lower linear bounding functions for an arbitrary polynomial $u_h(x)$ defined by control nodes $\boldsymbol{\eta}$ and control values $\mathbf{q}$.
• Figure 1: Predicted extrema for a polynomial interpolating a step function $[-0.5, 0.5]$ on the Gauss--Lobatto nodes at varying approximation orders. Error reduction in the bounds between the present approach (with $M = N$) and Bernstein approach shown on bottom.
• Figure 2: Example of the bounding approach directly applied to a high-order polynomial (left), a visualization of the projection step (middle), and the bounding approach with the projection step (right).
• Figure 2: Minimum/maximum values of the solution (computed via brute force sampling) for the solid body rotation problem at $t=1$ as computed by a $\mathbb P_3$ DG approximation with global maximum principle bounds $u_h(\mathbf{x},t) \in [0,1]$ on a varying number of quadrilateral mesh elements $N_e$.
• Figure 3: Examples of optimal bounding boxes for $N=4$ Gauss--Lobatto (left), Gauss--Legendre (middle), and Bernstein (right) basis functions with $M=5$ equispaced control nodes.