- The paper presents a novel unfitted finite element framework using cubic spline boundary tracking to optimize PDE-constrained problems with dynamic domain evolution.
- It employs ghost-penalized CutFEM discretization that avoids remeshing and provides robust convergence and optimal error bounds.
- Numerical experiments validate the method's efficiency in transforming elliptical and complex boundaries into optimal configurations with quantifiable error estimates.
Unfitted Finite Element Methods for PDE-Constrained Shape Optimization via Shape Gradient Flow
The paper focuses on PDE-constrained shape optimization, particularly on the minimization of shape-dependent functionals where the admissible set is comprised of domain boundaries, and the state variable is governed by an elliptic PDE. The authors adopt the framework of shape gradient flow, formulating the problem as a coupled system of PDEs: the state, adjoint, and velocity equations, together with a flow map tracking the time evolution of the domain boundary via the shape gradient. The core minimization prototype is: Γ=∂ΩminJ(Γ)=∫Ωj(⋅,u)dx
subject to an elliptic state equation (with, e.g., Neumann or Dirichlet conditions). The shape calculus approach is used to derive first-order shape gradients, which may be in boundary or distributed form.
Distinct from prior works relying on evolving or body-fitted meshes, this study embeds the moving, evolving domain within a fixed background mesh. The boundary is tracked by a cubic spline representation, whose control points are updated according to the velocity field solution to the regularized shape gradient equation. The PDEs on the physical domain are discretized using a ghost-penalized CutFEM framework, which robustly handles the cut elements resulting from the evolving, unfitted domain.
CutFEM Discretization, Cubic Spline Boundary Evolution, and Algorithm Design
The CutFEM implementation is characterized by:
- The physical domain at each timestep is embedded in a fixed background Cartesian grid.
- At each optimization step, the control points for the boundary spline are advected according to the computed velocity, giving the next domain boundary.
- For each new domain, the elliptic state/adjoint/velocity equations are solved with a CutFEM discretization, stabilized by a ghost-penalty term for faces in the boundary zone.
The update procedure is thus: solve the discrete state and adjoint PDEs (CutFEM, ghost penalty), evaluate the shape gradient (using either boundary or distributed forms), solve for the H1-regularized shape velocity field, advect the spline control points, and repeat. The regularization is critical for stability and for controlling boundary smoothness, yielding physically meaningful deformations.
A notable algorithmic detail is the addition or removal of control points on the spline boundary during the evolution, required to maintain shape quality and parametrization regularity—an advantage of the spline-based approach over level-set or parametric representations with fixed discretizations.
Mathematical Analysis: Error Bounds and Convergence
Theoretical analysis is provided with detailed focus on error estimates for the discrete solutions to the state, adjoint, and velocity fields, and for the flow map representing the evolving domain. The main result (Theorem 4.1) demonstrates that, under standard assumptions (smooth initial data, sufficient regularity of the exact solution, and mesh-step choice τ=O(hk), η=O(hk+21)), the error for the interpolated flow, state, adjoint, and velocity variables satisfies: ∥ ϕn∥H1(Ω0)+∥eun∥H1(Ωhn)+∥epn∥H1(Ωhn)+∥wn∥H1(Ωhn)≤C(τ+hk)
for BDF-1 time stepping, cubic or higher-degree elements (k≥2), and cubic spline geometries. The estimates are proven using a combination of geometrically sharp change-of-domain (pullback) lemmas, stability inequalities for unfitted CutFEM, and detailed treatment of material derivatives through the flow map.
Strong numerical evidence is provided for the optimality of these convergence rates.

Figure 1: Snapshots illustrate the evolution of the domain boundary governed by the shape gradient flow, starting from an elliptic initial guess and converging towards the optimal circular configuration.
Numerical Results
Two numerical experiments demonstrate the performance and convergence properties:
- Elliptical-to-circular test: The method recovers the unit circle from an initial ellipse, with the volume error and H1 errors for the state/velocity variables aligning with the theoretically predicted rates (Table 1).
- Complex star-shaped boundary: Starting from a more complicated geometrical form, the method drives the shape evolution towards the known optimum with similarly strong numerical convergence observed (Table 2).
In both cases, the method successfully avoids remeshing, robustly handles large deformations, and exhibits error scaling as h2 for the spline boundary and state variables for quadratic elements.

Figure 2: Evolution of a nontrivial, non-star-shaped domain to the optimal configuration, demonstrating the robustness of the CutFEM-unfitted approach.
Implications and Perspectives
This work constitutes a significant development at the intersection of shape optimization, unfitted FE discretization, and computational PDE-constrained optimization:
- The method achieves optimal convergence for evolved boundaries without complex remeshing or substantial loss in accuracy, even under substantial shape evolution.
- Stability and accuracy are maintained via the ghost-penalty CutFEM and H1 shape gradient regularization, which is crucial for practical algorithms due to the low regularity of shape gradient fields in PDE-constrained optimization.
- The spline-based boundary representation admits dynamic control point management, conferring an implementation advantage and flexibility in adapting to high curvature or topological changes.
Practically, this framework enables shape optimization for elliptic or parabolic PDE-constrained problems where repeated remeshing is prohibitive or where solution-driven boundary motion is required. Theoretical implications include a robust foundation for error analysis of coupled PDE-ODE moving boundary problems and for further development of adaptive strategies (automatic control point insertion/removal, mesh refinement) for long-term shape evolution.
Conclusion
The proposed unfitted finite element method, with cubic spline boundary tracking and CutFEM for moving domain PDEs, provides a mathematically rigorous and computationally robust tool for PDE-constrained shape optimization. Theoretical results are quantitatively supported by detailed error bounds that are numerically realized in canonical tests. Future work will extend to three-dimensional domains, highly nonlinear PDE constraints, and adaptivity in control point and mesh management for large-deformation or topology-changing shape flows.