🔢 mathematics

An $H^{-1}$ least-squares UnCut FEM on domains defined by a level set function

This paper proposes a novel UnCut finite element method for solving Poisson and Stokes equations on level-set-defined domains that utilizes an $H^{-1}$ least-squares formulation to ensure stability without user-tuned penalty parameters while achieving optimal-order convergence.

Original authors: Jiashun Hu, Buyang Li, Han Yang

Published 2026-08-04

📖 4 min read🧠 Deep dive

CC0 1.0

Original authors: Jiashun Hu, Buyang Li, Han Yang

Original paper dedicated to the public domain under CC0 1.0 (http://creativecommons.org/publicdomain/zero/1.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are an architect trying to build a model of a futuristic city, but the city is built on a perfectly smooth, curved hill, while your building blocks are all rigid, square cubes. This is the daily struggle of scientists who use computers to simulate how fluids flow or how heat spreads through objects with curved edges, like a blood vessel or a turbine blade. The standard way to do this is to chop up the space into a grid of tiny squares (a mesh) and try to make the squares fit the curve perfectly. But if the curve is wiggly, the squares get sliced into weird, tiny shards. These shards are a computational nightmare: they are hard to calculate, can make the math unstable, and require complex, error-prone tricks to handle. It's like trying to pave a winding garden path with square bricks; you end up with a lot of tiny, jagged pieces that don't fit well and might collapse the whole path if you aren't careful.

To solve this, scientists developed "unfitted" methods, where they just lay the square grid over the curve and ignore the fact that some bricks are cut in half. The challenge is making the math work on these cut bricks without the simulation blowing up. A recent method called $\phi$ -FEM was a big step forward because it avoided the messy math of the cut pieces, but it still required a "tuning knob" (a stabilization parameter) that had to be set just right. If you turned the knob too low, the simulation crashed; too high, and the results became inaccurate. Finding the perfect setting was like trying to balance a broom on your finger while blindfolded—you had to guess, and it was frustrating.

This paper introduces a new, smarter way to handle these curved boundaries, called the H $^{-1}$ least-squares UnCut FEM. Think of this new method as a self-correcting robot that doesn't need you to fiddle with a tuning knob at all. Instead of guessing how to stabilize the cut bricks, the new method uses a clever "least-squares" strategy. Imagine you are trying to balance a stack of books on a wobbly table. The old way was to add heavy weights (stabilization) to the bottom, but you had to guess exactly how heavy they needed to be. The new method is like having a smart sensor that automatically adjusts the balance by measuring how "off" the stack is and correcting it instantly, no matter how wobbly the table gets.

The authors, Jiashun Hu, Buyang Li, and Han Yang, applied this idea to two classic physics problems: the Poisson equation (which models things like heat distribution or electric potential) and the Stokes equations (which model slow-moving fluids like honey or blood). They proved mathematically that their new method is rock-solid stable without needing any user-tuned parameters. They showed that the solutions converge to the correct answer at the best possible speed (optimal-order convergence) for both problems. In their computer experiments, they tested the method on circles and spheres in both 2D and 3D. The results were impressive: the method worked perfectly even when the grid was cut in very strange ways, and unlike the older method, it didn't care if you changed the "stabilization weight" (they tested it with a weight of 1, and it worked great).

The paper explicitly argues against the idea that you must use large, carefully tuned penalty parameters to keep these simulations stable. They demonstrate that by minimizing the "residual" (the error) in a specific mathematical way (the H $^{-1}$ norm), they can achieve stability with fixed, simple weights. This means the method is "parameter-free" in practice, removing the guesswork that has plagued engineers and scientists for years. While the paper proves the math works and shows it in simulations, it notes that future work is needed to extend these ideas to other types of boundaries and time-dependent problems. But for now, this new approach offers a robust, simpler, and more reliable way to simulate the curved world using a grid of squares.

Technical Summary: An H−1 Least-Squares UnCut FEM on Domains Defined by a Level Set Function

Problem Statement
The paper addresses the numerical solution of the Poisson and Stokes equations on smooth domains $\Omega \subset \mathbb{R}^d$ with curved boundaries, where the domain is implicitly defined by a level set function $\phi$ (i.e., $\Omega = \{x : \phi(x) < 0\}$ ). A primary challenge in solving these problems on such domains is the "small cut cell" issue inherent in unfitted finite element methods (FEMs), where the boundary intersects background mesh elements arbitrarily, creating small or highly anisotropic sub-elements. While existing methods like CutFEM and XFEM address stability through ghost penalties or cell aggregation, they often require accurate numerical integration over curved cut cells, complicating implementation. Furthermore, traditional stabilized unfitted methods, such as the recently proposed $\phi$ -FEM, rely on stabilization parameters that must be chosen sufficiently large to ensure stability, yet the precise minimal values are often unknown in practice, requiring user tuning.

Methodology
The authors propose a novel UnCut finite element method that avoids numerical integration over cut subregions while introducing a least-squares formulation minimizing an $H^{-1}$ -residual of the governing equations.

Domain Extension and Function Spaces: Similar to $\phi$ -FEM, the method extends the computational domain to $\Omega_h$ , consisting of all full background elements intersecting the physical domain. The approximate solution is constructed as the product of a finite element function and the level set function $\phi_h$ , ensuring the Dirichlet boundary condition is satisfied on the discrete interface $\Gamma_h = \{\phi_h = 0\}$ .
Least-Squares Formulation: Instead of standard Galerkin formulations, the method minimizes a functional involving the squared discrete $H^{-1}(\Omega_h)$ $H^{- 1} (Ω_{h})$ -norm of the residual. This is achieved by introducing an auxiliary variable $w_h$ $w_{h}$ defined via the discrete Riesz map $(I - \Delta_h)^{-1}$ $(I - Δ_{h})^{- 1}$ , where $\Delta_h$ $Δ_{h}$ is the discrete Laplacian with homogeneous Neumann boundary conditions.
- For the Poisson equation, the scheme seeks $(w_h, u_h)$ such that the residual of the extended operator is minimized in the $H^1$ norm of the auxiliary variable, augmented by least-squares terms on cut elements and jump penalties on edges.
- For the Stokes equations, the formulation similarly minimizes the $H^{-1}$ -residual of the momentum equation and includes terms for the divergence constraint and least-squares stabilization on cut elements.
Stabilization: The method employs penalty terms on cut elements ( $T \in \mathcal{T}_h^\Gamma$ ) and edges ( $F \in \mathcal{F}_h^\Gamma$ ) to ensure stability. Crucially, the least-squares structure allows these stabilization terms to be effective with fixed unit weights, eliminating the need for large, user-tuned penalty parameters.

Key Contributions

Parameter-Free Robustness: The primary contribution is a formulation that ensures stability without requiring stabilization parameters to be "sufficiently large." The method remains robust with fixed unit weights, addressing a practical limitation of $\phi$ -FEM where optimal parameter values are often a priori unknown.
Avoidance of Curved Integration: Like $\phi$ -FEM, the method avoids complex numerical integration over curved cut cells by performing all volume integrations on full background elements.
Theoretical Analysis: The authors rigorously establish optimal-order convergence rates:
- For the Poisson equation, optimal $H^1$ -norm convergence is proven.
- For the Stokes equations, optimal $H^1$ -norm convergence for velocity and optimal $L^2$ -norm convergence for pressure (up to a constant) are derived.
Novel Stability Mechanism: The paper demonstrates that the intrinsic stability structure of the least-squares UnCut FEM allows stabilization terms to be strong enough for the proof without adversely affecting pressure stability in the Stokes problem, a common difficulty in unfitted discretizations.

Results

Theoretical: Under mild assumptions on the level set function and mesh regularity, the existence and uniqueness of the discrete solution are proven. Error estimates are derived in terms of the mesh size $h$ and the approximation error of the level set function $\|\phi - \phi_h\|_{H^1}$ .
Numerical: Experiments on 2D and 3D Poisson and Stokes problems using FEniCS confirm the theoretical convergence rates.
- The method achieves optimal convergence in both $H^1$ and $L^2$ norms for Poisson problems.
- For Stokes problems, the velocity converges optimally in $H^1$ and pressure in $L^2$ .
- Robustness Comparison: A direct comparison with the original $\phi$ -FEM shows that the proposed method maintains accuracy across a wide range of stabilization parameters $\sigma$ . In contrast, the original $\phi$ -FEM exhibits significant error growth when $\sigma$ is small or large, particularly for the pressure variable.

Significance and Claims
The paper claims that the proposed $H^{-1}$ least-squares UnCut FEM offers a significant improvement in the robustness and simplicity of unfitted finite element methods. By leveraging the least-squares formulation to minimize the $H^{-1}$ -residual, the method successfully eliminates the dependency on large, user-tuned stabilization parameters while retaining the implementation advantages of the $\phi$ -FEM (i.e., avoiding integration on cut cells). The authors assert that this approach provides a parameter-free scheme that is robust with respect to the cut configuration, maintaining accuracy even when mesh elements are arbitrarily small or irregularly cut. The work establishes a rigorous theoretical foundation for these claims, providing optimal-order error estimates for both Poisson and Stokes problems.

The authors note that future research directions include deriving optimal-order $L^2$ -error estimates theoretically, extending the method to Neumann/Robin boundary conditions, and applying it to time-dependent problems.

More like this