A reduced order model for a stable embedded boundary parametrized Cahn-Hilliard phase-field system based on cut finite elements
AA REDUCED ORDER MODEL FOR A STABLE EMBEDDED BOUNDARYPARAMETRIZED CAHN-HILLIARD PHASE-FIELD SYSTEM BASED ON CUTFINITE ELEMENTS
EFTHYMIOS N. KARATZAS
AND GIANLUIGI ROZZA Abstract.
In the present work, we investigate for the first time with a cut finite element method,a parameterized fourth order nonlinear geometrical PDE, namely the Cahn-Hilliard system. Wemanage to tackle the instability issues of such methods whenever strong nonlinearities appear and toutilize their flexibility of the fixed background geometry –and mesh– characteristic, through which,one can avoid e.g. in parametrized geometries the remeshing on the full order level, as well as,transformations to reference geometries on the reduced level. As a final goal, we manage to findan efficient global, concerning the geometrical manifold, and independent of geometrical changes,reduced order basis. The POD-Galerkin approach exhibits its strength even with pseudo-randomdiscontinuous initial data verified by numerical experiments. Introduction and Motivation
The Cahn–Hilliard model (CH), named after John Cahn and John Hilliard who suggested the sys-tem in 1958, describes a prototype of the process of phase separation, by which the two componentsof a binary fluid or material impulsively separate and form domains, pure in each component. Thesephase-field systems, are very interesting in the scientific community due to their great conservationproperties. They simulate many industrial systems, as the two-phase fluid flows for capturing theinterface location between two immiscible fluids, [3, 4, 47], spinodal decomposition in binary alloys–a process in which a mixture of two fluids or materials evolves– [38, 49], and the phase diagram formicrophase multiscale separation for diblock copolymer-linear chain molecule consisting of two sub-chains joined covalently to each other, [21, 48]. Additionally, we mention the image inpainting, i.e., thefilling in of damaged or missing regions of an image with the use of information from surrounding ar-eas, [10, 18, 65], micro-structure with elastic inhomogeneity which determines the transformation pathand the corresponding microstructure evolution, [74] and references therein, tumor growth simulationin order to provide optimal strategies for treatments, [1, 77], and topology optimization phase-fieldapproach to the problem of minimizing the mean compliance of a multi-material structure, see e.g.[83, 11, 64].The investigation of the behavior of such systems, started from the pioneer work [16] and theearly works of [62, 31, 29]. Thereafter, a lower solution space regularity and a second order splittingmethod have been investigated in [30], the solution existence and error analysis in [28, 51], for higherorder finite element we refer to [34] and for Discontinuous Galerkin in space approach to [76, 60].Optimal control for non-convective or optimal control for the convective case have been studied in[66, 81, 35, 45, 25, 44], stochastic partial differential equations and stochastic analysis in the very newwork of [33], Navier-Stokes/Cahn-Hilliard systems in [46, 43], and Cahn-Hilliard/Allen-Cahn systemsin [6, 5, 80, 2].Throughout this work, an efficient methodology for solving nonlinear systems governed by Cahn–Hilliard equations is studied within cut finite elements and reduced-order modeling. A stable full orderscheme for this geometrically parameterized nonlinear fourth-order diffusion system is introduced. Inthis embedded geometry framework, we propose a model order reduction technique using the advan-tages of a shape regular background mesh, recently investigated in [52, 55, 54]. The combination ofunfitted mesh finite element methods and reduced order modeling allows us to obtain a fast evaluation, SISSA, International School for Advanced Studies, Mathematics Area, mathLab Trieste, Italy. Department of Mathematics, School of Applied Mathematical and Physical Sciences, National Techni-cal University of Athens, Zografou, Greece.
E-mail addresses : [email protected] & [email protected], [email protected] .2010 Mathematics Subject Classification.
Key words and phrases.
Cut Finite Element Method, Cahn Hilliard, Reduced Order Model, POD, Stabilization. a r X i v : . [ m a t h . NA ] S e p ROM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS by considering the geometrical parametrized system, while we avoid remeshing, as well as the referencedomain formulation, often used in boundary fitted finite element formulations.This contribution is organized as follows: In Section 2 we define the continuous strong formulation ofthe mathematical problem. The semidiscrete cut elements Nitsche weak formulation and the implicitexplicit Euler method (IMEX) fully discrete problem under consideration is introduced, as well as,the incremental scheme used to solve the full order problem during the offline stage. In Section 3, wedemonstrate the reduced-order model formulation, based on the Proper Orthogonal Decomposition,and its main aspects. In section 4 the IMEX method and the high fidelity solver efficiency is validated.Subsequently, the proposed ROM technique is tested on a geometrical parametrized problem of a two-phase field problem starting from pseudo-random initial data around an embedded circular domain.Convergence results, errors and reduced execution times are introduced and analyzed. Finally inSection 5, conclusions and perspectives for future improvements and developments are demonstrated.To our best knowledge, the results of this work are original and applicable in many cases of nonlineartime-depended partial differential equation problems in terms of the prescribed methodology and inthe spirit of the geometrical parametrization.2.
The Model problem and the full-order approximation
Strong formulation of the Cahn-Hilliard problem.
We consider the model problem describ-ing a phase flow time evolution. An unknown function u indicates the perturbation of the concentrationof one of the phases of fluid components constituting a liquid mixture that contains a binary fluid. Letus consider an open domain Ω in R d , with d = 2 , ∂Ω . Let a k − dimensional parameter space P with a parameter vector µ ∈ P ⊂ R k . We statebelow, in a time interval [0 , T ], the strong form of the evolutionary Cahn-Hilliard phase flow systemof equations with open boundary conditions on ∂Ω , geometrically parametrized by µ . We denote by Ω ( µ ) and ∂Ω ( µ ) the parametrized domain and boundary, respectively.As first suggested by [26], and thereafter extended in [16], if we assume that the mobility is equalto 1 and ε is a measure of the size of the interface of two fluids, the mass flux is given by(1) J ( µ ) = −∇ (cid:18) ε F ( u ( µ )) − ε ∆u ( µ ) (cid:19) , where F denotes the chemical potential difference between the two species. Due to [16], the Ginzburg–Landau energy becomes E ( u ( µ )) = (cid:90) Ω (cid:18) F ( u ( µ )) + ε |∇ u ( µ ) | (cid:19) d x . (2)An equilibrium state of the considered mixture minimizes the above Ginzburg–Landau energy, subjectto the mass conservation ∂∂t (cid:90) Ω ( µ ) u ( µ ) d x = 0 . (3)Hence, the parametrized Cahn–Hilliard system can be described as: ∂u ( µ ) ∂t = − ε ∆ u ( µ ) + 1 ε ∆F (cid:48) ( u ( µ )) , in Ω ( µ ) × [0 , T ] , (4) ∂ n u ( µ ) = ∂ n ( − ε ∆u ( µ ) + 1 (cid:15) F (cid:48) ( u ( µ ))) = g N ( µ ) , on ∂Ω ( µ ) × [0 , T ] , (5) u ( · ,
0) = u ( · ) , in Ω ( µ ) , (6)where n is the unit outer normal vector of ∂Ω , and F is a double-well usually taken as a polynomialfunction of u of fourth power:(7) F ( u ( µ )) = γ u ( µ )4 + γ u ( µ )3 + γ u ( µ )2 with γ > . For more details, the interested reader is referred to [31] as well as for Dirichlet boundary conditionsin [61]. We remark that the equation (4) represents the equation of conservation of mass, with massflux as described in (1).2.2.
Full-order parametrized Nitsche cut elements weak variational formulation.
OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 3
A proper continuous weak formulation.
The Cahn-Hilliard equation as it is expressed in Equa-tions (4)–(6) is a fourth-order diffusion equation, involving first-order time derivatives, second andfourth-order spatial derivatives. Casting it in a weak form results in second-order spatial derivativesavoiding the fourth-order ones. Nevertheless, the problem still can not be solved using the standardLagrange finite element basis. The setback also is that if someone employs the Nitsche weak bound-ary enforcement, several integrals in the cut geometry should be calculated including various orderderivatives and normal derivatives, [76, 32, 82, 41], which in our case of unfitted mesh are avoided dueto time expensive integration. Moreover, a special stabilization term for the nonlinearity would beneeded, [24, 13].To overcome all these difficulties arising from this fourth-order nonlinear system, firstly we employa splitting method deriving a strong formulation system that requires H space regularity. Afterward,we multiply with a test function and we integrate by parts over Ω in order to drive the systemto a weak form that requires only H space regularity. In particular the pair solution ( u ( µ ) , w ( µ ))solves the following problem: find u ( µ ) ∈ L [0 , T ; H ( Ω ( µ ))] ∩ H [0 , T ; ( H ( Ω ( µ ))) (cid:48) ] and w ( µ ) ∈ L [0 , T ; H ( Ω ( µ ))] for all weight functions v ( µ ) ∈ H ( Ω ( µ )) such that( ∂u ( µ ) ∂t , v ( µ )) + ( ∇ w ( µ ) , ∇ v ( µ )) = ( g N ( µ ) , v ( µ )) ∂Ω ( µ ) , (8) − ( w ( µ ) , v ( µ )) + ε ( ∇ u ( µ ) , ∇ v ( µ )) + ε − ( F (cid:48) ( u ( µ )) , v ( µ )) = ( g N ( µ ) , v ( µ )) ∂Ω ( µ ) , (9) u ( · , µ ) = u ( · ; µ ) , (10)where w ( µ ) is depended on the geometry parameter µ and usually it is identified as a chemical potential, ∀ v ∈ H ( Ω ( µ )) and u ∈ L ( Ω ( µ )). In the following, we are focused on the concentration component u ( µ ), and we handle w ( µ ) as auxiliary function. Similarly, Dirichlet boundaries will be examinedfollowing [61]. Remark . We point out that we are taking into account nonsmooth initial data L ( Ω ( µ )) while theaforementioned regularity H [0 , T ; ( H ( Ω ( µ ))) (cid:48) ] of u ( µ ) is needed for the ∂u ( µ ) /∂t term , see [22, 23].2.2.2. Discretization, unfitted mesh and stability issues.
The system (8)–(9) in a discrete unfittedmesh formulation needs extra attention. Even if we use classical finite element methods and we applythe boundary conditions in a strong way, stability issues appear, and to achieve a stable solution aquite small time-stepping size seems necessary. We highlight that the evolution of the physics of theproblem during the first time steps is very fast, although, as time passes it slows down and finally itequilibrates, see also e.g. [78, 12], and references therein. An adaptive time-stepping approach wouldappear beneficial, although we will investigate it in a future work, as well as, the way it affects theaccuracy and efficiency of the reduced model detailed in Sections 3 and 4. Considering the Nitscheterms one may apply the simple formulation needed for linear systems. Under these considerations,we used cut finite element methods to solve the system applying a jump stabilization procedure in theboundary elements interface area. In the next paragraph, we derive the semi-discrete formulation ofthe system, while in paragraph 2.2.4 the fully discrete system is exploited.2.2.3.
Semidiscrete variational formulation.
We denote by T the background domain, and by T h itscorresponding mesh, see e.g. Figure 1. Considering the standard notation ( · , · ), (cid:104)· , ·(cid:105) ∂Ω for the L ( Ω ),and L ( ∂Ω ) inner products onto the truth geometry Ω , and ∂Ω , respectively. For every element K ∈ T h , we associate a parameter h K , denoting the diameter of the set K . The size of the mesh isdenoted by h = max K ∈T h h K .The continuous boundary value problem is next formulated on a domain ˜ Ω ( µ ) that contains Ω ( µ ) ⊂ ˜ Ω ( µ ), while its mesh ˜ Ω T ( µ ) := T h ∩ ˜ Ω ( µ ) is not fitted to the domain boundary and ˜ Ω T ( µ ) ⊂ T forall µ ∈ K . Let also Ω T ( µ ) := T h ∩ Ω ( µ ) and G h ( µ ) := { K ∈ T h ( µ ) : K ∩ ∂Ω ( µ ) (cid:54) = ∅} be the set ofelements that are intersected by the interface. We remark that G h ( µ ) and ˜ Ω T ( µ ) depend on µ through˜ Ω ( µ ) (or its boundary), while the background domain T and its mesh T h do not depend on µ .Furthermore, the set of element faces F G ( µ ) associated with G h ( µ ), is defined as follows: for eachface F ∈ F G ( µ ), there exist two simplices K (cid:54) = K (cid:48) , such that F = K ∩ K (cid:48) , and at least one of the twois a member of G h ( µ ). Note that the boundary faces of ˜ Ω T ( µ ) are excluded from F G ( µ ). On a face F ∈ F G ( µ ), F = K ∩ K (cid:48) , the jump of v is defined by [[ v ]] = v | K − v | K (cid:48) and the jump of the gradient[[ n F · ∇ v ]] = n F · ∇ v | K − n F · ∇ v | K (cid:48) , where n F denotes the outward pointing unit normal vector to F . namely, ∂u ( µ ) /∂t ∈ L [0 , T ; ( H ( Ω ( µ ))) (cid:48) ]. ROM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS (i) (ii)(iii) (iv)
Figure 1. (i) The geometry of an embedded disk, (ii) the cut finite element methodgeometry, (iii) zoom on the background mesh together with the surrogate cut dis-cretized geometry, and (iv) extended mesh and elements intersected by the true bound-ary.Let the space of continuous piecewise-linear functions V h ( Ω ( µ ))(11) V h ( Ω ( µ )) = (cid:8) υ ∈ C ( Ω T ( µ )) : υ | K ∈ P ( K ) , ∀ K ∈ T h ( µ ) (cid:9) . For the sake of simplicity, in the next set of equations, we will omit the parameter dependency notationwith respect to µ and the CutFEM discretization is as follows. We seek u h , w h ∈ V h ( ˜ Ω ( µ )) such that,for all the test functions v h , q h ∈ V h ( ˜ Ω ( µ )),( u ht , v h ) + ( ∇ w h , ∇ v h ) + ε ( ∇ u h , ∇ q h ) − ( w h , q h ) + 1 ε ( γ u h + γ u h + γ u h , q h )+ (cid:104) α N h n Ω · ∇ u h , n Γ · ∇ v h (cid:105) ∂Ω + (cid:88) F ∈F G (cid:0) α h − [[ u h ]] , [[ v h ]] (cid:1) F = (cid:104) g N , v h + q h + α N h n Γ · ∇ v h (cid:11) ∂Ω , (12)where α N , and α are positive penalty parameters related to Nitsche weak imposition of boundaryconditions and the boundary interface stabilization term respectively, see for instance [14]. Remark . For the sake of completeness, it would be convenient to emphasize that both the afore-mentioned kind of “jumps”, can be used to apply different type of ghost penalty stabilizations, namely (cid:80) F ∈F G ( α h [[ n F · ∇ u h ]] , [[ n F · ∇ v h ]]) F and (cid:80) F ∈F G (cid:0) α h − [[ u h ]] , [[ v h ]] (cid:1) F . Although, we prefer to employthe “projection based”-jump and not the “derivative”-jump one. All experiments consider first orderpolynomials. Finally, we highlight that the preferred ghost penalty is less computational expensivesince it does not involve any kind of derivatives. Concerning literature, most theoretical estimates arerelated to the “derivative-jump”, nevertheless, they can easily be extended to the “projection based”-jump ghost penalty since the jump can be bounded from above with the derivative jump, for moredetails we refer to [13, 15, 73, 59].2.2.4. IMEX type time discretization.
The goal here is to find concentration u h = u h ( x, t ) satisfyingthe equations (12) for all time instances t in a time interval [0 , T ] ⊂ R and space positions x ∈ Ω ⊂ R d .We fully discretize the system by an IMEX approach, see e.g. [67] and references therein, as wellas [57] for similar types of nonlinearities. Approximations will be constructed on a time partition0 = t < t < . . . < t N = T on which, each interval I n := ( t n , t n +1 ] is of length τ n = t n +1 − t n , n = 1 , ..., N −
1, starting from initial condition u . Therefore, we apply a splitting method onto timeintegration level –often called operator splitting method– which means that the differential operatoris rewritten as the sum of two complementary operators. The latter extensively attracted attentionIMEX method technique treats on the nonlinear term explicitly, allowing to act as a forcing term inthe w -equation (8). In this way, we avoid stability issues caused by the nonlinearity. All the other OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 5 terms have been handled implicitly for better accuracy. Hence, the IMEX method results in the set offully-discrete Cahn-Hilliard equations:Find ( u n +1 h , w n +1 h ) ∈ V h × V h , s.t. for all ( v nh , q nh ) ∈ V h × V h A ( u n +1 h ) + τ n · L ( u n +1 h , w n +1 h ) = A ( u nh ) − τ n · N ( u nh ) + B ( g N ) · τ n , where A ( u n +1 h ) = ( u n +1 h , v n +1 h ) = A u n +1 h ,L ( u n +1 h , w n +1 h ) = ( ∇ w n +1 h , ∇ v n +1 h ) + ε ( ∇ u n +1 h , ∇ q n +1 h ) − ( w n +1 h , q n +1 h )+ (cid:104) α N h n Ω · ∇ u n +1 h , n Γ · ∇ v n +1 h (cid:105) ∂Ω + (cid:88) F ∈F G (cid:0) α h − [[ u h ]] , [[ v h ]] (cid:1) F = L u n +1 h + L w n +1 h + C F G u n +1 h ,B ( g N ) = (cid:104) g N , v nh + q nh + α N h n Γ · ∇ v nh (cid:11) ∂Ω = B v + B q ,N ( u nh ) = 1 ε ( γ u nh + γ u nh + γ u nh , v nh ) = N ( u nh ) u nh . Next, we derive the related to that system matrix form. We define the parameter-depended Cahn-Hilliard operator G ( U nh ( µ )) := G (cid:18)(cid:20) u nh ( µ ) w nh ( µ ) (cid:21)(cid:19) = (cid:20) L + N ( u nh ( µ )) L C F G (cid:21) (cid:20) u nh ( µ ) w nh ( µ ) (cid:21) , and the right hand side consists of the forcing boundary data related to stabilization and Nitsche weakenforcement boundary terms F N ( µ ) := (cid:20) B v ( µ ) B q ( µ ) (cid:21) . These definitions result in the following residual R ( U nh ( µ )) = G ( U nh ( µ )) − F N ( µ ) , which yield the following algebraic system of equations for the increment δU n +1 h ( µ ) = (cid:20) δu nh ( µ ) δw nh ( µ ) (cid:21) = U n +1 h ( µ ) − U nh ( µ ): (cid:20) A + τ n L τ n L τ n C F G (cid:21) δU n +1 h ( µ ) = − τ n R ( U nh ( µ )) , (13)noting that in the above formulation the nonlinear term is treated only explicitly and the remainingpart implicitly. Remark . IMEX method approach is beneficial since we avoid the extra computation of iterativeapproaches, e.g. Newton-Raphson method, and simultaneously we avoid the small-time stepping ofan explicit time integration method, needed for a stable solution in order to minimize the dispersionerror associated with these schemes. In the above system of equations, it is important to underlinethat the discretized differential operators A , L , L and C F G are parameter-dependent, a feature ofgreat importance as we will see in Section 3 and the ROM basis construction. Also, a pre-assemblingtechnique for the involved matrices will be employed in Section 4, by minimizing the time-consumingintegration of several involved inner products.3. Reduced order model with a POD-Galerkin method
In this paragraph, a POD-Galerkin approach is briefly recalled as in [42, 68]. We emulate thehigh fidelity model system with a reproduced one, which allows predictive errors, within the aim of areduced computational cost and solution time in a way adjusted to embedded-immersed boundary finiteelement methods. This reduced-order system has been approved advantageous when geometricallydeformed systems appear and in comparison with traditional finite element methods and/or reduced-order modeling, see for instance [7, 52, 54, 55]. In particular, we employ a projection-based reducedorder model which consists of the projection of the governing equations onto the reduced basis spaceconstructed on a fixed background mesh.Following the literature, one could see reduced basis (RB) methods applied to linear elliptic equationsin [70], to linear parabolic equations in [37] and to non-linear problems in [75, 36]. Although the numberof works on reduced-order models with classical FEM are now significant big including Cahn-Hilliardsystems, see e.g. [35, 42, 68, 70, 37, 75, 36] and references therein, to the best of the authors’ knowledge,only very few research works [7, 52, 54, 55, 56] can be found concerning embedded boundary methods onlinear systems and ROMs and much fewer for nonlinear, [56, 53]. So, for the first time we investigate
ROM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS and focus onto how we can achieve for a fourth order evolutionary Cahn–Hilliard PDE system, in aFull Order Method (FOM) and in a Reduced Order Model (ROM) framework, a stable solution withinan embedded finite element method namely in a cut finite elements setting. The new techniques of[52, 54, 55], based on a proper orthogonal decomposition (POD) strategy will be employed. The keyfeature of our approach is that we avoid the remeshing effort or/and the need of a map of all thedeformed geometries to reference geometries often used in fitted mesh finite element methods, see e.g.[42, 72, 69, 8, 71, 70, 9].Regarding the reduced-order modeling, we investigate how ROMs can be applied to time dependedcut finite element methods and generally, to embedded boundary methods simulations considering timedepended nonlinear systems and Cahn–Hilliard systems. The main interest is to generate ROMs onparametrized geometries. The cut elements unfitted mesh finite element method with levelset geometrydescription is used to apply parametrization and the reduced order techniques (offline-online). Animportant aspect is also to test the efficiency of a geometrically parametrized reduced-order nonlinearmodel without the usage of the transformation to reference domains, which is an important advantageof embedded methods relying on fixed background meshes.Before going into thorough, we specify some basics for reduced basis modeling. We always start bythe generation of a set of full order solutions of the parametrized problem under a parameter valuesrandom choice. The final objective of RB methods is to emulate any member of this solution set witha low number of basis functions and this is based on a two-stage procedure, the offline and the onlinestage, [63, 72, 39].
Offline stage.
In order to derive reduced-order solutions emulating the full order system, a low dimen-sional reduced basis is constructed based on a specific number of full order solves. This reduced basiswill be able to approximate any member of the solution set to a predictive error accuracy. Hence, it ispossible to project the FOM differential operators, describing the governing equations, onto the reducedbasis space, applying a Galerkin projection technique and to create a reduced system of equations.The offline stage is computationally very expensive, nevertheless, it is executed only once.
Online stage.
Thereafter, during this stage, one can, even with very few computational resources,calculate a reproduced reduced system of equations involving any new value of the input parameters.For more details and applications, we refer to [9, 19] and the references therein.3.1.
POD.
The full-order model, as illustrated in paragraph 2.2.4, is solved for each µ k ∈ K = { µ , . . . , µ N k } ⊂ P where K is a finite-dimensional training set of parameters chosen inside the pa-rameter space P . The considered problem can be simultaneously parameter and time-dependent. Inorder to collect snapshots for the generation of the reduced basis spaces, one needs to consider boththe time and parameter dependency. For this reason, discrete-time instants t k ∈ { t , . . . , t N t } ⊂ [0 , T ]belong in a finite-dimensional training set, which is a subset of the simulation time interval and areconsidered as parameters. The total number of the full order snapshots is then equal to N s = N k · N t .The snapshots matrices S u and S w , are then given by N s full-order snapshots: S u = [ u ( µ , t ) , . . . , u ( µ N r , t N t )] ∈ R N hu × N s , (14) S w = [ w ( µ , t ) , . . . , w ( µ N r , t N t )] ∈ R N hw × N s , (15)where N hu and N hw are the number of degrees of freedom for the discrete full-order solution for theconcentration u and the auxiliary variable w , respectively. The reduced-order problem can be efficientlysolved for both sets of parameters and time instants. In order to generate the reduced basis spaces,for the projection of the governing equations, one can find in literature several techniques such asthe Proper Orthogonal Decomposition (POD), the Proper Generalized Decomposition (PGD) and theReduced Basis (RB) with a greedy sampling strategy. For more details about the different strategiesthe reader may see [42, 70, 19, 50, 63, 20, 27]. In this work, the POD strategy is applied onto the fullsnapshots matrices. The aforementioned procedure includes both time and parameter dependency. Inthe case of parametric and time-dependent problems also other approaches are available such as thePOD-Greedy approach [39] or the nested POD approach, where the POD is applied first in the timedomain and then on the parameter space. Given a general scalar u ( t ) in R d , with a certain number ofrealizations u , . . . , u N s , the POD problem consists in finding, for each value of the dimension of PODspace N P OD = 1 , . . . , N s , the scalar coefficients a , . . . , a N s , . . . , a N s , . . . , a N s N s and functions ϕ , . . . , ϕ N s OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 7 minimizing the quantity: E N POD = N s (cid:88) i =1 || u i − N POD (cid:88) k =1 a ki ϕ k || ∀ N P OD = 1 , . . . , N, (16) with (cid:104) ϕ i , ϕ j (cid:105) L ( Ω ) = δ ij ∀ i, j = 1 , . . . , N s . (17)The detailed results for the auxiliary variable w are omitted for shortness of space. It can be shown,[58], that the minimization problem of Equation (16) is equivalent of solving the following eigenvalueproblem: C u Q u = Q u λ u , (18) C uij = (cid:104) u i , u j (cid:105) L ( Ω ) for i, j = 1 , . . . , N s , (19)where C u is the correlation matrix obtained starting from the snapshots S u , Q u is a square matrix ofeigenvectors and λ u is a vector of eigenvalues. The basis functions can then be obtained with:(20) ϕ i = 1 N s λ ui N s (cid:88) j =1 u j Q uij . The POD spaces are constructed using the aforementioned methodology resulting in the spaces: B u = span { ϕ , . . . , ϕ N ru } ∈ R N hu × N ru , (21)and similarly for the auxiliary variable: B w = span { χ , . . . , χ N rw } ∈ R N hp × N rw , (22)where N ru , N rw < N s are chosen according to the eigenvalue decay of the vectors of eigenvalues λ u and λ w .Once the POD functional spaces are set, the reduced quantities fields can be approximated with:(23) u r ≈ N ru (cid:88) i =1 a i ( t, µ ) ϕ i ( x ) , w r ≈ N rw (cid:88) i =1 b i ( t, µ ) χ i ( x ) , where the coefficients a i and b i depend only on the time and parameter spaces and the basis functions ϕ i and χ i depend only on the physical space and not on the parametrized geometry.By denoting B = (cid:20) B u B w (cid:21) and B T = (cid:20) B Tu B Tw (cid:21) , the unknown vectors of coefficients V = (cid:20) ab (cid:21) then can be obtained through a Galerkin projection of the full-order system of equations onto thePOD reduced basis spaces with the resolution of a consequent reduced iterative algebraic system ofequations for the increment δV n +1 h ( µ ) = V n +1 h ( µ ) − V nh ( µ ),(24) B T (cid:20) A + τ n L τ n L τ n C F G (cid:21) B δV n +1 h ( µ ) = − τ n B T R ( B V nh ( µ )) , which leads to the following algebraic reduced system:(25) (cid:20) A + τ n L τ n L τ n C F G (cid:21) r δV n +1 h,r ( µ ) = − τ n R r ( V nh,r ( µ )) . Remark . The initial conditions for the ROM system of equation (25) are obtained performing aGalerkin projection of the initial full-order condition u ( · , µ ), w ( · , µ ) onto the POD basis spaces.For an efficient ROM basis construction, we follow some ideas of the authors as demonstrated in [52].In particular, concerning the full order method snapshots extension, and the extension of the solutionto the surrogate domain into the ghost area, we use the solution values as they have been computedusing the cut finite element method smooth mapping from the true to the unfitted mesh domain.This approach allows a smooth extension of the boundary solution to the neighboring ghost elementswith values which are decreasing smoothly to zero, see for instance the zoomed image in Figure 3.This approach provides a regular solution in the background domain and permits, therefore, theconstruction of a “functional” reduced basis. ROM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS (i) (ii)
Figure 2.
The background mesh (i), and a sketch of the embedded domain and theparameters considered in the numerical examples (ii).
Figure 3.
A zoom onto the embedded circular domain of the numerical exampleshows the smoothing natural smooth extension procedure employed by the cut finiteelement method inside the ghost area.4.
Numerical experiments
In the present section, we test the presented methodology considering numerical experiments for theevolutionary Cahn-Hilliard system while natural homogeneous Neumann boundary conditions, and/orzero Dirichlet ones are present. We start by testing the robustness of the full order model (FOM) fortwo classical benchmark test cases and we continue with a numerical example in which geometricalparametrization for the embedded domain is considered. The background domain in all experiments isthe rectangle [ − . , . × [ − . , . γ , γ and γ according to formula (7), are taken as γ = 2, γ = 9 and γ = 4, with the two fluids interfacesize ε = 10 − . The results for all test problems have been obtained with mesh size h = 1 /
48 unlessotherwise stated, Nitsche parameter a N = 10 and jump stabilization parameter a = 0 . .4.1. Robustness of the FOM solver.
In this first numerical example, we test the validity of the fullorder solver and the IMEX method. Two benchmark tests are examined, [43, 34, 79], a cross-shapedand an ellipse-shaped phase-field interface applying time step τ n = O (10 − ), for 10000 and 7000 timesteps respectively and without any embedded geometry, see Figure 2 (i). For the first experiment, weconsider the initial state u ( x, y ) = . , if 5 | ( y − . − ( x − . | + | ( x − . − ( y − . | ≤ , . , if 5 | ( x − . − ( y − . | + | ( y − . − ( x − . | ≤ , , otherwise , and for the second ellipse-shaped interface the initial data u ( x, y ) = − tanh((( x − . / . ) + (( y − . / . ) − . In both cases, with elliptic, and with the additional difficulty of presenting sharp corners in a cross-shaped interface case, we observe evolution toward to a circular interface, see Figures 4, 5, that validatesthe efficiency of the FOM methodology. For these benchmark tests, we also refer to the work [17].
Remark . Moreover, we emphasize that the cut elements FEM stabilization robustness is addition-ally verified using the conservation of mass test, as it is illustrated in the next paragraph numericalexperiment. In particular, in Figure 12, we can easily notice that the FOM solution and in particular
OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 9
Figure 4.
Testing the IMEX method: Cross type initial data and the phase fieldresults for times t = [0 , , , , , τ n which visualizes that it evolvesefficiently to a circular interface. Figure 5.
Testing the IMEX method: Ellipse type initial data and the phase-fieldresults for times t = [0 , , τ n which visualizes that it evolves efficiently to acircular interface.the mass preserves very well, as time passes fulfilling the conservation of mass property (see Equation3).4.2. Geometrical parametrization.
The numerical examples consider a geometrical parametriza-tion on the embedded domain. The embedded domain is in fact parametrized through µ according tothe expression: x + y ≤ µ / , where the parameter µ describes the diameter size of the circular embedded domain located on thecenter (0 ,
0) of the background domain, see e.g. Figure 2 (ii). The experiments focus on the initialevolution period when the phase-field changes fast, namely in the time interval 0 − τ n for time stepsize τ n = O (10 − ). We tested pseudo-random initial concentration u data and zero initial data for theauxiliary variable w . For all parameters, we experimented the same pseudo-random initial state appliedon the whole background mesh and afterward restricted onto the active parametrized geometry. In thefirst experiment the boundaries are free (open boundary condition) everywhere, while in the secondone we consider Dirichlet only for the embedded geometry, [61]. Though, only the cylinder is treatedas embedded. Linear P × P , T ] = [0 , τ n ], someone can notice fast phase-field changes startingfrom pseudo-random initial data, both challenging for the ROM construction. Nevertheless, as timepasses these changes are weakened, see e.g. Figure 6 for a visualization of the evolution of an openand a Dirichlet type embedded boundary.For the reduced basis solution, the ROM has been trained onto 900 parameter samples and testedonto 30 samples in a parameter range µ test ∈ [0 . , .
48] chosen randomly inside the parameter space.From these snapshots selection, a ROM basis has been derived with a POD procedure. The visual-ization of the first six basis components can be seen in Figure 7. To test the accuracy of the ROMwe compared its results on 30 additional samples which were not used to create the ROM and wereselected randomly within the aforementioned range.
Figure 6.
The full order solver and the concentration field (whole background ge-ometry) for a fixed parameter µ test = 0 .
436 at times t = [1 , , , , , τ n andmesh size h = 1 / OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 11
Figure 7.
The first 6 basis components for concentration u considering Neumannand Dirichlet geometrical parametrization for µ ∈ [0 . , . Figure 8.
Results for the geometrical parametrized embedded circle with diameter µ test = 0 .
436 and open boundary . In the first, second and third column we reportthe full-order solution, the reduced order solution and the absolute error plots for theconcentration field. Each row corresponds to a different time t = [10 , , , , τ n . OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 13
Figure 9.
Results for the geometrical parametrized embedded circle with diame-ter µ test = 0 .
436 and
Dirichlet embedded boundary . In the first, second and thirdcolumn we report the full-order solution, the reduced order solution and the abso-lute error plots for the concentration field. Each row corresponds to a different time t = [10 , , , , τ n . Snapshots for: 900 train parametersModes Relative error for Neumann case Relative error for Dirichlet case( N ) Concentration u Potential w Concentration u Potential w Table 1.
Geometrical parametrization and the relative error between the full-ordersolution and the reduced basis solution for the concentration and potential component.Results are reported for different dimensions of the reduced basis spaces.(i) N N / m a x Potential w (Neumann)Concentration u (Neumann)
Number of Modes N M e a n r e l a t i v e e rr o r Potential w (Neumann)Concentration u (Neumann) (ii)(iii) N N / m a x Potential w (Dirichlet)Concentration u (Dirichlet) 0 10 20 30 40
Number of Modes N M e a n r e l a t i v e e rr o r Potential w (Dirichlet)Concentration u (Dirichlet) (iv)
Figure 10.
Geometrical parametrization and visualization of the results. The leftplots depict the eigenvalues decay with respect to the numbers of the modes for theconcentration and potential variable. On the right plots, we visualize the relativeerrors of the reduced-order problem for one random parameter value and for variousnumber of modes.As it is displayed in Figures 8 and 9 –Neumann and Dirichlet embedded boundary experimentsrespectively– and with a minds eye comparison of the FOM and ROM solutions in first and secondcolumn, with a brief look they seem identical for every row associated to the time instances t =[10 , , , , τ n . In the third column and looking from a more detailed point of view, the absoluteerror in each point of the geometry domain is visualized for the open boundary and the Dirichlet case.In Table 1 and again for both types of boundaries, the relative errors for the concentration phasefield for various number of basis components are reported and their graph can be seen in Figure 10 aswell as the eigenvalues and their decay which have been used for the ROM. Thereafter, for a better OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 15
Snapshots for: 900 train parametersModes execution times Savings( N ) (t) ( t FOM − t RB ) /t FOM
Table 2.
Execution time at the reduced order level and per cent (%) savings. Thecomputational time includes the projection of the full order matrices, the executiontime of the online solver and the resolution of the reduced problem. Times are forthe resolution of one random value of the input parameter. The time execution at fullorder method level (FOM) is equal to ≈ . Time r e l a t i v e e rr o r ( t ) f o r t h e N e u mm a n B C e x p . Potential wConcentration u 0 20 40 60 80 100
Time r e l a t i v e e rr o r ( t ) f o r t h e D i r i c h l e t B C e x p . Potential wConcentration u (iii)(ii) N e x e c u t i o n t i m e s reducedtruth Figure 11.
Time related results, and visualization for the Neumann and Dirichletcase. The left and right plots depict the relative errors with respect the time for theconcentration and potential variable and 45 modes. On the bottom plot, we visualizethe execution times of the reduced-order problem for one random parameter value andfor various number of modes.understanding, the relative errors evolution with respect to time, for a 45 modes test, is demonstratedin Table 3 and visualized in Figure 11 (i), (iii). Finally in Table 2, and in Figure 11 (ii) we reportthe execution times for various number of basis components including the projection of the full ordermatrices, the execution time of the online solver and the determination of the reduced problem andwe compare them with the time execution at the full order level, namely 9 . Some comments.
Numerical experiments consider open boundary conditions as well as Dirichletones. The tests clearly indicate that an efficient orthogonal decomposition projection-based reduced-order model can be derived over a full-order cut finite element method solver for the challenging (forboth the full and the reduced level, see also the uncomfort basis components in Figure 7) nonlinear
Snapshots for: 900 train parameterstime ( i × τ n ) Relative error (Neumann) Relative error (Dirichlet) i Concentration u Potential w Concentration u Potential w1 0.00342 0.01692 0.00357 0.0094410 0.00472 0.01016 0.00699 0.0151720 0.00510 0.01613 0.00894 0.0225030 0.00952 0.02374 0.01330 0.0375140 0.01634 0.03580 0.01938 0.0507850 0.02184 0.05944 0.02752 0.0711760 0.03639 0.12172 0.03086 0.1045270 0.04636 0.11315 0.03381 0.1075080 0.04611 0.11615 0.03606 0.1258590 0.04787 0.10791 0.03856 0.13514100 0.05181 0.11652 0.04169 0.11201
Table 3.
The L relative error norm is reported overtime for the concentration andpotential field for the Neumann and Dirichlet type of boundaries experiments. TheROMs have been obtained with 45 modes for the concentration for both cases.Cahn–Hilliard system. Another aspect is that we avoid the costly traditional approach of geometrictransformation in a reference domain and we rely only on an appropriate smooth enough backgroundgrid. Last Section tests have clearly shown that accurate results can be obtained for this phase-fieldsystem at a reduced level.5. Concluding remarks and future developments
We conclude this work by noting that the above approach and the combination of unfitted meshfinite element methods with an embedded POD basis and IMEX type discretization, by using linearpolynomials, imply good reduced basis approximation properties for a Cahn–Hilliard fourth-orderdiffusion nonlinear PDE system. Considering the reduced-order approximation, the background meshapproach appears beneficial even for nonsmooth pseudo-random initial data. As expected, and referringto previous related authors’ works, [55, 54], an increased error is noticed onto the Nitsche embeddedboundary interface in the reduced level, which will be studied elsewhere. We underline the significantexecution time reduction considering the projection of the full order matrices, through the reducedexecution time of the online solver and the resolution of the reduced problem, by capturing efficientlythe full order solution information in a reduced level solution. As a perspective, we mention theconstruction of higher-order IMEX methods, more efficient methodologies for the affine decompositionof the discretization differential operator, and the investigation of the applicability of well-knownhyper reduction techniques, such as the empirical interpolation method, [40], in the context of theaforementioned embedded reduced-order basis method. Of future interest also parabolic nonlinearpartial differential equations set in a more general framework and to test snapshots transportationtechniques as presented in [52].
Acknowledgments
This work is supported by the European Research Council Executive Agency by means of the H2020ERC Consolidator Grant project AROMA-CFD “Advanced Reduced Order Methods with Applicationsin Computational Fluid Dynamics” - GA 681447, (PI: Prof. G. Rozza), FARE-X-AROMA-CFD projectby MIUR, INdAM-GNCS 2018 and 2019 and by project FSE - European Social Fund - HEaD “HigherEducation and Development” SISSA operazione 1, Regione Autonoma Friuli - Venezia Giulia, theHellenic Foundation for Research and Innovation (HFRI) and the General Secretariat for Researchand Technology (GSRT), under grant agreement No[1115] and National Infrastructures for Researchand Technology S.A. (GRNET S.A.) in the National HPC facility - ARIS - under project ID pa190902.The authors would like also to thank Dr Andrea Mola for useful instructions regarding the cutfemstabilization, and Dr Francesco Ballarin for fruitful discussions for the reduced order part. Numericalsimulations have been obtained with ngsxfem/ngsolve [85, 86, 84] and
RBniCS [87]; we acknowledgedevelopers and contributors of each of the aforementioned software packages.
OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 17
Time M a ss ( t ) f o r t h e N e u mm a n B C e x p . Time M a ss ( t ) f o r t h e D i r i c h l e t B C e x p . Figure 12.
The conservation of mass for Neumann and Dirichlet embedded bound-aries, evolutionary with respect to time, together with its RB approximation is re-ported and verified in FOM and ROM level for various number of modes as consid-ered in the numerical examples, for the parameter value µ = 0 .
436 and time instances t = nτ n , for n = 1 , ..., References [1] Agosti, A., Antonietti, P.F., Ciarletta, P., Grasselli, M., Verani, M.: A Cahn-Hilliard–type equation with applicationto tumor growth dynamics. Mathematical Methods in the Applied Sciences (18), 7598–7626 (2017)[2] Alikakos, N., Fusco, G., Smyrnelis, P.: Elliptic Systems of Phase Transition Type. 91. Monograph in the seriesProgress in Nonlinear Differential Equations and Their Applications, Birkhauser (2018)[3] Alpak, F.O., Riviere, B., Frank, F.: A phase-field method for the direct simulation of two-phase flows in pore-scalemedia using a non-equilibrium wetting boundary condition. Computational Geosciences (5), 881–908 (2016)[4] Anderson, D.M., McFadden, G.B., Wheeler, A.A.: Diffuse-interface methods in fluid mechanics. Annual Review ofFluid Mechanics (1), 139–165 (1998)[5] Antonopoulou, D.C., Karali, G., Millet, A.: Existence and regularity of solution for a stochastic Cahn–Hilliard/Allen-Cahn equation with unbounded noise diffusion (2013)[6] Antonopoulou, D.C., Karali, G., Millet, A.: Existence and regularity of solution for a stochastic Cahn–Hilliard/Allen–Cahn equation with unbounded noise diffusion. Journal of Differential Equations (3), 2383 –2417 (2016)[7] Balajewicz, M., Farhat, C.: Reduction of nonlinear embedded boundary models for problems with evolving inter-faces. Journal of Computational Physics , 489–504 (2014) [8] Ballarin, F., Manzoni, A., Quarteroni, A., Rozza, G.: Supremizer stabilization of POD-Galerkin approximationof parametrized steady incompressible Navier–Stokes equations. International Journal for Numerical Methods inEngineering (5), 1136–1161 (2015)[9] Benner, P., Ohlberger, M., Patera, A., Rozza, G., Urban, K.: Model Reduction of Parametrized Systems, MS&Aseries , vol. 17. Springer (2017)[10] Bertozzi, A.L., Esedoglu, S., Gillette, A.: Inpainting of Binary Images Using the Cahn–Hilliard Equation. Trans.Img. Proc. (1), 285–291 (2007)[11] Blank, L., Garcke, H., Sarbu, L., Srisupattarawanit, T., Styles, V., Voigt, A.: Phase-field Approaches to StructuralTopology Optimization, pp. 245–256. Springer Basel, Basel (2012)[12] Bosch, J.: Fast Iterative Solvers for Cahn-Hilliard Problems. Ph.D. thesis, Otto-von-Guericke Universitat, Magde-burg (2016)[13] Burman, E.: Ghost penalty. Comptes Rendus Mathematique (21), 1217 – 1220 (2010)[14] Burman, E., Hansbo, P.: Fictitious domain finite element methods using cut elements: II. A stabilized Nitschemethod. Applied Numerical Mathematics (6), 2837–2862 (2011)[15] Burman, E., Hansbo, P.: Fictitious domain methods using cut elements: III. A stabilized Nitsche method for Stokes’problem. ESAIM: M2AN (5-8), 859–874 (2014)[16] Cahn, J.W., Hilliard, J.E.: Free Energy of a Nonuniform System. I. Interfacial Free Energy. The Journal of ChemicalPhysics (2), 258–267 (1958)[17] Chave, F., Di Pietro, D., Marche, F., Pigeonneau, F.: A Hybrid High-Order Method for the Cahn–Hilliard problemin Mixed Form. SIAM Journal on Numerical Analysis (3), 1873–1898 (2016)[18] Cherfils, L., Fakih, H., Miranville, A.: A Complex Version of the Cahn–Hilliard Equation for Grayscale ImageInpainting. Multiscale Modeling & Simulation (1), 575–605 (2017)[19] Chinesta, F., Huerta, A., Rozza, G., Willcox, K.: chap. Model Reduction Methods, Encyclopedia of ComputationalMechanics, Second Edition, pp. 1–36. John Wiley & Sons (2017)[20] Chinesta, F., Ladeveze, P., Cueto, E.: A Short Review on Model Order Reduction Based on Proper GeneralizedDecomposition. Archives of Computational Methods in Engineering (4), 395 (2011)[21] Choksi, R., Peletier, M., Williams, J.: On the Phase Diagram for Microphase Separation of Diblock Copolymers:An Approach via a Nonlocal Cahn–Hilliard Functional. SIAM Journal on Applied Mathematics (6), 1712–1738(2009)[22] Chrysafinos, K., Karatzas, E.N.: Error Estimates for Discontinuous Galerkin Time-Stepping Schemes for RobinBoundary Control Problems Constrained to Parabolic PDEs. SIAM Journal on Numerical Analysis (6), 2837–2862 (2014)[23] Chrysafinos, K., Karatzas, E.N.: Symmetric error estimates for discontinuous Galerkin time-stepping schemes foroptimal control problems constrained to evolutionary Stokes equations. Computational Optimization and Applica-tions (3), 719–751 (2015)[24] Claus, S., Kerfriden, P.: A cutfem method for two-phase flow problems. Computer Methods in Applied Mechanicsand Engineering , 185 – 206 (2019)[25] Colli, P., Farshbaf-Shaker, M., Gilardi, G., Sprekels, J.: Optimal Boundary Control of a Viscous Cahn–HilliardSystem with Dynamic Boundary Condition and Double Obstacle Potentials. SIAM Journal on Control and Opti-mization (4), 2696–2721 (2015)[26] De Groot, S., Mazur, P.: Non-equilibrium thermodynamics. (1962), Dover edition (2013)[27] Dumon, A., Allery, C., Ammar, A.: Proper general decomposition (PGD) for the resolution of Navier–Stokesequations. Journal of Computational Physics (4), 1387–1407 (2011)[28] Elliott, C., Larsson, S.: Error estimates with smooth and nonsmooth data for a finite element method for theCahn-Hilliard equation. Math. Comp. (S33-S36), 603–630 (1992)[29] Elliott, C.M.: The Cahn-Hilliard Model for the Kinetics of Phase Separation, pp. 35–73. Birkh¨auser Basel, Basel(1989)[30] Elliott, C.M., French, D.A., Milner, F.A.: A second order splitting method for the Cahn-Hilliard equation. Nu-merische Mathematik (5), 575–590 (1989)[31] Elliott, C.M., Songmu, Z.: On the Cahn-Hilliard equation. Archive for Rational Mechanics and Analysis (4),339–357 (1986)[32] Embar, A., Dolbow, J., Harari, I.: Imposing dirichlet boundary conditions with Nitsche’s method and spline-basedfinite elements. International Journal for Numerical Methods in Engineering (7), 877–898 (2010)[33] Furihata, D., Kov`acs, M., Larsson, S., Lindgren, F.: Strong Convergence of a Fully Discrete Finite Element Ap-proximation of the Stochastic Cahn–Hilliard Equation. SIAM Journal on Numerical Analysis (2), 708–731 (2018)[34] Gouden`ege, L., Martin, D., Vial, G.: High Order Finite Element Calculations for the Cahn-Hilliard Equation.Journal of Scientific Computing (2), 294–321 (2012)[35] Gr¨aßle, C., Hinze, M., Scharmacher, N.: POD for Optimal Control of the Cahn-Hilliard System Using SpatiallyAdapted Snapshots. In: F.A. Radu, K. Kumar, I. Berre, J.M. Nordbotten, I.S. Pop (eds.) Numerical Mathematicsand Advanced Applications ENUMATH 2017, pp. 703–711. Springer International Publishing (2019)[36] Grepl, M., Maday, Y., Nguyen, N., Patera, A.: Efficient reduced-basis treatment of nonaffine and nonlinear partialdifferential equations. ESAIM: M2AN (3), 575–605 (2007)[37] Grepl, M., Patera, A.: A posteriori error bounds for reduced-basis approximations of parametrized parabolic partialdifferential equations. ESAIM: M2AN (1), 157–181 (2005)[38] Gurtin, M.E., Polignone, D., Vinals, J.: Two-phase binary fluids and immiscible fluids described by an orderparameter. Mathematical Models and Methods in Applied Sciences (06), 815–831 (1996)[39] Haasdonk, B., Ohlberger, M.: Reduced basis method for finite volume approximations of parametrized linearevolution equations. Mathematical Modelling and Numerical Analysis (2), 277–302 (2008) OM FOR GEOMETRICALLY PARAMETRIZED CAHN-HILLIARD SYSTEM BASED ON CUT ELEMENTS 19 [40] Haasdonk, B., Ohlberger, M., Rozza, G.: A reduced basis method for evolution schemes with parameter-dependentexplicit operators. Electron. Trans. Numer. Anal. , 145–161 (2008)[41] Harari, I., Grosu, E.: A unified approach for embedded boundary conditions for fourth-order elliptic problems.International Journal for Numerical Methods in Engineering (7), 655–675 (2015)[42] Hesthaven, J., Rozza, G., Stamm, B.: Certified Reduced Basis Methods for Parametrized Partial DifferentialEquations. SpringerBriefs in Mathematics. Springer International Publishing (2016)[43] Hinterm¨uller, M., Hinze, M., Kahle, C.: An Adaptive Finite Element Moreau-Yosida-based Solver for a CoupledCahn-Hilliard/Navier-Stokes System. J. Comput. Phys. (C), 810–827 (2013)[44] Hinterm¨uller, M., Keil, T., Wegner, D.: Optimal Control of a Semidiscrete Cahn–Hilliard–Navier–Stokes Systemwith Nonmatched Fluid Densities. SIAM Journal on Control and Optimization (3), 1954–1989 (2017)[45] Hintermuller, M., Wegner, D.: Distributed Optimal Control of the Cahn–Hilliard System Including the Case of aDouble-Obstacle Homogeneous Free Energy Density. SIAM Journal on Control and Optimization (1), 388–418(2012)[46] Hinze, M., Kahle, C.: A Nonlinear Model Predictive Concept for Control of Two-Phase Flows Governed by theCahn-Hilliard Navier-Stokes System. In: D. H¨omberg, F. Tr¨oltzsch (eds.) System Modeling and Optimization, pp.348–357. Springer Berlin Heidelberg, Berlin, Heidelberg (2013)[47] Israelachvili, J.N.: Intermolecular and Surface Forces. Elsevier (2011)[48] Jeong, D., Kim, J.: Microphase separation patterns in diblock copolymers on curved surfaces using a nonlocalCahn-Hilliard equation. The European Physical Journal E (11), 117 (2015)[49] Junseok, K., Seunggyu, L., Yongho, C., Seok-Min, L., Darae, J.: Basic Principles and Practical Applications of theCahn–Hilliard Equation. Mathematical Problems in Engineering (1), 79–141 (2016)[50] Kalashnikova, I., Barone, M.F.: On the stability and convergence of a Galerkin reduced order model (ROM) ofcompressible flow with solid wall and far-field boundary treatment. International Journal for Numerical Methods inEngineering (10), 1345–1375 (2010)[51] Karali, G., Nagase, Y.: On the existence of solution for a Cahn–Hilliard/Allen–Cahn equation. Discrete and Con-tinuous Dynamical Systems - S (1937-1632 2014 1 127), 127 (2014)[52] Karatzas, E.N., Ballarin, F., Rozza, G.: Projection-based reduced order models for a cut finite element method inparametrized domains. Computers & Mathematics with Applications (3), 833 – 851 (2020)[53] Karatzas, E.N., Nonino, M., Ballarin, F., Rozza, G.: A Reduced order cut finite element basis for stationary andevolutionary geometrically parameterized Navier–Stokes systems (2020). In preparation.[54] Karatzas, E.N., Stabile, G., Atallah, N., Scovazzi, G., Rozza, G.: A Reduced Order Approach for the EmbeddedShifted Boundary FEM and a Heat Exchange System on Parametrized Geometries In: Fehr J., Haasdonk B. (eds)IUTAM Symposium on Model Order Reduction of Coupled Systems, Stuttgart, Germany, May 22–25, 2018. IUTAMBookseries, vol 36. Springer, Cham (2020)[55] Karatzas, E.N., Stabile, G., Nouveau, L., Scovazzi, G., Rozza, G.: A reduced basis approach for PDEs onparametrized geometries based on the shifted boundary finite element method and application to a Stokes flow.Computer Methods in Applied Mechanics and Engineering , 568 – 587 (2019)[56] Karatzas, E.N., Stabile, G., Nouveau, L., Scovazzi, G., Rozza, G.: A reduced-order shifted boundary method forparametrized incompressible Navier–Stokes equations. Computer Methods in Applied Mechanics and Engineering , 113–273 (2020)[57] Katsouleas, G., Karatzas, E.N., Travlopanos, F.: Cut finite element error estimates for a class of nonlinear ellipticPDEs. pp. 1–6. Loughborough University, doi: 10.17028/rd.lboro.12154854.v1, extended version at arXiv:2003.06489(2020)[58] Kunisch, K., Volkwein, S.: Galerkin proper orthogonal decomposition methods for a general equation in fluiddynamics. SIAM Journal on Numerical Analysis (2), 492–515 (2002)[59] Lehrenfeld, C., Reusken, A.: L2-error analysis of an isoparametric unfitted finite element method for elliptic interfaceproblems. pp. 85–99. Journal of Numerical Mathematics (2019)[60] Li, C., Qin, R., Ming, J., Wang, Z.: A discontinuous Galerkin method for stochastic Cahn–Hilliard equations.Computers & Mathematics with Applications (6), 2100 – 2114 (2018). 2nd Annual Meeting of SIAM CentralStates Section, September 30-October 2, 2016[61] Li, Y., Jeong, D., Shin, J., Kim, J.: A conservative numerical method for the Cahn–Hilliard equation with Dirichletboundary conditions in complex domains. Computers & Mathematics with Applications (1), 102 – 115 (2013)[62] Novick-Cohen, A., Segel, L.A.: Nonlinear aspects of the Cahn-Hilliard equation. Physica D: Nonlinear Phenomena (3), 277 – 298 (1984)[63] Quarteroni, A., Manzoni, A., Negri, F.: Reduced Basis Methods for Partial Differential Equations, vol. 92. UNI-TEXT/La Matematica per il 3+2 book series, Springer International Publishing (2016)[64] Regazzoni, F., Parolini, N., Verani, M.: Topology optimization of multiple anisotropic materials, with applicationto self-assembling diblock copolymers. Computer Methods in Applied Mechanics and Engineering , 562 – 596(2018)[65] Reshma, S., Hansa J Thattil, H.: Inpainting of Binary Images Using the Cahn-Hilliard Equation. InternationalJournal of Computer Science Engineering and Technology (11), 296–300 (2014)[66] Rocca, E., Sprekels, J.: Optimal Distributed Control of a Nonlocal Convective Cahn–Hilliard Equation by theVelocity in Three Dimensions. SIAM Journal on Control and Optimization (3), 1654–1680 (2015)[67] Rokhzadi, A.: IMEX and Semi-Implicit Runge-Kutta Schemes for CFD Simulations. Ph.D. thesis, Civil EngineeringDepartment, Faculty of Engineering, University of Ottawa (2018)[68] Rozza, G.: Reduced Basis Methods for Elliptic Equations in subdomains with A-Posteriori Error Bounds andAdaptivity. App. Num. Math. (4), 403–424 (2005)[69] Rozza, G.: Reduced basis methods for Stokes equations in domains with non-affine parameter dependence. Com-puting and Visualization in Science (1), 23–35 (2009) [70] Rozza, G., Huynh, D., Patera, A.: Reduced basis approximation and a posteriori error estimation for affinelyparametrized elliptic coercive partial differential equations: Application to transport and continuum mechanics.Archives of Computational Methods in Engineering (3), 229–275 (2008)[71] Rozza, G., Huynh, D.B.P., Manzoni, A.: Reduced basis approximation and a posteriori error estimation for Stokesflows in parametrized geometries: Roles of the inf-sup stability constants. Numerische Mathematik (1), 115–152(2013)[72] Rozza, G., Veroy, K.: On the stability of the reduced basis method for Stokes equations in parametrized domains.Computer Methods in Applied Mechanics and Engineering (7), 1244–1260 (2007)[73] Schott, B.: Stabilized Cut Finite Element Methods for Complex Interface Coupled Flow Problems. Ph.D. thesis,Technische Universit¨at M¨unchen (TUM) (2016)[74] Shenyang, H.: Phase-field Models of Microstructure Evolution in a System with Elastic Inhomogeneity and Defects.Ph.D. thesis, Pennsylvania State University, Department of Materials Science and Engineering (2004)[75] Veroy, K., Prud’homme, C., Patera, A.: Reduced-basis approximation of the viscous Burgers equation: rigorous aposteriori error bounds. Comptes Rendus Mathematique (9), 619–624 (2003)[76] Wells, G.N., Kuhl, E., Garikipati, K.: A discontinuous Galerkin method for the Cahn–Hilliard equation. Journal ofComputational Physics (2), 860 – 877 (2006)[77] Welper, G.: Optimal treatment for a phase field system of Cahn-Hilliard type modeling tumor growth by asymptoticscheme (2019). ArXiv:1902.01079v2[78] Wodo, O., Ganapathysubramanian, B.: Computationally efficient solution to the Cahn–Hilliard equation: Adaptiveimplicit time schemes, mesh sensitivity analysis and the 3D isoperimetric problem. Journal of Computational Physics (15), 6037 – 6060 (2011)[79] Xu, M., Guo, H., Zou, Q.: Hessian recovery based finite element methods for the Cahn-Hilliard equation. Journalof Computational Physics , 524 – 540 (2019)[80] Zhang, X., Li, H., Liu, C.: Optimal Control Problem for the Cahn–Hilliard/Allen–Cahn Equation with StateConstraint. Applied Mathematics & Optimization (2018)[81] Zhao, X., Liu, C.: Optimal Control for the Convective Cahn–Hilliard Equation in 2D Case. Applied Mathematics& Optimization (1), 61–82 (2014)[82] Zhao, Y., Schillinger, D., Xu, B.X.: Variational boundary conditions based on the Nitsche method for fitted andunfitted isogeometric discretizations of the mechanically coupled Cahn-Hilliard equation. Journal of ComputationalPhysics , 177 – 199 (2017)[83] Zhou, S., Wang, M.Y.: Multimaterial structural topology optimization with a generalized Cahn–Hilliard model ofmultiphase transition. Structural and Multidisciplinary Optimization (2), 89 (2006)[84] Sch¨oberl, J.: C++11 implementation of finite elements in NGSolve. Institute for Analysis and Scientific Computing,Vienna University of Technology (2014)[85] ngsxfem – Add-On to NGSolve for unfitted finite element discretizations. https://github.com/ngsxfem/ngsxfemhttps://github.com/ngsxfem/ngsxfem