Level-set topology optimisation with unfitted finite elements and automatic shape differentiation

Technical Summary: Level-set Topology Optimisation with Unfitted Finite Elements and Automatic Shape Differentiation

Problem Statement
Topology optimisation is a critical computational technique with broad industrial applications, traditionally relying on density methods or level-set methods. Conventional level-set methods often employ "smudging" techniques, where the domain boundary is smoothed over an interface using a Heaviside function to allow integration over a static background domain. While computationally efficient, this approach can introduce undesirable artefacts in nonlinear or fluid-structure interaction problems and fails to fully capture the shape derivative of functionals under boundary perturbations.

Unfitted finite element methods (such as CutFEM) offer an alternative by embedding complex domains into simple bounding domains without additional degrees of freedom, providing a precise description of the boundary. However, these methods face challenges:

1. Analytic Derivation: Deriving first- and second-order shape derivatives for unfitted discretisations is mathematically difficult and error-prone.
2. Boundary Intersections: Existing analytic shape calculus results (e.g., Berggren [16]) often assume the domain boundary is entirely enclosed within the background domain, excluding cases where the domain boundary intersects the background boundary.
3. Isolated Volumes: In unfitted methods, disconnected regions of the design domain (isolated volumes) can lead to singular system matrices because they are not connected to Dirichlet boundaries via a weak "ersatz" material, as is common in density methods.
4. Precision and Generality: It remains unclear whether automatic differentiation techniques can recover analytic shape derivatives to machine precision in the context of unfitted discretisations, particularly for complex integrals like shape Hessians.

Methodology
The authors develop a framework combining unfitted finite element methods with automatic shape differentiation, implemented within the Julia package GridapTopOpt and the broader Gridap ecosystem.

1. Extension of Analytic Shape Calculus:
   The paper extends existing analytic results (Berggren [16]) to cases where the domain boundary $\partial\Omega$ intersects the background domain boundary $\partial D$. The authors derive directional shape derivatives for general volume and surface integrals, as well as flux integrals, under these intersection conditions. They prove that the semiderivatives remain consistent even when the domain boundary touches the background boundary, provided specific assumptions regarding mesh intersections are met.

2. Automatic Shape Differentiation (ASD):
   The core methodological contribution is the application of forward-mode automatic differentiation to unfitted discretisations.

   • Dual Numbers: The method utilizes dual numbers ($a + b\epsilon$) to propagate derivatives through the quadrature rules of cut cells.
   • Implementation: The level-set function values at nodes are replaced with dual numbers. The geometry of cut cells (defined by the intersection of the level-set zero-isosurface and the background mesh) and the associated Jacobians are computed using these dual numbers.
   • Verification: The authors demonstrate that ASD recovers the analytic directional derivatives derived in Section 2 to machine precision, regardless of mesh size. This eliminates the need for manual derivation of complex mesh-related quantities. Furthermore, the technique allows for the computation of shape Hessians and derivatives of functionals for which rigorous analytic counterparts (within the Berggren framework) are not yet established.

3. Level-Set Evolution and Reinitialisation:
   The level-set function is evolved using a transport equation with a velocity field derived from a Hilbertian extension-regularisation of the sensitivity. The authors propose a modified interior penalty term that includes the velocity magnitude $|n_F \cdot \beta|$ to prevent non-designable regions from being over-stabilised. Reinitialisation to a signed distance function is performed using a steady-state approach with artificial viscosity or continuous interior penalty terms.

4. Isolated Volume Detection:
   A novel graph-based approach is proposed to detect isolated volumes. Unlike previous methods relying on thermal convection problems (which fail to distinguish volumes separated by a single cell layer), this method constructs a conforming connectivity graph of the cut mesh (using polytopal cuts to ensure conformity) and applies a Breadth-First-Search (BFS) colouring algorithm. This identifies disconnected regions, allowing for the application of penalty terms to constrain the average displacement or pressure in these volumes, thereby regularising the system matrices.

5. Computational Framework:
   The implementation supports both serial and memory-distributed computing using GridapDistributed. It leverages GridapEmbedded for computation on subtriangulations and GridapTopOpt for the optimisation loop, which uses C´ea's method (adjoint method) combined with ASD for PDE-constrained optimisation.

Key Results
The framework is validated through two three-dimensional examples:

1. Minimum Compliance of a Linear Elastic Wheel: A topology optimisation problem involving a 3D elastic wheel with a volume constraint. The problem was solved on an unstructured background mesh with approximately 1.32 million elements using 48 processors in a distributed memory framework. The results demonstrated successful convergence and the ability to handle complex geometries without analytic derivation of shape derivatives.
2. Fluid-Structure Interaction (FSI): A topology optimisation of a linear elastic structure interacting with Stokes flow. This problem required handling one-way coupling between the fluid and solid phases. The authors introduced a penalty term to the objective function to prevent the formation of isolated solid volumes, which would otherwise artificially improve the compliance objective by blocking fluid flow. The problem was solved using staggered (sequential) solution of the Stokes and elasticity equations.

In both cases, the automatic differentiation techniques successfully computed the required gradients, and the graph-based isolated volume detection effectively managed disconnected regions.

Significance and Claims
The paper claims to provide a general and robust implementation for level-set topology optimisation on unstructured meshes that avoids the burdens of analytic shape derivative derivation. Specifically:

• Precision: The authors assert that their automatic differentiation approach recovers analytic derivatives to machine precision, validating the theoretical link between ASD and the discretise-then-differentiate approach of Berggren [16].
• Generality: By extending the analytic framework to include intersections between the domain and background boundaries, the method addresses a gap in existing literature relevant to general topology optimisation problems.
• Flexibility: The use of ASD allows for the optimisation of problems involving complex physical phenomena (e.g., FSI, non-standard function spaces like $H(\text{div})$ and $H(\text{curl})$) without requiring new mathematical derivations for every new functional.
• Scalability: The implementation in GridapTopOpt supports distributed computing, enabling the solution of large-scale problems on unstructured triangulations.

The authors position this work as a step towards making unfitted topology optimisation more accessible and applicable to a wider range of engineering problems by automating the most mathematically demanding aspects of the process. The software is released as open source.