# An Unfitted Discontinuous Galerkin Method for Elliptic Interface Problems

Authors: Qiuliang Wang; Jinru Chen
Journal: Journal of Applied Mathematics
 (2014)

---

## Abstract

An unfitted discontinuous Galerkin method is proposed for the elliptic interface problems. Based on a variant of the local discontinuous Galerkin method, we obtain the optimal
convergence for the exact solutionu in the energy norm and its flux p in the L2 norm. These results are the same as those in the case of elliptic problems without interface. Finally, some numerical experiments are presented to verify our theoretical results.


---

## Body

## 1. Introduction

Elliptic interface problems are often encountered in many multiphysics and multiphase applications in science computing and engineering. For example, second order elliptic equations with discontinuous coefficients are often used to model problems in material sciences and fluid dynamics when two or more distinct materials or fluids with different conductivities, densities, or permeability are involved. It is well known that, when the interface is smooth enough, the solution of elliptic interface problems has higher regularity in individual material or fluid region than in the entire physical domain.To numerically solve such interface problems, first we need to generate a mesh. One approach is to use a body fitted mesh. However, for those problems where the interface moves with time, repeated remeshing of the domain to obtain a fitted mesh is very costly. Another one is to use an unfitted grid independent of the location of the interface. This technique is particularly preferred to simulate time-dependent problems with moving interfaces. The major advantage for using an unfitted mesh is that it avoids repeatedly remeshing the domain for fitting the moving interfaces.As for fitted mesh method for elliptic problems with interface, Chen and Zou in [1] considered the finite element method for solving elliptic and parabolic interface problems, and almost-optimal error estimates in theL2 norm and energy norm were obtained. In [2], the authors studied a class of discontinuous Galerkin method for elliptic interface problems, which was shown to be optimally convergent inL2 norm. Recently, a high-order HDG method was presented to solve elliptic interface problems by Huynh et al. in [3], which was extended to solve Stokes interface flow in [4].Various unfitted grid methods for interface problems have been proposed in the literature. Finite difference methods are very popular unfitted grid methods due to their simplicity, for example, the immersed interface method [5, 6], the immersed boundary method [7], the boundary condition capturing method [8], and many others. And there exist many works for finite element methods on unfitted grid as well. Babuška in [9] studied the elliptic interface problem on unfitted mesh and derived suboptimal convergence behavior. Li et al. proposed an immersed interface finite element method in [10], which modified the basis functions near interface to satisfy the homogeneous jump conditions. Later, this method was applied to elliptic and elasticity interface problems with nonhomogeneous jump conditions in [11, 12]. Recently, an unfitted finite element method based on Nitsche’s method was presented by A. Hansbo and P. Hansbo in [13] and optimal order of convergence was proved without restrictions on the location of the interface relative to the mesh, which was used to solve incompressible elasticity with discontinuous modulus in [14]. More recently, Massjung in [15] considered anhp unfitted discontinuous Galerkin method which was viewed as a generalization of Hansbo’s method in [13]. An optimal convergence rate with respect tohand a suboptimal convergence rate with respect topin energy norm were proved. Later, Wu and Xiao also presented an unfittedhp interface penalty finite element method, which was extended to the three dimensional case in [16].The local discontinuous Galerkin (LDG) method was proposed by Cockburn and Shu in [17] to solve general time-dependent convection-diffusion problems. Later, the method was carried to elliptic problems for mixed discontinuous Galerkin formulation by Castillo et al. in [18]. The purpose of this paper is to extend the LDG method to a class of elliptic problems with a smooth interface. However, employing an unfitted mesh method, the interface can divide regular grid cells into degenerated subcells. If this situation happens, the standard inverse estimates can no longer be valid. In this paper, we use the weighted average instead of the arithmetic average in the classic LDG method to retrieve the inverse estimates (see Lemmas 8 and 10). Thus, we propose an unfitted discontinuous Galerkin method, based on a variant of LDG method. We prove the optimal convergence rate of the method for the exact solutionuin the energy norm and its fluxpin theL2norm, respectively.The rest of this paper is organized as follows. In Section2 we propose our DG method and present some necessary preliminaries. We prove optimal order error estimates for our DG method in Section 3. In Section 4, some numerical experiments are presented to justify our theoretical results. Finally, conclusions are given in Section 5.Let us now end this section with some notation to be used in this paper. We will use the standard notations for Sobolev spaces and norms in this paper (see [19, 20]). In particular, for a bounded open setΩinℝ2, if Ω=⋃i=1mΩi and Ωi∩Ωj=∅(i≠j), we denote byHk(⋃i=1mΩi)the Sobolev space of functionsusuch thatu|Ωi∈Hk(Ωi), whereHk(Ωi)denotes the standard Sobolev space with norm ∥·∥k,Ωi. As usual we define the broken norm:∥·∥k,Ω=(∑i=1m∥·∥k,Ωi2)1/2. Throughout the paper, the generic constantCis always independent of the mesh parameterh.


## 2. Discontinuous Galerkin Method and Preliminaries

LetΩbe a bounded domain inℝ2 with convex polygonal boundary∂ΩandΩ+⊂Ωan open domain withC2boundary Γ=∂Ω+⊂Ω. Let Ω-=Ω∖Ω+ (see Figure 1). We consider the following elliptic interface problem:
(1)-∇·(β∇u)=finΩ+∪Ω-,u=0on⁡∂Ω,[u]=0on⁡Γ,[β∂u∂n]=gon⁡Γ,
wherenis the outward pointing unit normal toΩ+and[u]:=u+|Γ-u-|Γis the jump ofuacross the interfaceΓ, whereu±is the restrictions ofuonΩ±. For the sake of simplicity, we assume that the coefficientβis a positive and piecewise constant; that is, β|Ω±=β±>0.Figure 1
DomainΩ, its subdomainsΩ+, Ω-, and interfaceΓ.Regarding the regularity for the solution of the interface problem (1), we state without proof the following theorem.Theorem 1 (cf. [1]).
Assume thatf∈L2(Ω)and g∈H1/2(Γ). Then problem (1) has a unique solution u∈H2(Ω+∪Ω-), and the following a priori estimate holds:
(2)∥u∥2,Ω+∪Ω-≤C(∥f∥0,Ω+∥g∥1/2,Γ).By introducing the fluxp=β∇u, the interface problem (1) can be rewritten into a first order system as
(3)1βp=∇uinΩ+∪Ω-,-∇·p=finΩ+∪Ω-,[u]=0on⁡Γ,[p·n]=gon⁡Γ,u=0on⁡∂Ω.Let𝒯hbe a shape regular and locally quasi-uniform simplicial triangulation of Ω, generated independently of the location of the interfaceΓ. For the definition of shape regular and locally quasi-uniform, we refer to [20, 21]. Suppose 𝒯hto be made of straight triangles T with diameter hT. As usual, let h∶=max⁡T∈𝒯hhT. The set of edges of the triangulation 𝒯h is denoted by ℰh, and ℰho is the set of interior edges of 𝒯h. For any elementT∈𝒯h, denote the part ofTin Ω± byT±; that is, T±=T∩Ω±. For any edge e∈ℰh, let e±=e∩Ω±. We call the elements whose interiors are cut through by Γ “interface elements”, and denote the set of the interface elements by 𝒯ho. For an interface element T∈𝒯ho, assume that ΓT is the part of interface Γ intersecting T. For the geometrical features of the interface Γ, we give the following plausible assumptions (cf. [13, 15]).Assumption 2.
We assume thatΓintersects the boundary ∂T of an element T∈𝒯ho exactly twice and each (open) edge at most once.Assumption 3.
LetΓT,h be the straight line segment connecting the points of intersection between Γ and ∂T. We assume that ΓT is a function of length on ΓT,h, in local coordinates:
(4)ΓT,h={(ξ,η)∶0<ξ<|ΓT,h|,η=0},ΓT={(ξ,η)∶0<ξ<|ΓT,h|,η=δ(ξ)}.Assumption 4.
Suppose thate=T1∩T2 with T1,T2∈𝒯h; then there exist triangles S1±⊂T1± and S2±⊂T2± such that e±=S1±∩S2± and
(5)|S1±|+|S2±|≥C|e±|2.To formulate our numerical scheme, first we define two usual discontinuous finite element spaces as(6)V~h={vh∈L2(Ω)∶vh|T∈Pk(T),∀T∈𝒯h},W~h={qh∈(L2(Ω))2∶qh|T∈Pk(T)2,∀T∈𝒯h},
where Pk(T) denotes the space of polynomials of degree less than or equal to k(k≥1) on each elementT.We define our discontinuous finite element spaces as(7)Vh={vh∶vh|Ω±=vh±,wherevh±∈V~h},Wh={qh∶qh|Ω±=qh±,whereqh±∈W~h}.Following the notation of [22], let e be an interior edge shared by two triangles T1 and T2 in 𝒯h. For a scalar valued function v, piecewise smooth on 𝒯h with vj=v|Tj, we define the jump and the weighted average ofvas
(8)⟦v±⟧=v1±n1+v2±n2,{{v±}}w=μ1±v1±+μ2±v2±,{{v±}}w=μ2±v1±+μ1±v2±on⁡e⊂ℰho.
Similarly, for a vector valued function w, piecewise smooth on 𝒯h with wj=w|Tj, we set
(9)⟦w±⟧=w1±·n1+w2±·n2,{{w±}}w=μ1±w1±+μ2±w2±,{{w±}}w=μ2±w1±+μ1±w2±on⁡e⊂ℰho,
where nj is the unit normal of e pointing towards the outside of Tj and μj±=|Sj±|/(|S1±|+|S2±|), j=1,2. If e is an edge on the boundary of Ω, we define on e(10)⟦v⟧=vn,{{w}}w=w,
where n denotes the unit outer normal of e pointing towards the outside of Ω.For the weight average across interfaceΓof any piecewise smooth function v discontinuous on Γ, we set
(11){v}w=λ+v++λ-v-,{v}w=λ-v++λ+v-,
where λ++λ-=1 whose specific definitions will be given in Lemma 10.For simplicity, forT±=T∩Ω± and e±=e∩Ω±, we define
(12)〈v,w〉T:=(v+,w+)T++(v-,w-)T-,〈v,w〉e:=(v+,w+)e++(v-,w-)e-.Following [18], testing the problem (3) by qh∈Wh and vh∈Vh, respectively, using integration by parts and noting the identities ⟦qv⟧=⟦v⟧{{q}}w+{{v}}w⟦q⟧ and [vw]={w}w[v]+{v}w[w], we obtain our DG method: find (qh,uh)∈Wh×Vh, such that
(13)∑T∈𝒯h‍〈1βph,qh〉T+∑T∈𝒯h‍〈uh,∇·qh〉T-∑e∈ℰho‍〈{{uh}}w,⟦qh⟧〉e-∑T∈𝒯ho‍({uh}w,[qh·n])ΓT=0,∑T∈𝒯h‍〈ph,∇vh〉T-∑e∈ℰh‍〈{{ph}}w,⟦vh⟧〉e-∑T∈𝒯ho‍({ph·n}w,[vh])ΓT+∑e∈ℰh‍〈γ|e±|⟦uh⟧,⟦vh⟧〉e+∑T∈𝒯ho‍γhT([uh],[vh])ΓT=∑T∈𝒯h‍〈f,vh〉T+∑T∈𝒯ho‍(g,{vh}w)ΓT,
for all (qh,vh)∈Wh×Vh, where γ>0 is the stabilization parameter.We define the bilinear and linear forms(14)a(ph,qh)=∑T∈𝒯h‍〈1βph,qh〉T,b(vh,ph)=∑T∈𝒯h‍〈ph,∇vh〉T-∑e∈ℰh‍〈{{ph}}w,⟦vh⟧〉emmmmmm-∑T∈𝒯ho‍({ph·n}w,[vh])ΓT,c(uh,vh)=∑e∈ℰh‍〈γ|e±|⟦uh⟧,⟦vh⟧〉emmmmmm+∑T∈𝒯ho‍γhT([uh],[vh])ΓT,F(vh)=∑T∈𝒯h‍〈f,vh〉T+∑T∈𝒯ho‍(g,{vh}w)ΓT.
And integral by parts yields that
(15)b(uh,qh)=-∑T∈𝒯h‍〈uh,∇·qh〉T+∑e∈ℰho‍〈{{uh}}w,⟦qh⟧〉e+∑T∈𝒯ho‍({uh}w,[qh·n])ΓT.Hence, our DG approximation can be written as the following mixed variational problem: find(qh,uh)∈Wh×Vh, such that
(16)a(ph,qh)-b(uh,qh)=0∀qh∈Wh,b(ph,vh)+c(uh,vh)=F(vh)∀uh∈Vh.
For the exactuof the interface problem (1) andp=β∇u, using Theorem 1, we have ⟦u⟧=0, ⟦p⟧=0  on  e∈ℰho, and[u]=0, [p·n]=g  across Γ. Then the following consistency property holds:
(17)a(p,qh)-b(u,qh)=0,b(p,vh)+c(u,vh)=F(vh).Let the mesh-dependent norm|∥·∥|hbe defined by
(18)|∥vh∥|h2=∑T∈𝒯h‍∥∇vh∥0,T±2+∑e∈ℰh‍1|e±|∥⟦vh⟧∥0,e±2+∑T∈𝒯ho‍1hT∥[vh]∥0,ΓT2.Theorem 5.
Suppose that the stabilization parameterγ is positive; then the DG method (16) defines a unique approximate solution (ph,uh)∈Wh×Vh.Proof.
Since (16) is a square system, it is enough to show uniqueness. Let f=0, g=0. Setting qh=ph andvh=uh, adding the two equations of (16), we have
(19)a(ph,ph)+c(uh,uh)=0,
which deducesph≡0,⟦uh⟧=0 on ℰhand[uh]=0 across Γ. As a consequence, the first equation of (16) becomes
(20)∑T∈𝒯h‍〈∇uh,qh〉=0∀qh∈Wh.
Hence, taking qh=∇uh implies ∇uh=0. Since⟦uh⟧=0 on ℰh and [uh]=0 across Γ, we conclude that uh≡0. This completes the proof.The following lemma comes from the famous Stein’s extension theorem.Lemma 6 (cf. [19]).
There exist two extension operatorsE±∶Hk(Ω±)→Hk(Ω)for all nonnegative integers k such that
(21)(E±v)|Ω±=v,∥E±v∥Hk(Ω)≤C∥v∥Hk(Ω±),
where v∈Hk(Ω±).Next, we state a standard approximation lemma.Lemma 7 (cf. [20, 21]).
Letu∈Hk+1(T). Then form=0,1there exists a linear continuous operatorΠT∶Hk+1(T)→Pk(T)such that
(22)∥u-ΠTu∥m,T≤ChTk+1-m∥u∥k+1,T,(23)∥u-ΠTu∥L∞(T)≤ChTk∥u∥k+1,T,(24)∥u-ΠTu∥0,∂T≤ChTk+1/2∥u∥k+1,T.The following two lemmas are variant inverse estimates involving interfaceΓ which play important role in our analysis.Lemma 8 (cf. [15]).
Fore∈ℰh let e=T1∩T2; then for q∈Wh the following inverse inequality holds:
(25)∥{{q}}w∥0,e±2≤C|e±|-1∑i=12‍∥q∥0,Ti±2.Lemma 9 (cf. [15, 16]).
The following estimate holds for eitheri=+ or i=-, for any v∈Pk(T):
(26)∥v∥0,ΓT≤ChT-1/2∥v∥0,Ti.By Lemma9, we can immediately obtain the following result.Lemma 10.
LetT∈𝒯ho; then there exists a positive constant C such that
(27)λ±∥v∥0,ΓT≤ChT-1/2∥v∥0,T±,∀v∈Pk(T),
where
(28)λ±={12if(26)holdsforbothT±andT∓,1if(26)holdsonlyforT±,0if(26)holdsonlyforT∓.Lemma 11.
LetT∈𝒯hoandv∈H1(T); then we have
(29)∥v∥0,ΓT≤C(hT-1/2∥v∥0,T+hT1/2|v|1,T).Proof.
Under Assumptions2 and 3, the trace inequality (29) follows from Lemma  3 in [13] and a scaling argument.


## 3. Error Estimate of Our DG Method

Now, we define interpolation operatorΠ by Π|T∶=ΠT and Π∶=(Π,Π). Let ξp±=p±-ΠE±p±,ηp±=ΠE±p±-ph± and ξu±=u±-ΠE±u±,ηu±=ΠE±u±-uh±. To obtain the convergence result, we need to show the following approximate error estimates.Lemma 12.
Suppose thatu∈Hk+1(Ω+∪Ω-) and p∈(Hk+1(Ω+∪Ω-))2; then the following approximate estimates hold:
(30)∑T∈𝒯h‍∥ξu∥0,T±2+∑e∈ℰh‍1|e±|∥⟦ξu⟧∥0,e±2+∑T∈𝒯ho‍1hT∥[ξu]∥0,ΓT2≤Ch2k∥u∥k+1,Ω+∪Ω-2,(31)∑T∈𝒯h‍∥ξp∥0,T±2+∑e∈ℰh‍|e±|∥{{ξp}}w∥0,e±2+∑T∈𝒯ho‍hT∥{ξp·n}w∥0,ΓT2≤Ch2k∥p∥k,Ω+∪Ω-2.Proof.
We only need to show inequality (30); inequality (31) can be shown similarly. For the first term on the left-hand side of (30), by the inequality (22) from Lemma 7 we obtain
(32)∥ξu∥0,T±2=∥E±u±-ΠE±u±∥0,T±2≤∥E±u±-ΠE±u±∥0,T2≤Ch2k+2∥E±u±∥k+1,T2.
Summing over all triangles, it follows by Lemma 6 that
(33)∑T∈𝒯h‍∥ξu∥0,T±2≤Ch2k+2∥u∥k+1,Ω+∪Ω-2.
Next, we estimate the second term on the left-hand side of (30) as follows. Suppose that e=T1∩T2 with T1,T2∈𝒯h; using the inequality (23) from Lemma 7 yields
(34)1|e±|∥⟦ξu⟧∥0,e±2=1|e±|∥⟦E±u±-ΠE±u±⟧∥0,e±2≤∥⟦E±u±-ΠE±u±⟧∥L∞(e±)2≤Ch2k∥E±u±∥k+1,T1∪T22.
Summing over all edges, using Lemma 6 implies that
(35)∑e∈ℰh‍1|e±|∥⟦ξu⟧∥0,e±2≤Ch2k∥u∥k+1,Ω+∪Ω-2.
Similarly, due to the inequality (24) from Lemma 7 and Lemmas 6 and 11, we find
(36)∑T∈𝒯ho‍1hT∥[ξu]∥0,ΓT2≤Ch2k∥u∥k+1,Ω+∪Ω-2.
Thus the inequality (30) follows combining (33)–(36), which completes the proof.Next we present a priori error estimate of the exact solutionu in the energy norm |∥·∥|h and its fluxpin theL2-norm.Theorem 13.
Let(p,u) be the solution of (3) and (ph,uh) the solution of (16), respectively. Then, for (p,u)∈(Hk(Ω+∪Ω-))2×Hk+1(Ω+∪Ω-), the following error estimate holds:
(37)∥p-ph∥0,Ω+∪Ω-+|∥u-uh∥|h≤Chk(∥p∥k,Ω+∪Ω-+∥u∥k+1,Ω+∪Ω-).Proof.
By using the consistency property (17), we have
(38)a(ηp,qh)-b(ηu,qh)=a(ξp,qh)-b(ξu,qh),b(ηp,vh)+c(ηu,vh)=b(ξp,vh)+c(ξu,vh).
Setqh=ηp, vh=ηu in (38) to give
(39)a(ηp,ηp)+c(ηu,ηu)=a(ξp,ηp)-b(ξu,ηp)+b(ξp,ηu)+c(ξu,ηu).
Using (31) in Lemma 12 and ε-Cauchy-Schwartz inequality, we estimate the first term on the right-hand side of (39) as
(40)a(ξp,ηp)=∑T∈𝒯h‍〈ξp,ηp〉T=∑T∈𝒯h‍(ξp,ηp)T±≤(∑T∈𝒯h‍∥ξp∥0,T±2)1/2(∑T∈𝒯h‍∥ηp∥0,T±2)1/2≤Ch2k∥p∥k,Ω+∪Ω-2+ε∑T∈𝒯h‍∥ηp∥0,T±2.
For the second term on the right-hand side of (39), by using (30) in Lemma 12 and trace inequalities from Lemmas 8 and 10, we obtain
(41)b(ξu,ηp)=∑T∈𝒯h‍〈ηp,∇ξu〉T-∑e∈ℰh‍〈{{ηp}}w,⟦ξu⟧〉e-∑T∈𝒯ho‍({ηp·n}w,[ξu])ΓT=∑T∈𝒯h‍(ηp,∇ξu)T±-∑e∈ℰh‍({{ηp}}w,⟦ξu⟧)e±-∑T∈𝒯ho‍({ηp·n}w,[ξu])ΓT≤(∑T∈𝒯h‍∥ηp∥0,T±2)1/2(∑T∈𝒯h‍∥∇ξu∥0,T±2)1/2+(∑e∈ℰh‍|e±|∥{{ηp}}w∥0,e±2)1/2×(∑e∈ℰh‍1|e±|∥⟦ξu⟧∥0,e±2)1/2+(∑T∈𝒯ho‍hT∥{ηp·n}w∥0,ΓT2)1/2(∑T∈𝒯ho‍1hT∥[ξu]∥0,ΓT2)1/2≤Ch2k∥u∥k+1,Ω+∪Ω-2+ε∑T∈𝒯h‍∥ηp∥0,T±2.
Similarly, by using the approximate estimation (31), the third term on (39) can be estimated as
(42)b(ξp,ηu)=∑T∈𝒯h‍〈ξp,∇ηu〉T-∑e∈ℰh‍〈{{ξp}}w,⟦ηu⟧〉e-∑T∈𝒯ho‍({ξp·n}w,[ηu])ΓT=∑T∈𝒯h‍(ξp,∇ηu)T±-∑e∈ℰh‍({{ξp}}w,⟦ηu⟧)e±-∑T∈𝒯ho‍({ξp·n}w,[ηu])ΓT≤(∑T∈𝒯h‍∥ξp∥0,T±2)1/2(∑T∈𝒯h‍∥∇ηu∥0,T±2)1/2+(∑e∈ℰh‍|e±|∥{{ξp}}w∥0,e±2)1/2×(∑e∈ℰh‍1|e±|∥⟦ηu⟧∥0,e±2)1/2+(∑T∈𝒯ho‍hT∥{ξp·n}w∥0,ΓT2)1/2×(∑T∈𝒯ho‍1hT∥[ηu]∥0,ΓT2)1/2≤Ch2k∥p∥k,Ω+∪Ω-2+ε(∑e∈ℰh‍1|e±|0,e±∥⟦ηu⟧∥0,e±2+∑T∈𝒯ho‍1hT∥[ηu]∥0,ΓT2mmmmm+∑T∈𝒯h‍∥∇ηu∥0,T±2).
Then, by virtue of inequality (30) in Lemma 12, we bound the fourth term on (39) as
(43)c(ξu,ηu)=∑e∈ℰh‍〈γ|e±|⟦ξu⟧,⟦ηu⟧〉e+∑T∈𝒯ho‍γhT([ξu],[ηu])ΓT=∑e∈ℰh‍(γ|e±|⟦ξu⟧,⟦ηu⟧)e±+∑T∈𝒯ho‍γhT([ξu],[ηu])ΓT≤(∑e∈ℰh‍γ|e±|0,e±∥⟦ξu⟧∥0,e±2)1/2×(∑e∈ℰh‍γ|e±|∥⟦ηu⟧∥0,e±2)1/2+(∑T∈𝒯ho‍γhT∥[ξu]∥0,ΓT2)1/2×(∑T∈𝒯ho‍γhT∥[ηu]∥0,ΓT2)1/2≤Ch2k∥u∥2k+1,Ω+∪Ω-+εc(ηu,ηu).
Thus, combining (39)–(43) yields
(44)a(ηp,ηp)+c(ηu,ηu)≤Ch2k(∥p∥k,Ω+∪Ω-2+∥u∥k+1,Ω+∪Ω-2)+ε∑T∈𝒯h‍∥∇ηu∥0,T±2.
By using the triangle inequality, we have
(45)∥p-ph∥0,Ω+∪Ω-2+c(ηu,ηu)≤Ch2k(∥p∥k,Ω+∪Ω-2+∥u∥k+1,Ω+∪Ω-2)+ε∑T∈𝒯h‍∥∇ηu∥0,T±2.
At the other hand, settingqh=∇ηu in the first equality of (38), we obtain
(46)a(ηp,∇ηu)-b(ηu,∇ηu)=a(ξp,∇ηu)-b(ξu,∇ηu).
By the definition of b(·,·), an integration by parts implies that
(47)∑T∈𝒯h‍∥∇ηu∥0,T±2=a(p-ph,∇ηu)+b(ξu,∇ηu)+∑e∈ℰh‍〈{{∇ηu}}w,⟦ηu⟧〉e+∑T∈𝒯ho‍({∇ηu·n}w,[ηu])ΓT.
Using Lemmas 8 and 10 yields
(48)∑T∈𝒯h‍∥∇ηu∥0,T±2≤C∥p-ph∥0,Ω+Ω-2+Ch2k∥u∥k+1,Ω+∪Ω-+c(ηu,ηu)+ε∑T∈𝒯h‍∥∇ηu∥0,T±2.
Combining (45) and (48), choosingεenough small, from the definition |∥·∥|h and the triangle inequality, we can arrive at
(49)∥p-ph∥0,Ω+∪Ω-2+|∥u-uh∥|h2≤Ch2k(∥p∥k,Ω+∪Ω-2+∥u∥k+1,Ω+∪Ω-2),
which completes the proof.Using the standard duality argument, we can obtain the following error estimate in theL2norm.Theorem 14.
Under the condition of Theorem13, we have
(50)∥u-uh∥0,Ω+∪Ω-≤Chk+1(∥p∥k,Ω+∪Ω-+∥u∥k+1,Ω+∪Ω-).Proof.
Consider the following so-called adjoint problem(51)-∇·(β∇ϕ)=u-uhinΩ1∪Ω2,ϕ=0on⁡∂Ω,[ϕ]=0on⁡Γ,[β∂ϕ∂n]=0on⁡Γ.
Using Theorem 1, we get
(52)∥ϕ∥2,Ω+∪Ω-≤C∥u-uh∥0,Ω+∪Ω-.
By introducing an auxiliary variable q=β∇ϕ, we obtain
(53)1βq=∇ϕinΩ+∪Ω-,-∇·q=u-uhinΩ+∪Ω-,ϕ=0on∂Ω,[ϕ]=0on⁡Γ,[q·n]=0on⁡Γ.
Since [ϕ]=0  on  Γ and ⟦ϕ⟧=0  on  ℰh, using integration by parts and the consistency property (17), we deduce
(54)∥u-uh∥0,Ω+∪Ω-2=(u-uh,u-uh)=a(q,p-ph)-b(ϕ,p-ph)+b(u-uh,q)+c(ϕ,u-uh)=∑T∈𝒯h‍〈1βξq,p-ph〉T-∑T∈𝒯h‍〈∇ξϕ,p-ph〉T+∑e∈ℰh‍〈{{p-ph}}w,⟦ξϕ⟧〉e+∑T∈𝒯ho‍({(p-ph)·n}w,[ξϕ])ΓT+∑T∈𝒯h‍〈ξq,∇(u-uh)〉T-∑e∈ℰh‍〈{{ξq}}w,⟦u-uh⟧〉e-∑T∈𝒯ho‍({ξq·n}w,[u-uh])ΓT+∑e∈ℰh‍〈1|e±|⟦ξϕ⟧,⟦u-uh⟧〉e+∑T∈𝒯ho‍1hT([ξϕ],[u-uh])ΓT=∑i=19‍Ei.
Due to q=β∇ϕ, as in the proof of Theorem 13, an application of Theorem 1 implies that
(55)|E1|+|E2|+∑i=59‍|Ei|≤Ch∥ϕ∥2,Ω+∪Ω-(∥p-ph∥0,Ω+∪Ω-+|∥u-uh∥|h).
Now we estimate E3 as follows:
(56)E3=∑e∈ℰh‍〈{{ξp}}w,⟦ξϕ⟧〉e+∑e∈ℰh‍〈{{ηp}}w,⟦ξϕ⟧〉e=(∑e∈ℰh‍|e±|∥{{ξp}}w∥0,e±2)1/2(∑e∈ℰh‍1|e±|∥⟦ξϕ⟧∥0,e±2)1/2+(∑e∈ℰh‍|e±|∥{{ηp}}w∥0,e±2)1/2(∑e∈ℰh‍1|e±|∥⟦ξϕ⟧∥0,e±2)1/2≤Ch(hk∥p∥k,Ω+∪Ω-+∥p-ph∥0,Ω+∪Ω-)∥ϕ∥2,Ω+∪Ω-.
Similarly,
(57)|E4|≤Ch(hk∥p∥k,Ω+∪Ω-+∥p-ph∥0,Ω+∪Ω-)∥ϕ∥2,Ω+∪Ω-.
We now combine (54)–(57) and use (52) and Theorem 13 to obtain
(58)∥u-uh∥0,Ω+∪Ω-≤Chk+1(∥p∥k,Ω+∪Ω-+∥u∥k+1,Ω+∪Ω-),
which completes the proof.


## 4. Numerical Experiments

In this section, we consider the following nontrivial example with a homogeneous jump condition (from [10, 13]). The exact solution is given by
(59)u(x,y)={r2β1ifr≤r0,r2β2-r02β2+r02β1,ifr>r0,
where r=(x-0.5)2+(y-0.5)2, and on the domain Ω=(0,1)×(0,1) we choose r0=0.3. A simple calculation shows that p(x,y)=(x,y)T.We compute the order of convergence foreu=u-uh and ep=p-ph, when piecewise linear polynomials are used to approximateuand p, respectively. In Figure 2, we plot the numerical solution for this example on a 64×64 uniform mesh for β1:β2=1 : 1000. Figures 3, 4, 5, and 6 show the computed order of convergence for ∥eu∥0,Ω1∪Ω2 and ∥ep∥0,Ω1∪Ω2 when the jump in the coefficient is taken as β1:β2=1 : 10, 1 : 1000, 10 : 1, 1000 : 1, respectively, in the log-log scale. These computed results coincide with the theoretical results in Theorems 13 and 14.Figure 2
The numerical solution on a64×64uniform mesh forβ1∶β2= 1 : 1000.Figure 3
The convergence rates ofL2 norm in the exactuand its fluxpforβ1∶β2= 1 : 10.Figure 4
The convergence rates ofL2 norm in the exactuand its flux p for β1∶β2= 1 : 1000.Figure 5
The convergence rates ofL2 norm in the exactuand its flux p for β1∶β2= 10 : 1.Figure 6
The convergence rates ofL2 norm in the exactuand its flux p for β1∶β2= 1000 : 1.


## 5. Conclusions

In this paper, we have discussed an unfitted discontinuous Galerkin method for elliptic problems with a smooth interface. Based on a variant of local discontinuous Galerkin method, we have obtained the optimal order error estimates in the energy norm inuand its flux p. And by using the standard duality argument the optimal convergence rate in L2norm for uhas also been derived. These presented results are the same as that of elliptic problems without interface. Finally, numerical experiments are given to confirm our theoretical results. We note that the convergence behavior in most existing works concerning the elliptic interface problems depends on the jump in the discontinuous coefficients. It will be one of our future subjects to design an efficient numerical scheme that is robust with respect to the jump in the discontinuous coefficients.



---

*Source: 241890-2014-02-20.xml*