Unfitted schemes

Warning

This page is still under construction! Please refer to the tutorial for a discussion of functionality.

UnfittedEvolution

Level-set evolution

In the following, we use the approach for evolution proposed by Burman et al. (2018) to update the level-set function. Namely, we evolve the level-set function by solving a transport equation

\[\begin{equation} \begin{cases} \displaystyle\frac{\partial\phi(t,\boldsymbol{x})}{\partial t}+\boldsymbol{\beta}\cdot\boldsymbol{\nabla}\phi(t,\boldsymbol{x})=0,\\ \phi(0,\boldsymbol{x})=\phi_0(\boldsymbol{x}),\\ \boldsymbol{x}\in D,~t\in(0,T), \end{cases} \end{equation}\]

where $\boldsymbol{\beta}$ is a velocity field. In this work $\boldsymbol{\beta}$ is computed as $\boldsymbol{\beta}=\boldsymbol{n}g_\Omega$, where $g_\Omega$ is the regularised sensitivity resulting from the Hilbertian extension-regularisation approach. The transport equation above is solved using an interior penalty approach and Crank-Nicolson for the discretisation in time (see Burman et al. (2018) and references there in for further discussion). The weak formulation of this problem is: for all $t\in(0,T)$, find $\phi\in W$ such that

\[\int_D\left[v\frac{\partial\phi}{\partial t}+v\boldsymbol{\beta}\cdot\boldsymbol{\nabla}\phi\right]~\mathrm{d}\boldsymbol{x}+\sum_{F\in\mathscr{S}_h}\int_Fc_eh_F^2\lvert\boldsymbol{n}_F\cdot\boldsymbol{\beta}\rvert\llbracket \boldsymbol{n}_F\cdot\boldsymbol{\nabla}\phi\rrbracket\llbracket\boldsymbol{n}_F\cdot\boldsymbol{\nabla}v\rrbracket~\mathrm{d}s=0,~\forall v\in W^h,\]

where $\mathscr{S}_h$ is the set of interior mesh facets, $h_F$ is the average element diameters of the elements sharing a facet, $\boldsymbol{n}_F$ is the normal to the facet $F$, $\llbracket v\rrbracket = v^+-v^-$ is the jump in a function $v$ over the facet $F$, and $c_e$ is a stabilisation coefficient. Further discussion regarding solving this can be found at Wegert et al. (2025).

GridapTopOpt.CutFEMEvolveType
mutable struct CutFEMEvolve{V,M} <: Evolver

CutFEM method for level-set evolution based on method developed by

  • Burman et al. (2018). DOI: 10.1016/j.cma.2017.09.005.
  • Burman et al. (2017). DOI: 10.1016/j.cma.2016.12.021.
  • Burman and Fernández (2009). DOI: 10.1016/j.cma.2009.02.011

Parameters

  • ode_solver::ODESolver: ODE solver
  • Ωs::B: EmbeddedCollection holding updatable triangulation and measures from GridapEmbedded
  • dΩ_bg::C: Measure for integration
  • space::B: Level-set FE space
  • assembler::Assembler: FE assembler
  • params::D: Tuple of stabilisation parameter γg, mesh sizes h, and max steps max_steps, and background mesh skeleton parameters
  • cache: Cache for evolver, initially nothing.
Warning

Caching for the CutFEMEvolve method is currently disabled. This will be re-enabled in the future."

source
GridapTopOpt.CutFEMEvolveMethod
CutFEMEvolve(V_φ::B,Ωs::EmbeddedCollection,dΩ_bg::A,h;
  max_steps=10,
  γg = 0.1,
  ode_ls = LUSolver(),
  ode_nl = ode_ls,
  ode_solver = MutableRungeKutta(ode_nl, ode_ls, 0.1, :DIRK_CrankNicolson_2_2),
  assembler=SparseMatrixAssembler(V_φ,V_φ)
) where {A,B}

Create an instance of CutFEMEvolve with the space for the level-set V_φ, the EmbeddedCollection Ωs for the triangulation and measures, the measure dΩ_bg for the background mesh, and the mesh size h. The mesh size h can either be a scalar or a CellField object.

The optional arguments are:

  • max_steps: Maximum number of steps for the ODE solver.
  • γg: Stabilisation parameter for the continuous interior penalty term.
  • ode_ls: Linear solver for the ODE solver.
  • ode_nl: Non-linear solver for the ODE solver.
  • ode_solver: ODE solver, default is MutableRungeKutta(ode_nl, ode_ls, 0.1, :DIRK_CrankNicolson_2_2).
  • assembler: Assembler for the finite element space, default is SparseMatrixAssembler(V_φ,V_φ).

Note

  • The stepsize dt = 0.1 in MutableRungeKutta is a place-holder and is updated using the γ passed to solve!.
source

Level-set reinitialisation

We solve the reinitialisation equation

\[\left\lbrace\begin{aligned} \frac{\partial\phi(t,\boldsymbol{x})}{\partial t} + \operatorname{sign}(\phi_0(\boldsymbol{x}))\left(\lvert\boldsymbol{\nabla}\phi(t,\boldsymbol{x})\rvert-1\right) = 0,\\ \phi(0,\boldsymbol{x})=\phi_0(\boldsymbol{x}),\\ \boldsymbol{x}\in D, t>0. \end{aligned}\right.\]

To steady state using an approach based on the one proposed by Mallon et al. (2025). The weak formulation for this problem is given as: find $\phi\in W$ such that

\[ \int_Dv\boldsymbol{w}\cdot\boldsymbol{\nabla}\phi-v\operatorname{sign}(\phi_0)~\mathrm{d}\boldsymbol{x}+\int_\Gamma \frac{\gamma_d}{h}\phi v~\mathrm{d}s+j(\phi,v)=0,~\forall v\in W^h,\]

where $\boldsymbol{w}=\operatorname{sign}(\phi_0)\frac{\boldsymbol{\nabla}\phi}{\lVert\boldsymbol{\nabla}\phi\rVert}$. Further discussion regarding solving this can be found at Wegert et al. (2025).

GridapTopOpt.StabilisedReinitType
mutable struct StabilisedReinit{V,M} <: Reinitialiser

Stabilised FE method for level-set reinitialisation. Artificial viscosity approach (ArtificialViscosity) based on Mallon et al. (2023). DOI: 10.48550/arXiv.2303.13672. Interior jump penalty approach (InteriorPenalty) adapted from that work and replaces the artifical viscosity term with an interior jump penalty term.

Parameters

  • nls::NonlinearSolver: Nonlinear solver for solving the reinitialisation equation
  • stabilisation_method::A: A StabilisationMethod method for stabilising the problem
  • Ωs::B: EmbeddedCollection holding updatable triangulation and measures from GridapEmbedded
  • dΩ_bg::D: Background mesh measure
  • space::C: FE space for level-set function
  • assembler::Assembler: Assembler for LS FE space
  • params::E: Tuple of Nitsche parameter γd and mesh size h
  • cache: Cache for reinitialiser, initially nothing.

Note

  • We expect the EmbeddedCollection Ωs to contain :dΓ. If this is not available we add it to the recipe list in Ωs and a warning will appear.
source
GridapTopOpt.StabilisedReinitMethod
StabilisedReinit(V_φ::C,Ωs::EmbeddedCollection,dΩ_bg::B,h;
  γd = 20.0,
  nls = NewtonSolver(LUSolver();maxiter=20,rtol=1.e-14,verbose=true),
  assembler=SparseMatrixAssembler(V_φ,V_φ),
  stabilisation_method::A = InteriorPenalty(V_φ)
) where {A,B,C}

Create an instance of StabilisedReinit with the space for the level-set V_φ, the EmbeddedCollection Ωs for the triangulation and measures, the measure dΩ_bg for the background mesh, and the mesh size h. The mesh size h can either be a scalar or a CellField object.

The optional arguments are:

  • γd: Interface penalty parameter for the reinitialisation equation.
  • nls: Nonlinear solver for solving the reinitialisation equation.
  • assembler: Assembler for the finite element space.
  • stabilisation_method: A StabilisationMethod method for stabilising the problem.
source

Custom UnfittedEvolution

GridapTopOpt.EvolverType
abstract type Evolver

Your own unfitted level-set evolution method can be created by implementing concrete functionality for solve!.

source
GridapTopOpt.ReinitialiserType
abstract type Reinitialiser

Your own unfitted level-set reinitialisation method can be created by implementing concrete functionality for solve!.

source
Gridap.Algebra.solve!Method
solve!(::Reinitialiser,φ,args...)

Reinitialise the level set function φ according to an Reinitialiser method.

source