An efficient and convergent finite element scheme for Cahn--Hilliard equations with dynamic boundary conditions
aa r X i v : . [ m a t h . NA ] A ug AN EFFICIENT AND CONVERGENT FINITE ELEMENT SCHEMEFOR CAHN–HILLIARD EQUATIONS WITH DYNAMIC BOUNDARYCONDITIONS
STEFAN METZGER
Department of Applied Mathematics, Illinois Institute of Technology, Chicago IL, 60616,USA
Abstract.
The Cahn–Hilliard equation is a widely used model that describes amongstothers phase separation processes of binary mixtures or two-phase-flows. In the recentyears, different types of boundary conditions for the Cahn–Hilliard equation were pro-posed and analyzed. In this publication, we are concerned with the numerical treatmentof a recent model which introduces an additional Cahn–Hilliard type equation on theboundary as closure for the Cahn–Hilliard equation [C. Liu, H. Wu, Arch. Ration.Mech. An., 2019]. By identifying a mapping between the phase-field parameter and thechemical potential inside of the domain, we are able to postulate an efficient, uncondi-tionally energy stable finite element scheme. Furthermore, we establish the convergenceof discrete solutions towards suitable weak solutions of the original model, which servesalso as an additional pathway to establish existence of weak solutions. Introduction
Different approaches to model the hydrodynamics of fluid mixtures have been widelyused in literature.. In addition to the conventional sharp interface models which consistof separate hydrodynamic systems for each component of the mixture, there are diffuseinterface models. In these models, the hyper-surface description of the fluid-fluid interfaceis replaced by a small transition region, where mixing of the macroscopically immisciblefluids is allowed, which leads to a smooth transition between the pure phases. In its easiestform, a diffuse interface model for two phases in a domain Ω with boundary Γ = ∂ Ω reads ∂ t φ = m ∆ µ in Ω , (1.1a) µ = − δσ ∆ φ + δ − σF ′ ( φ ) in Ω , (1.1b) ∇ µ · n = 0 on Γ , (1.1c) ∇ φ · n = 0 on Γ , (1.1d) E-mail address : [email protected] . Date : August 19, 2019.2010
Mathematics Subject Classification.
Key words and phrases.
Cahn–Hilliard, dynamic boundary conditions, finite elements, convergence.THIS WORK WAS FUNDED BY THE NSF THROUGH GRANT NUMBER NSF-DMS1759536. in combination with suitable initial conditions for the phase-field parameter φ . Here, n is the outer normal vector of the domain Ω , m > is a mobility constant, σ > is aparameter related to surface tension, and the parameter δ > prescribes the width of theinterfacial area. The chemical potential µ given as the first variation of the free energy E Ω ( φ ) := σδ ˆ Ω 12 |∇ φ | + σδ − ˆ Ω F ( φ ) , (1.2)where F is a double-well potential with minima in φ = ± representing the pure fluidphases. Typical choices for F are F ( φ ) := θ (1 + φ ) log (1 + φ )+ θ (1 − φ ) log (1 − φ ) − θ c φ with < θ < θ c and F ( φ ) := ( φ − . The boundary condition (1.1c) states that thereis no flux across Γ , i.e. ´ Ω φ is conserved. The second boundary condition (1.1d) indicatesthat the fluid-fluid interface, i.e. the zero level set of φ , intersects the boundary Γ at astatic contact angle of π . This can be interpreted as neglecting the interactions betweenthe fluids and the walls of the surrounding container.Although (1.1) satisfies the energybalance equation E Ω ( φ ) (cid:12)(cid:12) T + ˆ T ˆ Ω m |∇ µ | = E Ω ( φ ) (cid:12)(cid:12) , (1.3)the boundary condition (1.1d) imposing a static contact angle is considered a majorflaw and there are several attempts to improve this boundary condition. In [16], it wassuggested to include the wetting energy ´ Γ γ fs ( φ ) to the system’s energy and to replace(1.1d) by σδ ∇ φ · n = − α∂ t φ − γ ′ fs ( φ ) . (1.4)Here, γ fs interpolates between liquid-solid interfacial energies of the two fluid phases andprescribes the stationary contact angle via Young’s formula. For α > , (1.4) allows alsofor hysteresis effects in the evolution of the contact angles and changes the energy balance(1.3) to E Ω ( φ ) (cid:12)(cid:12) T + ˆ ∂ Ω γ fs ( φ ) (cid:12)(cid:12)(cid:12)(cid:12) T + ˆ T ˆ Ω m |∇ µ | + α ˆ ∂ Ω | ∂ t φ | = E Ω ( φ ) (cid:12)(cid:12) + ˆ ∂ Ω γ fs ( φ ) (cid:12)(cid:12)(cid:12)(cid:12) . (1.5)Other approaches for improving (1.1d) include – among others – using an Allen–Cahn-typeequation to replace (1.1d). The resulting boundary condition reads δ Γ ∂ t φ = κδ Γ ∆ Γ φ − δ − G ′ ( φ ) − δσ ∇ φ · n on Γ (1.6)with ∆ Γ denoting the Laplace-Beltrami operator describing surface diffusion on Γ . Thistype of dynamic boundary condition was suggested in [11] and analyzed in [5]. Formally,it can be seen as an L -gradient flow of the surface free energy E Γ ( φ ) := κδ Γ ˆ Γ 12 |∇ Γ φ | + δ − ˆ Γ G ( φ ) , (1.7)where ∇ Γ denotes the tangential (surface) gradient operator defined on Γ and G is asurface potential function. In particular, the resulting model satisfies E Ω ( φ ) (cid:12)(cid:12) T + E Γ ( φ ) (cid:12)(cid:12) T + ˆ T ˆ Ω m |∇ µ | + ˆ T ˆ Γ | ∂ t φ | = E Ω ( φ ) (cid:12)(cid:12) + E Γ ( φ ) (cid:12)(cid:12) . (1.8) From a formal point of view, one may set δ Γ = 1 , κ = 0 , and G ≡ γ fs in order to recover(1.4).Recently, a new type of dynamic boundary condition was derived by C. Liu and H. Wu[13]. Using a variational approach with different flow maps for Ω and Γ , they deriveda model closing the Cahn–Hilliard equation in Ω with an additional Cahn–Hilliard-typeequation on Γ . The derived model assumes that φ is continuous on Ω and reads ∂ t φ = m ∆ µ in Ω , (1.9a) µ = − δσ ∆ φ + δ − σF ′ ( φ ) in Ω , (1.9b) ∇ µ · n = 0 on Γ , (1.9c) ∂ t φ = m Γ ∆ Γ µ Γ on Γ , (1.9d) µ Γ = − δ Γ κ ∆ Γ φ + δ − G ′ ( φ ) + δσ ∇ φ · n on Γ . (1.9e)This model is derived from physical principles and satisfies conservation of mass on Ω and Γ , balance of forces, and dissipation of free energy. In particular, model (1.9) satisfies theenergy equation E Ω ( φ ) (cid:12)(cid:12) T + E Γ ( φ ) (cid:12)(cid:12) T + ˆ T ˆ Ω m |∇ µ | + ˆ T ˆ Γ m Γ |∇ Γ µ Γ | = E Ω ( φ ) (cid:12)(cid:12) + E Γ ( φ ) (cid:12)(cid:12) . (1.10)Constructing solutions of a regularized system, where (1.9b) and (1.9e) are extended by α∂ t φ , and taking the limit α ց allowed them to prove existence and uniqueness of weakand strong solutions.The authors of [7] interpreted model (1.9) as a gradient flow equation of the total freeenergy E Ω ( φ ) + E Γ ( φ ) and used this structure for their proof of existence and uniquenessof weak solutions.The numerical treatment of the Cahn–Hilliard equation and its variants – often in com-bination with Navier–Stokes-equations – was intensely discussed through the last years.Consequently, there are various different discretization techniques at hand, which transferthe energy stability (1.3) to a discrete setting. These techniques include approaches basedon convex-concave splittings of the energy (cf. [20, 17]) or the double-well potential (cf.[10, 8, 6, 9]), stabilized linearly implicit approaches (cf. [21, 18]), the method of invariantenergy quadratization (cf. [4, 22]) and the recently developed scalar auxiliary variableapproach (see [12]).In this publication, we are interested in the numerical treatment of (1.9). A first finiteelement scheme was proposed in the Bachelor’s thesis [19] (see also [7] for numericalresults). In this thesis, a straightforward discretization based on continuous, piecewiselinear finite element functions was applied to model (1.9), and the arising nonlinear systemwas solved using Newton’s method. In this publication, we pursuit a different approachand investigate the connection between φ and the chemical potentials.The peculiarity of (1.9) is the coupling between the chemical potential defined inside ofthe domain and the one on the boundary. In the standard Cahn–Hilliard equation (1.1),the chemical potential is merely a definition in terms of φ . This allows us to write (1.1) asa sole, nonlinear, fourth-order equation. In (1.9), however, the chemical potentials µ and µ Γ are coupled via the normal derivative ∇ φ · n . Consequently, the chemical potentials S. METZGER have to be determined by solving a system consisting of (1.9b), (1.9e), and the additionalassumption that φ is continuous on Ω . The latter one translates to the constraint that(1.9a) and (1.9d) yield compatible results. Deducing a suitable expression for µ will be keyingredient for the derivation of an efficient numerical scheme, but also for the numericalanalysis, as the existence of a unique (discrete) µ for any given φ allows us to reusetechniques from the analysis of the standard Cahn–Hilliard equation.The outline of the paper is as follows. In Section 2, we introduce the discrete functionspaces and derive the discrete scheme. In Section 3, we will establish a first a prioriestimate which is discrete counterpart of (1.10), and use this estimate to prove the ex-istence of discrete solutions. The main convergence result, Theorem 4.4, can be foundin Section 4, where we establish improved regularity results and show the convergence ofdiscrete solutions towards weak solutions of (1.9). For the uniqueness results for theseweak solutions, we refer the reader to Section 5 in [7]. We will conclude this manuscript bybriefly discussing the case of Allen–Cahn-type boundary conditions (cf. Remark 4.5). Byshowing that the presented techniques are also applicable for Allen–Cahn-type boundaryconditions, we also cover (1.4) as a special case of (1.6).Notation. Given the spatial domain Ω ⊂ R d with d ∈ { , } and a time interval (0 , T ) ,we denote the space-time cylinder Ω × (0 , T ) by Ω T . By W k,p (Ω) we denote the spaceof k -times weakly differentiable functions with weak derivatives in L p (Ω) . The symbol W k,p (Ω) stands for the closure of C ∞ (Ω) in W k,p (Ω) . For p = 2 , we will denote W k, (Ω) by H k (Ω) and W k, (Ω) by H k (Ω) . The dual space of H (Ω) will be denoted by ( H (Ω)) ′ and the corresponding dual pairing by h ., . i . For a Banach space X and a time interval I ,the symbol L p ( I ; X ) stands for the parabolic space of L p -integrable functions on I withvalues in X . Denoting the boundary of Ω by Γ and the space-time cylinder Γ × (0 , T ) by Γ T , we use a notation similar to the one introduced above for the function spaces definedon Γ and denote the dual pairing between ( H (Γ)) ′ and H (Γ) by h ., . i Γ . In addition, wedefine the space X κ := (cid:26) { v ∈ H (Ω) : γ ( v ) ∈ H (Γ) } if κ > ,H (Ω) if κ = 0 , (1.11)where γ defines the trace operator. The trace operator is uniquely defined and lies in L ( W s,p (Ω) , W s − /p,p (Γ)) for ≤ p ≤ ∞ and s > /p , where s − /p is not an integer.For brevity, we will sometimes (in particular when the considered function is continuous)neglect the trace operator and write v instead of γ ( v ) .2. Derivation of an efficient numerical scheme
We start by introducing the general notation and the discretization techniques used inthe considered scheme. Concerning the discretization with respect to time, we assumethat (T) the time interval I := [0 , T ) is subdivided in intervals I n := [ t n , t n +1 ) with t n +1 = t n + τ n for time increments τ n > and n = 0 , ..., N − with t N = T . For simplicity,we take τ n ≡ τ = TN for n = 0 , ..., N − .The spatial domain Ω ⊂ R d in spatial dimensions d ∈ { , } is assumed to be bounded,convex, and polygonal (or polyhedral, respectively). We introduce partitions T h of Ω and T Γ h of Γ depending on a spatial discretization parameter h > satisfying the followingassumptions: (S1) Let {T h } h> a quasiuniform family (in the sense of [3]) of partitions of Ω intodisjoint, open, nonobtuse simplices K , so that Ω ≡ [ K ∈T h K with max K ∈T h diam ( K ) ≤ h . (S2) Let {T Γ h } h> a quasiuniform family of partitions of Γ into disjoint, open, nonobtusesimplices K Γ , so that ∀ K Γ ∈ T Γ h ∃ ! K ∈ T h such that K Γ = K ∩ Γ , and Γ ≡ [ K Γ ∈T Γ h K Γ with max K Γ ∈T Γ h diam ( K Γ ) ≤ h . (S2) implies that T Γ h is compatible to T h in the sense that all elements in T Γ h are edgesof elements in T h . For the approximation of the phase-field φ and the chemical potential µ we use continuous, piecewise linear finite element functions on T h . This space will bedenoted by U Ω h and is given by functions { χ h,k } k =1 ,..., dim U Ω h forming a dual basis to thevertices { x k } k =1 ,..., dim U Ω h of T h , i.e. χ h,k ( x k ) = δ k,l for k, l = 1 , ..., dim U Ω h . Analogously,we denote the space of continuous, piecewise linear finite element functions on T Γ h by U Γ h . This space is spanned by functions { χ Γ h,k } k =1 ,..., dim U Γ h forming a dual basis to thevertices { x Γ k } k =1 ,..., dim U Γ h of T Γ h , i.e. χ Γ h,k ( x Γ k ) = δ k,l for k, l = 1 , ..., dim U Ω h . Due to thecompatibility condition for T h and T Γ h , we have U Γ h = span { γ ( θ h ) : θ h ∈ U Ω h } . (2.1)Without loss of generality, we may assume that the first dim U Γ h vertices of T h are locatedon Γ , i.e. { x Γ k } k =1 ,..., dim U Γ h = { x k } k =1 ,..., dim U Γ h . We define the nodal interpolation operators I h : C (Ω) → U Ω h and I Γ h : C (Γ) → U Γ h by I h { a } := dim U Ω h X k =1 a ( x k ) χ h,k , and I Γ h { a } := dim U Γ h X k =1 a ( x k ) χ Γ h,k . (2.2)For future reference, we state the following estimate for the interpolation operators. Lemma 2.1.
Let T h and T Γ h satisfy (S1) and (S2). Furthermore, let p ∈ [1 , ∞ ) , ≤ q ≤∞ , and q ∗ = q − q for q < ∞ or q ∗ = 1 for q = ∞ . Then k ( I − I h ) { f h g h }k L p (Ω) ≤ Ch k∇ f h k L pq (Ω) k∇ g h k L pq ∗ (Ω) , (2.3) (cid:13)(cid:13) ( I − I Γ h ) { f h g h } (cid:13)(cid:13) L p (Γ) ≤ Ch k∇ Γ f h k L pq (Γ) k∇ Γ g h k L pq ∗ (Γ) . (2.4) holds true for all f h , g h ∈ U Ω h Proof.
Using the standard error estimates for the nodal interpolation operator (cf. [3])and Hölder’s inequality, we compute on each K ∈ T h : ˆ K | ( I − I h ) { f h g h }| p ≤ Ch p ˆ K | f h g h | pW , ∞ ( K ) ≤ Ch p d X i,j =1 ˆ K | ∂ i f h | p | ∂ j g h | p ≤ Ch p k∇ f h k pL pq ( K ) k∇ g h k pL pq ∗ ( K ) . (2.5) S. METZGER
Similar computations provide the result for I Γ h . (cid:3) Concerning the potentials F and G , we make the following assumptions: (P1) F, G ∈ C ( R ) are bounded from below, i.e. there exists a constant C > suchthat F ( s ) > − C and G ( s ) > − C for all s ∈ R . Furthermore, there exist convex,non-negative functions F + , G + ∈ C ( R ) and concave functions F − , G − ∈ C ( R ) such that F ≡ F + + F − and G ≡ G + + G − . (P2) The convex and concave parts of F and G can be further decomposed into apolynomial part of degree four and an additional part with a globally Lipschitz-continuous first derivative. Moreover, there exists β ≥ such that the concaveparts satisfy G ′− ( s )( s + s ) ≥ G − ( s ) − G − ( s ) + β | s − s | for all s , s ∈ R . In the case κ = 0 , we assume that the above assumption holdstrue for β > . Remark 2.2.
The Assumptions (P1) and (P2) are in particular satisfied by the (penal-ized) polynomial double-well potential W ( φ ) := (1 − φ ) + C pen max { ( | φ | − , } with C pen > , as well as by the typical fluid-solid interfacial energy γ fs with γ fs , + ( φ ) = sin ( π min { max { φ, − } , } ) + π φ and γ fs , − ( φ ) = − π φ . Remark 2.3.
In this publication, we consider only a convex-concave decomposition of thedouble-well potential. Other suitable, energy stable discretization techniques can be foundin [9] . For a comparison of these techniques, we refer to [14] . Defining the backward difference quotient ∂ − τ a n := τ − ( a n − a n − ) , we have all tools athand to propose a fully discrete finite element scheme. Due to the compatibility condition(2.1), we may write the scheme as ˆ Ω I h { φ nh θ h } + τ m ˆ Ω ∇ µ nh · ∇ θ h = ˆ Ω I h (cid:8) φ n − h θ h (cid:9) , (2.6a) ˆ Γ I Γ h { φ nh θ h } + τ m Γ ˆ Γ ∇ Γ µ n Γ ,h · ∇ Γ θ h = ˆ Γ I Γ h (cid:8) φ n − h θ h (cid:9) , (2.6b) ˆ Ω I h { µ nh θ h } + ˆ Γ I Γ h (cid:8) µ n Γ ,h θ h (cid:9) = δσ ˆ Ω ∇ φ nh · ∇ θ h + δ − σ ˆ Ω I h (cid:8)(cid:0) F ′ + ( φ nh ) + F ′− ( φ n − h ) (cid:1) θ h (cid:9) + κδ Γ ˆ Γ ∇ Γ φ nh · ∇ Γ θ h + δ − ˆ Γ I Γ h (cid:8)(cid:0) G ′ + ( φ nh ) + G ′− ( φ n − h ) (cid:1) θ h (cid:9) (2.6c)for all θ h ∈ U Ω h . At this point, the main difference between (2.8) and the establishedschemes developed for the Cahn–Hilliard models without dynamic boundary conditionsbecomes evident. While schemes for the standard Cahn–Hilliard equation (see e.g. [10,9, 6]) allow us to compute µ nh for given φ nh and φ n − h in an easy way, model (1.9) with thedynamic boundary conditions provides only the expression (2.6c) for µ nh and µ n Γ ,h and theconstraint that (2.6a) and (2.6b) have to yield the same results for γ ( φ nh ) . Consequently, the goal for this section will be to derive an equivalent formulation for (2.6) which allowsto compute the chemical potentials in an efficient way.We define the lumped mass matrices M Ω and M Γ via ( M Ω ) ij := ˆ Ω I h { χ hj χ hi } ∀ i, j = 1 , ..., dim U Ω h , (2.7a) ( M Γ ) ij := ˆ Γ I Γ h (cid:8) χ Γ hj χ Γ hi (cid:9) ∀ i, j = 1 , ..., dim U Γ h , (2.7b)and the stiffness matrices L Ω and L Γ via ( L Ω ) ij := ˆ Ω ∇ χ hj · ∇ χ hi ∀ i, j = 1 , ..., dim U Ω h , (2.7c) ( L Γ ) ij := ˆ Γ ∇ Γ χ Γ hj · ∇ Γ χ Γ hi ∀ i, j = 1 , ..., dim U Γ h . (2.7d)Furthermore, we collect the nodal values of φ nh , φ n − h , µ nh , and µ n Γ ,h in the vectors Φ n , Φ n − , P n , and P n Γ . In a slight misuse of notation, we will write F (Φ n ) , when we apply afunction F to all components of Φ n . With this notation, we are able to rewrite (2.6) as M Ω Φ n + τ m L Ω P n = M Ω Φ n − , (2.8a) M Γ (Φ n (cid:12)(cid:12) Γ ) + τ m Γ L Γ P n Γ = M Γ (Φ n − (cid:12)(cid:12) Γ ) , (2.8b) M Ω P n + ( M Γ P n Γ ) (cid:12)(cid:12) Ω = δσ L Ω Φ n + δ − σ M Ω F ′ + (Φ n ) + δ − σ M Ω F ′− (Φ n − )+ (cid:0) δ Γ κ L Γ (Φ n (cid:12)(cid:12) Γ ) + δ − M Γ G ′ + (Φ n (cid:12)(cid:12) Γ ) + δ − M Γ G ′− (Φ n − (cid:12)(cid:12) Γ ) (cid:1)(cid:12)(cid:12) Ω . (2.8c)Here, we used the extension operator . (cid:12)(cid:12) Ω : R dim U Γ h → R dim U Ω h defined via R dim U Γ h ∋ A ( A ) ∈ R dim U Ω h (2.9)and the restriction operator . (cid:12)(cid:12) Γ : R dim U Ω h → R dim U Γ h , which restricts a vector its first dim U Γ h entries.In order to derive a scheme allowing to solve (2.8) efficiently, we define restriction operatorsfor matrices. In particular, we will split a matrix A ∈ R dim U Ω h × dim U Ω h into submatrices A (cid:12)(cid:12) Γ × Γ ∈ R dim U Γ h × dim U Γ h , A (cid:12)(cid:12) Γ × ◦ Ω ∈ R dim U Γ h × (dim U Ω h − dim U Γ h ) , A (cid:12)(cid:12) ◦ Ω × Γ ∈ R (dim U Ω h − dim U Γ h ) × dim U Γ h , A (cid:12)(cid:12) ◦ Ω × ◦ Ω ∈ R (dim U Ω h − dim U Γ h ) × (dim U Ω h − dim U Γ h ) , A (cid:12)(cid:12) Γ × Ω ∈ R dim U Γ h × dim U Ω h , A (cid:12)(cid:12) ◦ Ω × Ω ∈ R (dim U Ω h − dim U Γ h ) × dim U Ω h , A (cid:12)(cid:12) Ω × Γ ∈ R dim U Ω h × dim U Γ h , and A (cid:12)(cid:12) Ω × ◦ Ω ∈ R dim U Ω h × (dim U Ω h − dim U Γ h ) , (2.10)such that A = A (cid:12)(cid:12) Γ × Γ A (cid:12)(cid:12) Γ × ◦ Ω A (cid:12)(cid:12) ◦ Ω × Γ A (cid:12)(cid:12) ◦ Ω × ◦ Ω ! = A (cid:12)(cid:12) Γ × Ω A (cid:12)(cid:12) ◦ Ω × Ω ! = (cid:16) A (cid:12)(cid:12) Ω × Γ A (cid:12)(cid:12) Ω × ◦ Ω (cid:17) . (2.11) S. METZGER
Hence, the chemical potentials are given as solutions of the (dim U Ω h + dim U Γ h ) × (dim U Ω h + dim U Γ h ) -system M Ω (cid:12)(cid:12) Γ × Γ Γ Ω (cid:12)(cid:12) ◦ Ω × ◦ Ω m M − (cid:12)(cid:12) Γ × Ω L Ω (cid:12)(cid:12) Ω × Γ m M − (cid:12)(cid:12) Γ × Ω L Ω (cid:12)(cid:12) Ω × ◦ Ω − m Γ M − L Γ P n (cid:12)(cid:12) Γ P n (cid:12)(cid:12) ◦ Ω P n Γ = R Γ (Φ n )R ◦ Ω (Φ n )0 (2.12)with R Γ (Φ n ) := δσ L Ω (cid:12)(cid:12) Γ × Ω Φ n + δ − σ M Ω (cid:12)(cid:12) Γ × Ω F ′ + (Φ n ) + δ − σ M Ω (cid:12)(cid:12) Γ × Ω F ′− (Φ n − )+ κδ Γ L Γ Φ n (cid:12)(cid:12) Γ + δ − M Γ G ′ + (Φ n (cid:12)(cid:12) Γ ) + δ − M Γ G ′− (Φ n − (cid:12)(cid:12) Γ ) and R ◦ Ω (Φ n ) := δσ L Ω (cid:12)(cid:12) ◦ Ω × Ω Φ n + δ − σ M Ω (cid:12)(cid:12) ◦ Ω × Ω F ′ + (Φ n ) + δ − σ M Ω (cid:12)(cid:12) ◦ Ω × Ω F ′− (Φ n − ) . Here,the first two lines are a consequence of (2.8c) and the last line guarantees that (2.8a) and(2.8b) provide the same result for Φ n (cid:12)(cid:12) Γ .As the (2.8) is nonlinear in Φ n , computing a possible solution requires the application ofan iterative scheme (e.g. Newton’s method) and therefore solving (2.12) multiple timesper time step. Hence, solving a (dim U Ω h + dim U Γ h ) × (dim U Ω h + dim U Γ h ) -system eachtime is not desirable and we have to continue reducing the complexity of the system.From the second line in (2.12), we immediately get P n (cid:12)(cid:12) ◦ Ω = M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) , (2.13)while the first line provides P n Γ = − M − M Ω (cid:12)(cid:12) Γ × Γ P n (cid:12)(cid:12) Γ + M − R Γ (Φ n ) . (2.14)Plugging this into the last line, we obtain m M − (cid:12)(cid:12) Γ × Ω L Ω (cid:12)(cid:12) Ω × Γ P n (cid:12)(cid:12) Γ = − m M − (cid:12)(cid:12) Γ × Ω L Ω (cid:12)(cid:12) Ω × ◦ Ω P n (cid:12)(cid:12) ◦ Ω + m Γ M − L Γ P n Γ = − m M − (cid:12)(cid:12) Γ × Ω L Ω (cid:12)(cid:12) Ω × ◦ Ω M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) − m Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ P n (cid:12)(cid:12) Γ + m Γ M − L Γ M − R Γ (Φ n ) , (2.15)and therefore (cid:16) m M − (cid:12)(cid:12) Γ × Ω L Ω (cid:12)(cid:12) Ω × Γ + m Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) P n (cid:12)(cid:12) Γ = − m M − (cid:12)(cid:12) Γ × Ω L Ω (cid:12)(cid:12) Ω × ◦ Ω M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) + m Γ M − L Γ M − R Γ (Φ n ) . (2.16)As M − is a diagonal matrix, M − (cid:12)(cid:12) Γ × Ω L Ω (cid:12)(cid:12) Ω × Γ = M − (cid:12)(cid:12) Γ × Γ L Ω (cid:12)(cid:12) Γ × Γ holds true. This allowsus to multiply (2.16) by M Ω (cid:12)(cid:12) Γ × Γ to obtain (cid:16) m L Ω (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) P n (cid:12)(cid:12) Γ = − m L Ω (cid:12)(cid:12) Γ × ◦ Ω M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − R Γ (Φ n ) . (2.17)In order to show that (2.17) is well defined, we need to prove that the matrix on theleft-hand side is indeed invertible. Lemma 2.4.
The matrix (cid:16) m L Ω (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) , that is definedvia (2.7) and (2.10) , is symmetric, positive definite.Proof. It is obvious that m L Ω (cid:12)(cid:12) Γ × Γ and m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ are symmetric,positive semi-definite matrices. Therefore, it will be sufficient to show that A T L Ω (cid:12)(cid:12) Γ × Γ A > for all = A ∈ R dim U Γ h × U Γ h to complete the proof. This is equivalent to showing ˜ A T L Ω ˜ A > with ˜ A = A (cid:12)(cid:12) Ω = ( A ) for all = A ∈ R dim U Γ h × U Γ h . (2.18)From (2.7), we have that L Ω is symmetric, positive semi-definite with only constantvectors corresponding to the zero eigenvalue. As the restrictions in (2.18) do not allowfor constant vectors, the proof is complete. (cid:3) Combining (2.17) with (2.13), we obtain an expression for the chemical potential whichrequires to solve only a dim U Γ h by dim U Γ h linear system with a sparse, symmetric, positivedefinite matrix. Having an expression for the chemical potential, we propose the followingnonlinear equation for Φ n . Φ n + τ m M − L Ω (cid:16) m L Ω (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) − ! · − m L Ω (cid:12)(cid:12) Γ × ◦ Ω M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − R Γ (Φ n ) M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) ! = Φ n − . (2.19)This nonlinear equation can be tackled using e.g. Newton’s method. Thereby, a linearsystem with the matrix (cid:16) m L Ω (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) has to be solvedrepeatedly. As this matrix is symmetric, positive definite, it is predestined for the appli-cation of a conjugate gradient method or a Cholesky decomposition.3. Stability and existence of discrete solutions
In this section, we analyze the discrete scheme (2.19) proposed in the previous section.Although (2.19) is entirely written in terms of the unknown quantity Φ n , we will continueusing P n and P n Γ , which are defined in (2.17), (2.13), and (2.14), to simplify the notation.For the ease of representation, we will set σ = δ = δ Γ = 1 for the remainder of thispublication. As a first step, we shall verify that (2.19) indeed satisfies the compatibilityconstraint m M − (cid:12)(cid:12) Γ × Γ L Ω (cid:12)(cid:12) Γ × Ω P n = m Γ M − L Γ P n Γ . This auxiliary result allows us derivean a priori stability result for (2.19) which serves as the corner stone for proving theexistence of discrete solutions. Lemma 3.1.
Let P n and P n Γ be defined via (2.17) , (2.13) , and (2.14) . Then the identity m Γ M − L Γ P n Γ − m M − (cid:12)(cid:12) Γ × Γ L Ω (cid:12)(cid:12) Γ × Ω P n = 0 holds true. Proof.
Using (2.13), we compute m Γ M − L Γ P n Γ − m M − (cid:12)(cid:12) Γ × Γ L Ω (cid:12)(cid:12) Γ × Ω P n = − m Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ P n (cid:12)(cid:12) Γ + m Γ M − L Γ M − R Γ (Φ n ) − m M − (cid:12)(cid:12) Γ × Γ L Ω (cid:12)(cid:12) Γ × Ω P n = m Γ M − L Γ M − R Γ (Φ n ) − M − (cid:12)(cid:12) Γ × Γ (cid:16) m L Γ (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) P n (cid:12)(cid:12) Γ − m M − (cid:12)(cid:12) Γ × Γ L Ω (cid:12)(cid:12) Γ × ◦ Ω P n (cid:12)(cid:12) ◦ Ω =: I + II + III . (3.1)Recalling (2.17) and (2.14), we obtain II = m M − (cid:12)(cid:12) Γ × Γ L Ω (cid:12)(cid:12) Γ × ◦ Ω M − (cid:12)(cid:12) ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) − m Γ M − L Γ M − R Γ (Φ n )= m M − (cid:12)(cid:12) Γ × Γ L Ω (cid:12)(cid:12) Γ × ◦ Ω P n (cid:12)(cid:12) ◦ Ω − m Γ M − L Γ M − R Γ (Φ n ) = − III − I , (3.2)which completes the proof. (cid:3)
This result allows us to show that the phasefield parameter is conserved in Ω and on Γ . Testing (2.19) by T M Ω and by T Γ M Γ (cid:12)(cid:12) Ω proves the following corollary. Corollary 3.2.
Let Φ n be a discrete solution of (2.19) . Then T M Ω Φ n = T M Ω Φ n − T Γ M Γ Φ n (cid:12)(cid:12) Γ = T Γ M Γ Φ n − (cid:12)(cid:12) Γ with := (1 , ..., T ∈ R dim U Ω h and Γ := (cid:12)(cid:12) Γ . Using the above auxiliary results, we are now able to state a first stability result whichis a discrete version of the energy equality (1.10).
Lemma 3.3.
Let the assumptions (T), (S1), (S2), (P1), and (P2) hold true and let Φ n − ∈ R dim U Ω h be given. Then a solution Φ n ∈ R dim U Ω h to (2.19) , if it exists, satifies Φ nT L Ω Φ n + (Φ n − Φ n − ) T L Ω (Φ n − Φ n − ) + T M Ω F (Φ n ) + κ Φ n (cid:12)(cid:12) T Γ L Γ Φ n (cid:12)(cid:12) Γ + κ (Φ n − Φ n − ) (cid:12)(cid:12) T Γ L Γ (Φ n − Φ n − ) (cid:12)(cid:12) Γ + T Γ M Γ G (Φ n (cid:12)(cid:12) Γ )+ β (Φ n − Φ n − ) (cid:12)(cid:12) T Γ M Γ (Φ n − Φ n − ) (cid:12)(cid:12) Γ + τ mP nT L Ω P n + τ m Γ P n Γ T L Γ P n Γ ≤ Φ n − T L Ω Φ n − + T M Ω F (Φ n − ) + κ Φ n − (cid:12)(cid:12) T Γ L Γ Φ n − (cid:12)(cid:12) Γ + T Γ M Γ G (Φ n − (cid:12)(cid:12) Γ ) , with := (1 , ..., T ∈ R dim U Ω h , Γ := (cid:12)(cid:12) Γ , and P n and P n Γ defined in (2.13) , (2.17) , and (2.14) .Proof. We test (2.19) by (cid:0) M Ω P n + (cid:0) M Γ P n Γ (cid:1)(cid:1) and use Lemma 3.1 to obtain (cid:0) Φ n − Φ n − (cid:1) T M Ω P n + (cid:0) Φ n − Φ n − (cid:1)(cid:12)(cid:12) T Γ M Γ P n Γ + τ m ( P n ) T L Ω P n + τ m Γ ( P n Γ ) T L Γ P n Γ =: I + II + III + IV . (3.3)As
III and IV provide the dissipative parts of the desired estimate, we have show to that I and II yield the time difference of the energy. Recalling (2.14), we compute II = − (cid:0) Φ n − Φ n − (cid:1)(cid:12)(cid:12) T Γ M Ω (cid:12)(cid:12) Γ × Γ P n (cid:12)(cid:12) Γ + (cid:0) Φ n − Φ n − (cid:1)(cid:12)(cid:12) T Γ R Γ (Φ n ) . (3.4) Consequently, we obtain from (2.17) I + II = (Φ n − Φ n − ) (cid:12)(cid:12) T ◦ Ω R ◦ Ω (Φ n ) + (Φ n − Φ n − ) (cid:12)(cid:12) T Γ R Γ (Φ n )= (cid:0) Φ n − Φ n − (cid:1) T L Ω Φ n + (Φ n − Φ n − ) T M Ω ( F ′ + (Φ n ) + F ′− (Φ n − ))+ κ (Φ n − Φ n − ) (cid:12)(cid:12) T Γ L Γ Φ n (cid:12)(cid:12) Γ + (Φ n − Φ n − ) (cid:12)(cid:12) T Γ M Γ ( G ′ + (Φ n (cid:12)(cid:12) Γ ) + G ′− (Φ n − (cid:12)(cid:12) Γ )) . (3.5)As M Ω and M Γ are diagonal matrices, we may combine ( F ′ + (Φ n ) + F ′− (Φ n − )) and (Φ n − Φ n − ) ,and ( G ′ + (Φ n (cid:12)(cid:12) Γ ) + G ′− (Φ n − (cid:12)(cid:12) Γ )) and (Φ n − Φ n − ) (cid:12)(cid:12) Γ componentwise. In combination with s ( s − s ) = s + ( s − s ) − s , this provides the result. (cid:3) Using the a priori estimate from Lemma 3.3, we are able to prove the existence ofdiscrete solutions.
Lemma 3.4.
Let the assumptions (T), (S1), (S2), and (P1) hold true and let Φ n − ∈ R dim U Ω h be given. Then, there exists at least one vector Φ n ∈ R dim U Ω h solving (2.19) .Proof. We will prove the existence of discrete solutions by contradiction. Let ||| . ||| denotethe discrete L -norm which is derived from the inner product ( A, B ) := A T M Ω B . Ac-cording to Corollary 3.2, the mean-value of the phase-field is conserved in Ω . This allowsus to assume w.l.o.g. that T M Ω Φ n = T M Ω Φ n − = 0 . Therefore, p Φ nT L Ω Φ n is also anorm of Φ n . Under the assumption that (2.19) has no solution in B R := { A ∈ R dim U Ω h : T M Ω A = 0 and ||| A ||| ≤ R } (3.6)for any R > , the function H defined via Φ − Φ n − + τ m M − L Ω (cid:16) m L Ω (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) − ! · − m L Ω (cid:12)(cid:12) Γ × ◦ Ω M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ) + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − R Γ (Φ) M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ) ! =: H (Φ) (3.7)has no root and is continuous on B R . This allows us to define a function A : B R → ∂B R ⊂ B R as A (Φ) := − R H (Φ) |||H (Φ) ||| . (3.8)As A is continuous and maps a closed set onto itself, Brouwer’s fixed point theoremprovides the existence of at least one fixed point Φ ∗ . In the following, we will show < (Φ ∗ , Ψ) < (3.9)for a suitable Ψ ∈ B R and R large enough. This contradiction shows that the initialassumption of (2.19) not having solutions in B R is wrong. To prove the contradiction (3.9), we choose Ψ = ˜Ψ + ˜Ψ − T M Ω ( ˜Ψ + ˜Ψ )( T M Ω ) − with ˜Ψ := (cid:16) m L Ω (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) − ! · − m L Ω (cid:12)(cid:12) Γ × ◦ Ω M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ ∗ ) + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − R Γ (Φ ∗ ) M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ ∗ ) ! (3.10)and ˜Ψ := M − (cid:12)(cid:12) Γ × Γ M Γ (cid:16) − M − M Ω (cid:12)(cid:12) Γ × Γ ˜Ψ (cid:12)(cid:12)(cid:12) Γ + M − R Γ (Φ ∗ ) (cid:17) ! , (3.11)i.e. the test function is the sum of the chemical potentials deprived of their mean values.The computations from the proof of Lemma 3.3 provide ( H (Φ ∗ ) , Ψ) ≥ Φ ∗ T L Ω Φ ∗ − C , (3.12)where the constant C depends on Φ n − and the lower bound from (P1), but not on thefixed point Φ ∗ or R . Since all norms on finite dimensional spaces are equivalent, thereexists c > such that Φ ∗ T L Ω Φ ∗ ≥ c Φ ∗ T M Ω Φ ∗ and we obtain ( H (Φ ∗ ) , Ψ) ≥ c ||| Φ ∗ ||| − C = cR − C > (3.13)for R large enough. This provides the second inequality in (3.9). In order to establish thefirst inequality we again use the computations from the proof of Lemma 3.3 to show (Φ ∗ , Ψ) =Φ ∗ T L Ω Φ ∗ + Φ ∗ T M Ω (cid:0) F ′ + (Φ ∗ ) + F ′− ( ) (cid:1) + Φ ∗ T M Ω (cid:0) F ′− (Φ n − ) − F ′− ( ) (cid:1) + κ Φ ∗ (cid:12)(cid:12) T Γ L Γ Φ ∗ (cid:12)(cid:12) Γ + Φ ∗ (cid:12)(cid:12) T Γ M Γ (cid:0) G ′ + (Φ ∗ (cid:12)(cid:12) Γ ) + G ′− ( ) (cid:1) + Φ ∗ (cid:12)(cid:12) T Γ M Γ (cid:0) G ′− (Φ n − (cid:12)(cid:12) Γ ) − G ′− ( ) (cid:1) ≥ c ||| Φ ∗ ||| + T M Ω ( F (Φ ∗ ) − F ( )) − ε ||| Φ ∗ ||| − C ε (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) F ′− (Φ n − ) − F ′− ( ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) + κ Φ ∗ (cid:12)(cid:12) T Γ L Γ Φ ∗ (cid:12)(cid:12) Γ + T Γ M Γ (cid:0) G (Φ ∗ (cid:12)(cid:12) Γ ) − G ( ) (cid:1) − ˜ ε Φ ∗ (cid:12)(cid:12) T Γ M Γ Φ ∗ (cid:12)(cid:12) Γ − C ˜ ε (cid:0) G ′− (Φ n − (cid:12)(cid:12) Γ ) − G ′− ( ) (cid:1) T M Γ (cid:0) G ′− (Φ n − (cid:12)(cid:12) Γ ) − G ′− ( ) (cid:1) (3.14)with < ε, ˜ ε << . For every fixed h , there is a constant C h > such that Φ ∗ (cid:12)(cid:12) T Γ M Γ Φ ∗ (cid:12)(cid:12) Γ ≤ C h Φ ∗ T M Ω Φ ∗ . Hence, we have (Φ ∗ , Ψ) ≥ ( c − ε − C h ˜ ε ) ||| Φ ∗ ||| − C ε, ˜ ε = ( c − ε − C h ˜ ε ) R − C ε, ˜ ε with C ε, ˜ ε > independent of Φ ∗ and R . Choosing ε and ˜ ε small enough provides ( c − ε − C h ˜ ε ) > . Hence, we obtain the first inequality in (3.9) for R large enough,which completes the proof. (cid:3) Remark 3.5.
The existence result in Lemma 3.3 implies no constraints on the timeincrement τ . Therefore, we have the existence of discrete solutions for arbitrary timeincrements. Convergence of the discrete scheme
In this section, we show that the discrete solutions established in the last section con-verge towards suitable weak solutions of (1.9). This requires some assumptions on theinitial data. In particular, we will assume that (I) the initial data φ ∈ X κ and its projection φ h onto U Ω h satisfies ˆ Ω (cid:12)(cid:12) ∇ φ h (cid:12)(cid:12) + ˆ Ω I h (cid:8) F ( φ h ) (cid:9) + κ ˆ Γ (cid:12)(cid:12) ∇ Γ φ h (cid:12)(cid:12) + ˆ Γ I Γ h (cid:8) G ( φ h ) (cid:9) ≤ C with some C > independent of h and τ .Furthermore, the regularity results provided in this section require additional assumptionson h and τ . In particular, we will need (C) that h τ ց for ( h, τ ) ց when κ > and that h τ ց for ( h, τ ) ց when κ = 0 .Assumption (I) allows us to state our first regularity result. Corollary 4.1.
Let the assumptions (T), (S1), (S2), (P1), (P2), and (I) hold true andlet h > be small enough. Then a solution ( φ nh , µ nh , µ n Γ ,h ) n =1 ,...,N to (2.6) satisfies max n =0 ,...,N k φ nh k H (Ω) + max n =0 ,...,N ˆ Ω I h { F ( φ nh ) } + κ max n =0 ,...,N k φ nh k H (Γ) + max n =0 ,...,N ˆ Γ I Γ h { G ( φ nh ) } + N X n =1 ˆ Ω (cid:12)(cid:12) ∇ φ nh − ∇ φ n − h (cid:12)(cid:12) + κ N X n =1 ˆ Γ (cid:12)(cid:12) ∇ Γ φ nh − ∇ Γ φ n − h (cid:12)(cid:12) + β N X n =1 ˆ Γ (cid:12)(cid:12) φ nh − φ n − h (cid:12)(cid:12) + τ m N X n =1 k µ nh k H (Ω) + τ m Γ N X n =1 (cid:13)(cid:13) µ n Γ ,h (cid:13)(cid:13) H (Γ) ≤ C , with a constant
C > independent of h and τ .Proof. After summing the result of Lemma 3.3 over all time steps and recalling Corollary3.2 and (I), it remains to show that we have indeed control over the complete H normof µ nh and µ n Γ ,h . To establish this result, we will follow the lines of [7].Testing (2.6c) by I h { η } with η ∈ C ∞ (Ω; [0 , , which is not identically zero, we obtain ˆ Ω I h { µ nh η } = ˆ Ω ∇ φ nh · ∇I h { η } + ˆ Ω I h (cid:8) ( F ′ + ( φ nh ) + F ′− ( φ n − h )) η (cid:9) . (4.1)From (P2) and (I), we obtain (cid:12)(cid:12)(cid:12)(cid:12) ˆ Ω I h (cid:8) ( F ′ + ( φ nh ) + F ′− ( φ n − h )) η (cid:9)(cid:12)(cid:12)(cid:12)(cid:12) ≤ C k φ nh k L (Ω) + C k φ nh k L (Ω) + C (cid:13)(cid:13) φ n − h (cid:13)(cid:13) L (Ω) + C (cid:13)(cid:13) φ n − h (cid:13)(cid:13) L (Ω) + C ≤ C (4.2)Hence, there exists a constant ˜ C ( η ) independent of h and τ such that (cid:12)(cid:12) ´ Ω I h { µ nh η } (cid:12)(cid:12) ≤ ˜ C ( η ) .We now define M η := (cid:26) v ∈ H (Ω) : ˆ Ω I h { vη } ≤ ˜ C ( η ) (cid:27) . (4.3) From standard error estimates for the interpolation operator I h (cf. [3]), we derive theexistence of c ( η ) > such that ´ Ω I h { η } ≥ c ( η ) for h small enough.Therefore, we mayuse the generalized Poincaré inequality (cf. [1]), which we cite in the appendix as LemmaA.1, with u ≡ and C := ˜ C ( η ) /c ( η ) to obtain k µ nh k L (Ω) ≤ C (1 + k∇ µ nh k L (Ω) ) for all n ∈ { , ..., N } . (4.4)To obtain the L -bound for µ n Γ ,h , we test (2.6c) by θ h ≡ and obtain (cid:12)(cid:12)(cid:12)(cid:12) ˆ Γ µ n Γ ,h (cid:12)(cid:12)(cid:12)(cid:12) ≤ (cid:12)(cid:12)(cid:12)(cid:12) ˆ Ω µ nh (cid:12)(cid:12)(cid:12)(cid:12) + (cid:12)(cid:12)(cid:12)(cid:12) ˆ Ω I h (cid:8) F ′ + ( φ nh ) + F ′− ( φ n − h ) (cid:9)(cid:12)(cid:12)(cid:12)(cid:12) + (cid:12)(cid:12)(cid:12)(cid:12) ˆ Γ I Γ h (cid:8) G ′ + ( φ nh ) + G ′− ( φ n − h ) (cid:9)(cid:12)(cid:12)(cid:12)(cid:12) . Considerations similar to (4.2) show that the last term on the right-hand side is alsobounded by a constant independent of h and τ . Therefore, we may use Poincaré’s in-equality to complete the proof. (cid:3) In a second step, we derive uniform bounds for the time difference quotient of thephase-field parameter on Ω and Γ . Lemma 4.2.
Let the assumptions (T), (S1), (S2), (P1), (P2), (I), and (C) hold true.Furthermore, let h > be small enough such that Corollary 4.1 holds true. Then a solution ( φ nh ) n =1 ,...,N to (2.6) satisfies τ N X n =1 (cid:13)(cid:13) ∂ − τ φ nh (cid:13)(cid:13) H (Ω)) ′ ≤ C , and τ N X n =1 (cid:13)(cid:13) ∂ − τ φ nh (cid:13)(cid:13) H (Γ)) ′ ≤ C , (4.5) with
C > independent of h and τ .Proof. We take θ ∈ H (Ω) and test (2.6a) by θ h := P U Ω h θ , where P U Ω h is the orthogonal L -projection onto U Ω h . We decompose the first term in (2.6a) into ˆ Ω I h (cid:8) ∂ − τ φ nh θ h (cid:9) = ˆ Ω ∂ − τ φ nh θ − ˆ Ω ( I − I h ) (cid:8) ∂ − τ φ nh θ h (cid:9) . (4.6)The first term will be used to obtain a norm on the dual space of H (Ω) . The secondterm can be controlled via Lemma 2.1 and the H -stability of P U Ω h (cf. [2]). Using theseconsiderations and Hölder’s inequality, we obtain (cid:12)(cid:12)(cid:12)(cid:12) ˆ Ω ∂ − τ φ nh θ (cid:12)(cid:12)(cid:12)(cid:12) ≤ C h τ (cid:13)(cid:13) ∇ φ nh − ∇ φ n − h (cid:13)(cid:13) L (Ω) k θ k H (Ω) + k∇ µ nh k L (Ω) k θ k H (Ω) . (4.7)Dividing by k θ k H (Ω) , taking the second power on both sides, multiplying by τ , andsumming over all time steps provides τ N X n =1 (cid:13)(cid:13) ∂ − τ φ nh (cid:13)(cid:13) H (Ω)) ′ ≤ C h τ N X n =1 (cid:13)(cid:13) ∇ φ nh − ∇ φ n − h (cid:13)(cid:13) L (Ω) + Cτ N X n =1 k∇ µ nh k L (Ω) . (4.8)Applying the already established regularity results and (C) completes the proof of the leftinequality in (4.5). For the case κ > , the right inequality in (4.5) can be establishedusing similar computations. In the case κ = 0 , we combine Lemma 2.1 with an inverseestimate and obtain τ N X n =1 (cid:13)(cid:13) ∂ − τ φ nh (cid:13)(cid:13) H (Γ)) ′ ≤ C h τ N X n =1 (cid:13)(cid:13) φ nh − φ n − h (cid:13)(cid:13) L (Γ) + Cτ N X n =1 (cid:13)(cid:13) ∇ µ n Γ ,h (cid:13)(cid:13) L (Γ) . (4.9) Again, the already established regularity results and (C) complete the proof. (cid:3)
In order to pass to the limit ( h, τ ) ց , we define time-interpolants of time-discretefunctions a n , n = 0 , ..., N , and introduce some time-index-free notation as follows. a τ ( ., t ) := t − t n − τ a n ( . ) + t n − tτ a n − ( . ) t ∈ [ t n − , t n ] , n ≥ , (4.10a) a τ, + ( ., t ) := a n ( . ) , a τ, − ( ., t ) := a n − ( . ) t ∈ ( t n − , t n ] , n ≥ . (4.10b)We want to point out that the time derivative of a τ coincides with the difference quotient,i.e. ∂ t a τ = ∂ t (cid:16) t − t n − τ a n + t n − tτ a n − (cid:17) = τ a n − τ a n − = ∂ − τ a n . (4.11)If a statement is valid for a τ , a τ, + , and a τ, − , we use the abbreviation a τ, ( ± ) . With thisnotation, system (2.6) reads as follows. ˆ Ω T I h { ∂ t φ τh θ h } + m ˆ Ω T ∇ µ τ, + h · ∇ θ h =0 , (4.12a) ˆ Γ T I Γ h { ∂ t φ τh θ h } + m Γ ˆ Γ T ∇ Γ µ τ, +Γ ,h · ∇ Γ θ h =0 , (4.12b) ˆ Ω T I h (cid:8) µ τ, + h θ h (cid:9) + ˆ Γ T I Γ h (cid:8) µ τ, +Γ ,h θ h (cid:9) = ˆ Ω T ∇ φ τ, + h · ∇ θ h + ˆ Ω T I h (cid:8)(cid:0) F ′ + ( φ τ, + h ) + F ′− ( φ τ, − h ) (cid:1) θ h (cid:9) + κ ˆ Γ T ∇ Γ φ τ, + h · ∇ Γ θ h + ˆ Γ T I Γ h (cid:8)(cid:0) G ′ + ( φ τ, + h ) + G ′− ( φ τ, − h ) (cid:1) θ h (cid:9) (4.12c)for all θ h ∈ L (0 , T ; U Ω h ) . Similarly, we can rewrite the regularity results obtained inCorollary 4.1 and Lemma 4.2 as (cid:13)(cid:13)(cid:13) φ τ, ( ± ) h (cid:13)(cid:13)(cid:13) L ∞ (0 ,T ; H (Ω)) + κ (cid:13)(cid:13)(cid:13) φ τ, ( ± ) h (cid:13)(cid:13)(cid:13) L ∞ (0 ,T ; H (Ω)) + τ − (cid:13)(cid:13) ∇ φ τ, + h − ∇ φ τ, − h (cid:13)(cid:13) L (0 ,T ; L (Ω)) + κτ − (cid:13)(cid:13) ∇ Γ φ τ, + h − ∇ Γ φ τ, − h (cid:13)(cid:13) L (0 ,T ; L (Γ)) + βτ − (cid:13)(cid:13) φ τ, + h − φ τ, − h (cid:13)(cid:13) L (0 ,T ; L (Γ)) + (cid:13)(cid:13) µ τ, + h (cid:13)(cid:13) L (0 ,T ; H (Ω)) + (cid:13)(cid:13) µ τ, +Γ ,h (cid:13)(cid:13) L (0 ,T ; H (Γ)) ≤ C , (4.13a)as well as k ∂ t φ τh k L (0 ,T ;( H (Ω)) ′ ) ≤ C and k ∂ t φ τh k L (0 ,T ;( H (Γ)) ′ ) ≤ C . (4.13b)These regularity results can be used to identify converging subsequences.
Lemma 4.3.
Let the assumptions (T), (S1), (S2), (P1), (P2), (I), and (C) hold true.Furthermore, let ( φ τ, ( ± ) h , µ τ, + h , µ τ, +Γ ,h ) be a solution to (4.12) . Then there exists a subse-quence (again denoted by ( φ τ, ( ± ) h , µ τ, + h , µ τ, +Γ ,h ) ) and functions φ ∈ L ∞ (0 , T ; H (Ω)) ∩ H (0 , T ; ( H (Ω)) ′ ) , (4.14a) ψ ∈ (cid:26) L ∞ (0 , T ; H (Γ)) ∩ H (0 , T ; ( H (Γ)) ′ ) if κ > ,L ∞ (0 , T ; H / (Γ)) ∩ H (0 , T ; ( H (Γ)) ′ ) if κ = 0 , (4.14b) µ ∈ L (0 , T ; H (Ω)) , (4.14c) µ Γ ∈ L (0 , T ; H (Γ)) (4.14d) such that γ ( φ ) = ψ almost everywhere on Γ T and for ( h, τ ) ց φ τ, ( ± ) h ∗ ⇀ φ in L ∞ (0 , T ; H (Ω)) , (4.15a) ∂ t φ τh ⇀ ∂ t φ in L (0 , T ; ( H (Ω)) ′ ) , (4.15b) φ τ, ( ± ) h → φ in L p (0 , T, L s (Ω)) , (4.15c) γ ( φ τ, ( ± ) h ) ∗ ⇀ ψ in (cid:26) L ∞ (0 , T ; H (Γ)) if κ > ,L ∞ (0 , T ; H / (Γ)) if κ = 0 , (4.15d) ∂ t γ ( φ τh ) ⇀ ∂ t ψ in L (0 , T ; ( H (Γ)) ′ ) , (4.15e) γ ( φ τ, ( ± ) h ) → ψ in (cid:26) L p (0 , T ; L ∞ (Γ)) if κ > ,L p (0 , T ; L ˜ s (Γ)) if κ = 0 , (4.15f) µ τ, + h ⇀ µ in L (0 , T ; H (Ω)) , (4.15g) µ τ, +Γ ,h ⇀ µ Γ in L (0 , T ; H (Γ)) (4.15h) for all p < ∞ , s ∈ [1 , dd − ) , and ˜ s ∈ [1 , d − d − ) .Proof. The weak and weak ∗ convergence expressed in (4.15a), (4.15b), (4.15g), and (4.15h)follows directly from the bounds in (4.13a) and (4.13b). The strong convergence in(4.15c) then follows from the bounds for φ τ, ( ± ) h in L ∞ (0 , T ; H (Ω)) , the bounds on ∂ t φ τh in L (0 , T ; ( H (Ω)) ′ ) , the Aubin–Lions theorem, and the fact that φ τ, + h , φ τ, − h , and φ τh convergetowards the same limit function due to the bound on τ − (cid:13)(cid:13) ∇ φ τ, + h − ∇ φ τ, − h (cid:13)(cid:13) L (0 ,T ; L (Ω)) .Similar arguments provide (4.15d)-(4.15f) in the case κ > . In the case κ = 0 , weuse the uniform bound on k φ nh k H (Ω) to deduce a uniform bound for k γ ( φ nh ) k H (1 / (Γ) . As H / (Γ) is compactly embedded in L ˜ s (Γ) for ˜ s ∈ [1 , d − d − ) (cf. [15]), we verify (4.15d)-(4.15f) for κ = 0 . It remains to show that ψ can be identified with γ ( φ ) . We choose θ ∈ L (0 , T ; ( C ∞ (Ω)) d ) and compute ˆ Ω T φ div θ ← ˆ Ω T φ τ, ( ± ) h div θ = − ˆ Ω T ∇ φ τ, ( ± ) h · θ + ˆ Γ T γ ( φ τ, ( ± ) h ) θ · n → − ˆ Ω T ∇ φ · θ + ˆ Γ T ψ θ · n = ˆ Ω T φ div θ − ˆ Γ T γ ( φ ) θ · n + ˆ Γ T ψ θ · n . (4.16) (cid:3) Theorem 4.4.
Let d ∈ { , } and let the assumptions (T), (S1), (S2), (P1), (P2), (I),and (C) hold true. Then a tuple ( φ, µ, µ Γ ) can be obtained from discrete solutions to (2.6) by passing to the limit ( h, τ ) ց that solves (1.9) in the following weak sense: ˆ T h ∂ t φ, θ i + m ˆ Ω T ∇ µ · ∇ θ =0 ∀ θ ∈ L (0 , T ; H (Ω)) , (4.17a) ˆ T h ∂ t γ ( φ ) , θ i Γ + m Γ ˆ Γ T ∇ Γ µ Γ · ∇ Γ θ =0 ∀ θ ∈ L (0 , T ; H (Γ)) , (4.17b) ˆ Ω T µθ + ˆ Γ T µ Γ θ = ˆ Ω T ∇ φ · ∇ θ + ˆ Ω T F ′ ( φ ) θ + κ ˆ Γ T ∇ Γ γ ( φ ) · ∇ Γ θ + ˆ Γ T G ′ ( γ ( φ )) θ ∀ θ ∈ L (0 , T ; X κ ) . (4.17c) Proof.
We start by passing to the limit in (4.12a). Choosing θ h := I h { θ } for θ ∈ L (0 , T ; C ∞ (Ω)) , we have θ h → θ in L (0 , T ; H (Ω)) (cf. [3]). We decompose the firstterm as ˆ Ω T I h { ∂ t φ τh θ h } = ˆ Ω T ∂ t φ τh θ h − ˆ Ω T ( I − I h ) { ∂ t φ τh θ } . (4.18)This allows us to combine the results from (4.6) and (4.7) with (4.15b) and (4.15g)to derive (4.17a) for θ ∈ L (0 , T ; C ∞ (Ω)) . Noting that L (0 , T ; C ∞ (Ω)) is dense in L (0 , T ; H (Ω)) yields the result. Similar arguments allow us to pass to the limit in(4.12b) to obtain (4.17b).In order to pass to the limit in (4.12c), we choose θ h := I h { θ } with θ ∈ L (0 , T ; C ∞ (Ω)) and assume that γ ( φ τ, ( ± ) h ) ∈ L ∞ (0 , T ; H / (Γ)) , which is the case for κ > and κ = 0 .While the convergence of the left-hand side of (4.12c) and the gradient terms is straight-forward, the convergence of the terms including the derivative of the potential functions F and G require more finesse. We will showcase the convergence of ´ Γ T I Γ h (cid:8) G ′ + ( φ τ, + h ) θ h (cid:9) .Then, the convergence of the remaining parts can be obtained in an analogous manner.According to (P2), G ′ + can be written as the sum of a polynomial of degree three and aglobally Lipschitz-continuous component G L + ′ . We start with the decomposition ˆ Γ T I Γ h n ( φ τ, + h ) θ h o = ˆ Γ T ( φ τ, + h ) θ h − ˆ Γ T ( I − I Γ h ) n ( φ τ, + h ) o φ τ, + h θ h − ˆ Γ T ( I − I Γ h ) n I Γ h n ( φ τ, + h ) o φ τ, + h o θ h − ˆ Γ T ( I − I Γ h ) n I Γ h n ( φ τ, + h ) o θ h o . (4.19)The convergence of the first term on the right-hand side follows directly from (4.15f) andthe strong convergence of θ h → θ . Therefore, it remains to show that the remaining termsvanish when passing to the limit. Recalling that H / (Γ) is continuously embedded in L (Γ) (cf. [15]), the estimates in Lemma 2.1 and the standard inverse estimates (cf. [3]) provide (cid:12)(cid:12)(cid:12)(cid:12) ˆ Γ ( I − I Γ h ) n I Γ h n ( φ τ, + h ) o θ h o(cid:12)(cid:12)(cid:12)(cid:12) ≤ Ch (cid:13)(cid:13)(cid:13) I Γ h n ( φ τ, + h ) o(cid:13)(cid:13)(cid:13) L (Γ) k∇ Γ θ h k H (Γ) ≤ Ch (cid:13)(cid:13) φ τ, + h (cid:13)(cid:13) L (Γ) k∇ Γ θ h k H (Γ) ≤ Ch / (cid:13)(cid:13) φ τ, + h (cid:13)(cid:13) L (Γ) k∇ Γ θ h k H (Γ) . (4.20)Therefore, the last term in (4.19) vanishes. Furthermore, we derive the estimates (cid:12)(cid:12)(cid:12)(cid:12) ˆ Γ ( I − I Γ h ) n I Γ h n ( φ τ, + h ) o φ τ, + h o θ h (cid:12)(cid:12)(cid:12)(cid:12) ≤ (cid:13)(cid:13)(cid:13) ( I − I Γ h ) n I Γ h n ( φ τ, + h ) o φ τ, + h o(cid:13)(cid:13)(cid:13) L / (Γ) k θ h k H (Γ) ≤ Ch (cid:13)(cid:13) ∇ φ τ, + h (cid:13)(cid:13) L / (Γ) (cid:13)(cid:13)(cid:13) ∇I Γ h n ( φ τ, + h ) o(cid:13)(cid:13)(cid:13) L (Γ) k θ h k H (Γ) ≤ Ch / (cid:13)(cid:13) φ τ, + h (cid:13)(cid:13) L (Γ) k θ h k H (Γ) (4.21)and (cid:12)(cid:12)(cid:12)(cid:12) ˆ Γ ( I − I Γ h ) n(cid:0) φ τ, + h (cid:1) o φ τ, + h θ h (cid:12)(cid:12)(cid:12)(cid:12) ≤ (cid:13)(cid:13)(cid:13) ( I − I Γ h ) n(cid:0) φ τ, + h (cid:1) o(cid:13)(cid:13)(cid:13) L / (Γ) (cid:13)(cid:13) φ τ, + h (cid:13)(cid:13) L (Γ) k θ h k H (Γ) ≤ Ch (cid:13)(cid:13) ∇ Γ φ τ, + h (cid:13)(cid:13) L (Γ) (cid:13)(cid:13) φ τ, + h (cid:13)(cid:13) L (Γ) k θ h k H (Γ) ≤ Ch / (cid:13)(cid:13) φ τ, + h (cid:13)(cid:13) L (Γ) k θ h k H (Γ) . (4.22)As φ τ, + h ∈ L p (0 , T ; L (Γ)) , we obtain the convergence of the polynomial part of G ′ + . Todeal with the Lipschitz-continuous part G L + ′ , we start with the decomposition ˆ Γ T I Γ h n G L + ′ ( φ τ, + h ) θ h o = ˆ Γ T G L + ′ ( φ τ, + h ) θ h − ˆ Γ T ( I − I Γ h ) n G L + ′ ( φ τ, + h ) o θ h − ˆ Γ T ( I − I Γ h ) n I Γ h n G L + ′ ( φ τ, + h ) o θ h o := I + II + III . (4.23)Combining Lemma 2.1 with a standard inverse estimate, we compute | III | ≤ ˆ T Ch (cid:13)(cid:13)(cid:13) ∇ Γ I Γ h n G L + ′ ( φ τ, + h ) o(cid:13)(cid:13)(cid:13) L (Γ) k∇ Γ θ h k L (Γ) ≤ ˆ T Ch / (cid:13)(cid:13)(cid:13) I Γ h n G L + ′ ( φ τ, + h ) o(cid:13)(cid:13)(cid:13) L (Γ) k∇ Γ θ h k L (Γ) . (4.24)Using the Lipschitz-continuity of G L + ′ , we deduce (cid:13)(cid:13)(cid:13) G L + ′ ( φ τ, + h ) (cid:13)(cid:13)(cid:13) L ∞ (0 ,T ; L (Γ)) + (cid:13)(cid:13)(cid:13) I Γ h n G L + ′ ( φ τ, + h ) o(cid:13)(cid:13)(cid:13) L ∞ (0 ,T ; L (Γ)) ≤ C (cid:13)(cid:13) φ τ, + h (cid:13)(cid:13) L ∞ (0 ,T ; L (Γ)) + C , (4.25) with a constant C depending on the Lipschitz-constant of G L + ′ . Furthermore, the Lipschitz-continuity provides on every K Γ ∈ T Γ h ˆ K Γ (cid:12)(cid:12)(cid:12) I Γ h n G L + ′ ( φ τ, + h ) o − G L + ′ ( φ τ, + h ) (cid:12)(cid:12)(cid:12) ≤ C ˆ K Γ (cid:12)(cid:12)(cid:12)(cid:12) max K Γ { φ τ, + h } − min K Γ { φ τ, + h } (cid:12)(cid:12)(cid:12)(cid:12) ≤ Ch ˆ K Γ (cid:12)(cid:12) ∇ Γ φ τ, + h (cid:12)(cid:12) . (4.26)Consequently, an inverse estimate yields (cid:13)(cid:13)(cid:13) ( I − I Γ h ) n G L + ′ ( φ τ, + h ) o(cid:13)(cid:13)(cid:13) L (Γ) ≤ Ch (cid:13)(cid:13) ∇ Γ φ τ, + h (cid:13)(cid:13) L (Γ) ≤ Ch / (cid:13)(cid:13) φ τ, + h (cid:13)(cid:13) L (Γ) , (4.27)which proves that II will also vanish when passing to the limit. From the strongconvergence (4.15f), we deduce G L + ′ ( φ τ, + h ) → G L + ′ ( γ ( φ )) almost everywhere. Recalling G L + ′ ( φ τ, + h ) ∈ L ∞ (0 , T ; L (Γ)) , we may use Vitali’s convergence theorem (see e.g. [1]) toshow G L + ′ ( φ τ, + h ) → G L + ′ ( γ ( φ )) in L ∞ (0 , T ; L ˜ s (Γ)) for ˜ s < . The convergence of deriva-tives of the concave parts of G follows from the same arguments. The uniform boundsof φ τ, ( ± ) h in L ∞ (0 , T ; H (Ω)) provide enough regularity, to adapt the previously presentedarguments to three spatial dimensions, which proves the convergence of the remainingterms. Then, a denseness argument concludes the proof. (cid:3) Remark 4.5.
The results presented in the preceding sections carry over to the case ofAllen–Cahn-type dynamic boundary conditions (cf. (1.6) ), where we use ˆ Γ I Γ h (cid:8) ∂ − τ φ nh θ h (cid:9) = − m Γ ˆ Γ I Γ h (cid:8) µ n Γ ,h θ h (cid:9) for all θ h ∈ U Ω h . (4.28) instead of (2.6b) . The resulting scheme reads Φ n + τ m M − L Ω (cid:16) m L Ω (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) − ! · − m L Ω (cid:12)(cid:12) Γ × ◦ Ω M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − L Γ M − R Γ (Φ n ) M Ω (cid:12)(cid:12) − ◦ Ω × ◦ Ω R ◦ Ω (Φ n ) ! = Φ n − (4.29) and is well defined, as (cid:16) m L Ω (cid:12)(cid:12) Γ × Γ + m Γ M Ω (cid:12)(cid:12) Γ × Γ M − M Ω (cid:12)(cid:12) Γ × Γ (cid:17) is obviously a symmetric,positive definite matrix.Although, ´ Γ φ nh is not conserved when using Allen–Cahn-type boundary conditions, testing (4.28) by shows that (cid:12)(cid:12) ´ Γ φ nh (cid:12)(cid:12) is bounded. Consequently, the energy estimate still providescontrol over k φ nh k H (Γ) .Testing (4.28) by ∂ − τ φ nh shows τ P Nn =1 k ∂ − τ φ nh k L (Γ) ≤ C , i.e. we obtain a slightly betterregularity result for the discrete time derivative than we obtained for Cahn–Hilliard-typeboundary conditions. Using the time-index-free notation introduced in (4.10) , the bounds read (cid:13)(cid:13)(cid:13) φ τ, ( ± ) h (cid:13)(cid:13)(cid:13) L ∞ (0 ,T ; H (Ω)) + κ (cid:13)(cid:13)(cid:13) φ τ, ( ± ) h (cid:13)(cid:13)(cid:13) L ∞ (0 ,T ; H (Ω)) + τ − (cid:13)(cid:13) ∇ φ τ, + h − ∇ φ τ, − h (cid:13)(cid:13) L (0 ,T ; L (Ω)) + κτ − (cid:13)(cid:13) ∇ Γ φ τ, + h − ∇ Γ φ τ, − h (cid:13)(cid:13) L (0 ,T ; L (Γ)) + βτ − (cid:13)(cid:13) φ τ, + h − φ τ, − h (cid:13)(cid:13) L (0 ,T ; L (Γ)) + (cid:13)(cid:13) µ τ, + h (cid:13)(cid:13) L (0 ,T ; H (Ω)) + (cid:13)(cid:13) µ τ, +Γ ,h (cid:13)(cid:13) L (0 ,T ; L (Γ)) + k ∂ t φ τh k L (0 ,T ;( H (Ω)) ′ ) + k ∂ t φ τh k L (0 ,T ; L (Γ)) ≤ C (4.30a) with C > independent of h and τ . Based on these uniform bounds, we are able toidentify converging subsequences and pass to the limit. Appendix A. Appendix
For the reader’s convenience, we provide the generalized Poincaré inequality which canbe found in [1].
Lemma A.1.
Let Ω ⊂ R d be open, bounded and connected with Lipschitz boundary ∂ Ω .Moreover, let < p < ∞ and let M ⊂ W ,p (Ω) be nonempty, closed and convex. Thenthe following items are equivalent for every u ∈ M : (1) There exists a constant C < ∞ such that for all ξ ∈ R u + ξ ∈ M = ⇒ | ξ | ≤ C . (2) There exists a constant
C < ∞ with k u k L p (Ω) ≤ C (1 + k∇ u k L p (Ω) ) for all u ∈ M . References
1. H. W. Alt,
Linear functional analysis : An application-oriented introduction , Universitext„Birkhäuser, London, 2016.2. J.H. Bramble, J.E. Pasciak, and O. Steinbach,
On the stability of the L projection in H (Ω) , Math.Comp. (2002), no. 237, 147–156.3. S. C. Brenner and L. R. Scott, The Mathematical Theory of Finite Element Methods , Springer, 2002.4. Q. Cheng, X. Yang, and J. Shen,
Efficient and accurate numerical schemes for a hydro-dynamicallycoupled phase field diblock copolymer model , Journal of Computational Physics (2017), 44 – 60.5. L. Cherfils, S. Gatti, and A. Miranville,
A Variational Approach to a Cahn–Hilliard Model in aDomain with Nonpermeable Walls , Journal of Mathematical Sciences (2013), no. 4, 604–636.6. H. Garcke, M. Hinze, and C. Kahle,
A stable and linear time discretization for a thermodynamicallyconsistent model for two-phase incompressible flow , Applied Numerical Mathematics (2016), 151– 171.7. H. Garcke and P. Knopf, Weak solutions of the Cahn-Hilliard system with dynamic boundary condi-tions: A gradient flow approach , arXiv e-prints (2018), arXiv:1810.09817.8. G. Grün,
On convergent schemes for diffuse interface models for two-phase flow of incompressiblefluids with general mass densities , SIAM Journal on Numerical Analysis (2013), no. 6, 3036–3061.9. G. Grün, F. Guillén-González, and S. Metzger, On fully decoupled, convergent schemes for diffuseinterface models for two-phase flow with general mass densities , Communications in ComputationalPhysics (2016), no. 5, 1473–1502, DOI:10.4208/cicp.scpde14.39s.10. G. Grün and F. Klingbeil, Two-phase flow with mass density contrast: stable schemes for a thermo-dynamic consistent and frame-indifferent diffuse-interface model , Journal of Computational Physics (2014), 708–725.
11. R. Kenzler, F. Eurich, P. Maass, B. Rinn, J. Schropp, E. Bohl, and W. Dieterich,
Phase separation inconfined geometries: Solving the Cahn–Hilliard equation with generic boundary conditions , ComputerPhysics Communications (2001), no. 2, 139 – 157.12. X. Li and J. Shen,
On a SAV-MAC scheme for the Cahn-Hilliard-Navier-Stokes Phase Field Model ,arXiv e-prints (2019), arXiv:1905.08504.13. C. Liu and H. Wu,
An Energetic Variational Approach for the Cahn–Hilliard Equation with DynamicBoundary Condition: Model Derivation and Mathematical Analysis , Archive for Rational Mechanicsand Analysis (2019), no. 1, 167–247.14. S. Metzger,
On stable, dissipation reducing splitting schemes for two-phase flow of electrolyte solu-tions , Numerical Algorithms (2018), DOI: 10.1007/s11075-018-0530-2.15. E. Di Nezza, G. Palatucci, and E. Valdinoci,
Hitchhiker’s guide to the fractional Sobolev spaces ,Bulletin des Sciences Mathématiques (2012), no. 5, 521–573.16. T. Qian, X. Wang, and P. Sheng,
A variational approach to the moving contact line hydrodynamics ,Journal of Fluid Mechanics (2006), 333–360.17. J. Shen, C. Wang, X. Wang, and S. Wise,
Second-order Convex Splitting Schemes for Gradient Flowswith Ehrlich–Schwoebel Type Energy: Application to Thin Film Epitaxy , SIAM Journal on NumericalAnalysis (2012), no. 1, 105–125.18. J. Shen and X. Yang, Numerical approximations of Allen-Cahn and Cahn-Hilliard equations , Discrete& Continuous Dynamical Systems (2010), no. 4, 1669–1691.19. D. Trautwein, Finite-Elemente Approximation der Cahn–Hilliard-Gleichung mit Neumann- und dy-namischen Randbedingungen , Bachelor thesis, University of Regensburg, 2018.20. S. Wise, C. Wang, and J. Lowengrub,