[go: up one dir, main page]
More Web Proxy on the site http://driver.im/ Skip to content
Publicly Available Published by De Gruyter September 30, 2017

Combining the DPG Method with Finite Elements

  • Thomas Führer , Norbert Heuer EMAIL logo , Michael Karkulik and Rodolfo Rodríguez

Abstract

We propose and analyze a discretization scheme that combines the discontinuous Petrov–Galerkin and finite element methods. The underlying model problem is of general diffusion-advection-reaction type on bounded domains, with decomposition into two sub-domains. We propose a heterogeneous variational formulation that is of the ultra-weak (Petrov–Galerkin) form with broken test space in one part, and of Bubnov–Galerkin form in the other. A standard discretization with conforming approximation spaces and appropriate test spaces (optimal test functions for the ultra-weak part and standard test functions for the Bubnov–Galerkin part) gives rise to a coupled DPG-FEM scheme. We prove its well-posedness and quasi-optimal convergence. Numerical results confirm expected convergence orders.

MSC 2010: 65N30; 35J20

1 Introduction

The discontinuous Petrov–Galerkin method with optimal test functions (DPG method) is an approximation scheme that makes the use of optimal test functions, cf. [1, 6, 8], feasible by considering broken test norms [9]. Optimal test functions are those which maximize discrete inf-sup numbers, and the broken form of test spaces and norms allows for their local calculation or approximation. In this form, the DPG method has been developed by Demkowicz and Gopalakrishnan; see [8, 9].

The DPG method has been designed having in mind problems where standard methods suffer from locking phenomena (of small inf-sup numbers) or, otherwise, require specific stabilization techniques. This is particularly the case with singularly perturbed problems where DPG schemes have made some contributions [11, 7, 2, 3, 17]. Nevertheless, in the current form most of the schemes are not cheap to implement. On the one hand, corresponding formulations have several unknowns as is the case with mixed finite elements. On the other hand, the efficient approximation of optimal test functions for singularly perturbed problems is ongoing research. For these reasons, advanced DPG techniques are best used for specific problems whereas finite elements are hard to beat when solving uniformly well-posed problems. Though, it has to be said, that in the latter cases DPG schemes can also be efficient and are competitive in general; cf. the software package developed by Roberts [18].

In this paper, we develop a discretization method that combines DPG techniques with standard finite elements. In this way, one can restrict the use of more expensive DPG approximations to regions where they are beneficial. Examples are reaction-advection-diffusion problems with small diffusivity in a reduced area, or transmission problems that couple a singularly perturbed problem with an unperturbed problem. In a previous publication [14], we have proposed such a combination with boundary elements to solve transmission problems of the Laplacian in the full space, and studied a singularly perturbed case of reaction diffusion in [13]. In this paper, we follow the general framework from [14]. There, the basis is set by a heterogeneous variational formulation consisting of an ultra-weak one in a bounded domain and variational boundary integral equations for the exterior unbounded part. Here, we combine an ultra-weak formulation with a standard variational form. We remark that this approach of combining different variational formulations has been systematically analyzed in [12]. Indeed, it is not essential to use an ultra-weak formulation for the DPG scheme, any well-posed formulation would work. Though, the overall strategy in [12] is to employ DPG techniques throughout whereas we combine different discretization techniques.

An important aspect of the DPG method is the inherent guaranteed error control by the computed residuals. For early computational experiments with adaptivity based on these residuals, see [10]. In [4], an a posteriori error analysis is given that includes data approximation errors. When coupling the DPG method with other discretizations, the inherent residuals for the DPG approximation have to be combined with appropriate estimators. For the case of boundary elements, see [14, 13], and for the DPG method dealing with contact conditions, see [15]. In this paper, we do not specifically deal with a posteriori error estimation.

Having set our heterogeneous formulation, we proceed to rewrite it by using the so-called trial-to-test operator (which maps the test space to the ansatz space). This is only done for the ultra-weak formulation. The whole system then transforms into one where spaces on the ansatz and test sides are identical. In this way, our heterogeneous variational formulation fits the Lax–Milgram framework just as in [14]. We prove coercivity under the condition that the trial-to-test operator is weighted by a sufficiently large constant. Then quasi-optimal convergence of a discretized version follows by standard arguments. When proving coercivity we follow steps that are similar to the ones when studying the coupling of DPG with boundary elements. But whereas [14] analyzes only the Laplacian, here we set up the scheme and prove coercivity for a general second-order equation of reaction-advection-diffusion type. Throughout, we assume that our problem is uniformly well posed, i.e., we do not study variations for singularly perturbed cases as in [13]. Also note that, since coefficients are variable, transmission problems can be treated the same way by selecting the sub-domains accordingly. One only has to move the possibly non-homogeneous jump data to the right-hand side functional.

The remainder of this paper is as follows: In Section 2, we start by formulating the model problem. A heterogeneous variational formulation is given in Section 2.1. There, we also state its well-posedness and coercivity (Theorem 2.2) and briefly mention a simplified case where continuity across the sub-domain interface is incorporated strongly (Corollary 2.4). The corresponding discrete DPG-FEM scheme is presented in Section 2.2. Its quasi-optimal convergence is announced in Theorem 2.5. Most technical details and proofs are given in Section 3. In Section 4, we report on some numerical experiments.

Furthermore, throughout the paper, suprema are taken over sets excluding the null element, and the notation AB is used to say that ACB with a constant C>0 which does not depend on any quantity of interest. Correspondingly, the notation AB is used.

2 Mathematical Setting and Main Results

Let Ωd, d{2,3}, be a bounded, simply connected polygonal/polyhedral Lipschitz domain with boundary Ω, and normal vector 𝐧Ω on Ω pointing outside of Ω. We consider the following elliptic problem of diffusion-advection-reaction type: Given fL2(Ω), find uH01(Ω) such that

(2.1)Au:=div(-α¯u+𝜷u)+γu=fin Ω.

Here, L2(Ω) and H01(Ω) denote standard Sobolev spaces, the latter with zero trace on Ω. Furthermore, the coefficients α¯(x)d×d, 𝜷(x)d, γ(x) (xΩ¯) are assumed to be uniformly bounded. Furthermore, we require that the symmetric part of α¯ is positive definite and uniformly bounded from below, with minimum eigenvalue larger than or equal to α0>0, and that 12div𝜷+γ0 in Ω. These conditions imply that the operator A is bounded and coercive on H01(Ω).

2.1 Heterogeneous Variational Formulation

In order to solve (2.1) by a combination of DPG method and finite elements, we formulate the problem in a heterogeneous way, using different variational forms in different parts of the domain. For ease of illustration, we restrict ourselves to two Lipschitz sub-domains Ω1, Ω2 (again of polygonal/polyhedral form, each with one connected component) with boundaries Ω1, Ω2, as specified in Figure 1. There, also a notation for the boundary pieces is introduced. In particular, Γ denotes the interface between the sub-domains. The picture indicates that both sub-domains touch the boundary of Ω (where the homogeneous Dirichlet condition is imposed), but this is not essential. For instance, one sub-domain, Ω2, can be of annular type so that, in this case, ΩΩ2 and Γ=Ω1. Other combinations can be analyzed without difficulty, also including Neumann conditions. Nevertheless, since our analysis centers around proving coercivity of bilinear forms, we need positivity of the combined advection-reaction term on a sub-domain that does not touch the Dirichlet boundary.

Assumption 2.1.

For i=1,2, the following condition holds: If meas(Γi)=0, then there exists β>0 such that 12div𝜷+γβ a.e. in Ωi.

Figure 1

Decomposition of the domain Ω into sub-domains. Ω¯=Ω¯1Ω¯2, Ω1Ω2=, Ω¯1Ω¯2=Γ, Ω¯1Ω=Γ¯1, Ω¯2Ω=Γ¯2.

Standard and broken Sobolev spaces. Essential for the DPG method is to use broken test spaces. Therefore, at this early stage we consider a partitioning 𝒯1 of Ω1 into (regular non-intersecting) finite elements T such that Ω¯1={¯T;T𝒯1}, and with skeleton 𝒮:={T;T𝒯1}.

Before describing the variational formulation, we introduce the Sobolev spaces we need. For a domain ωΩ, we use the standard spaces L2(ω), H1(ω), H01(ω), and 𝐇(div,ω). The trace operator acting on H1(ω) will be denoted simply by ()|ω. Then we define the trace space H1/2(ω):=H1(ω)|ω and its dual space H-1/2(ω):=(H1/2(ω)) with canonical norms. The duality pairing on ω is ,ω and extends the L2(ω) bilinear form. Correspondingly, (,)ω is the L2(ω) bilinear form.

We also need HD1(Ωi) consisting of H1-functions with vanishing trace on Γi (i=1,2). Vector-valued spaces and functions will be denoted by bold symbols. Connected with 𝒯1 we use the product spaces H1(𝒯1) and 𝐇(div,𝒯1) with corresponding product norms.

Now, related with 𝒯1 are the skeleton trace spaces

H1/2(𝒮):={u^ΠT𝒯1H1/2(T); there exists wH1(Ω) such that u^|T=w|T for all T𝒯1},
H-1/2(𝒮):={σ^ΠT𝒯1H-1/2(T); there exists 𝒒𝐇(div,Ω) such that σ^|T=(𝒒𝐧T)|T for all T𝒯1}

and

H001/2(𝒮):={u^H1/2(𝒮);u^|Ω1=0},
HD1/2(𝒮):={u^H1/2(𝒮);u^|Γ1=0}.

Here, 𝐧T is the exterior unit normal vector on T, and (𝒒𝐧T)|T indicates the standard way of defining normal traces of 𝐇(div,T)-functions. The notation u^|Ω1=0 (resp. u^|Γ1=0) is to be understood in the sense that u^ is a 𝒯1-piecewise trace of an element of H01(Ω1) (resp. of HD1(Ω1)). These trace spaces are equipped with the norms

u^H1/2(𝒮):=inf{wH1(Ω);wH1(Ω) such that u^|T=w|T for all T𝒯1},
σ^H-1/2(𝒮):=inf{𝒒𝐇(div,Ω);𝒒𝐇(div,Ω) such that σ^|T=(𝒒𝐧T)|T for all T𝒯1},

and analogously for H001/2(𝒮) and HD1/2(𝒮). For functions u^H1/2(𝒮), σ^H-1/2(𝒮) (they are elements of product spaces) and 𝝉𝐇(div,𝒯1), vH1(𝒯1), we use the duality pairings

u^,𝝉𝐧𝒮:=T𝒯1u^|T,𝝉𝐧TT,σ^,v𝒮:=T𝒯1σ^|T,vT.

Heterogeneous formulation in Ω1Ω2. In Ω1, where the DPG method will be used, we consider an ultra-weak variational formulation. As mentioned before, this is just for illustration as any other formulation of primal, mixed, dual-mixed or strong type can be used and analyzed analogously to our case; cf. [12, Section 2.3].

The ultra-weak formulation requires additional independent unknowns

(2.2)𝝈:=α¯u-𝜷u on Ω1,u^:=ΠT𝒯1u|T,σ^:=ΠT𝒯1(𝝈𝐧T)|T.

Then we test the defining relation of 𝝈 with α¯-T𝝉 for 𝝉𝐇(div,𝒯1), and equation (2.1) with vH1(𝒯1). Integrating by parts element-wise and substituting the corresponding terms by 𝝈, u^, and σ^, we obtain

(2.3)(u,div𝒯𝝉+𝜷α¯-T𝝉+γv)Ω1+(𝝈,𝒯v+α¯-T𝝉)Ω1-u^,𝝉𝐧𝒮-σ^,v𝒮=(f,v)Ω1.

Here, div𝒯 and 𝒯 denote the 𝒯1-piecewise divergence and gradient operators, respectively.

In Ω2, we use the standard primal formulation

(2.4)(α¯u-𝜷u,w)Ω2+(γu,w)Ω2-𝐧Ω2(α¯u-𝜷u),wΩ2=(f,w)Ω2

for wHD1(Ω2).

Solving (2.1) in Ω is equivalent to solving (in appropriate spaces) (2.3) and (2.4) with homogeneous Dirichlet condition on Ω and transmission conditions on Γ. These transmission conditions will be imposed in variational form. For the time being, we replace 𝐧Ω2(α¯u-𝜷u)|Γ by -σ^|Γ in (2.4). Here, we slightly abuse the notation of σ^ noting that σ^,v𝒮=σ^,vΓ for vH1(Ω1) with v|Γ1=0; cf., e.g., [13, Section 2.2].

We formally distinguish between u1:=u|Ω1 and u2:=u|Ω2. Then our preliminary heterogeneous variational formulation consists in finding

(𝒖,u2)=(u1,𝝈,u^,σ^,u2)U:=U1×HD1(Ω2)withU1:=L2(Ω1)×𝐋2(Ω1)×HD1/2(𝒮)×H-1/2(𝒮)

such that u^|Γ=u2|Γ and

(u1,div𝒯𝝉+𝜷α¯-T𝝉+γv)Ω1+(𝝈,𝒯v+α¯-T𝝉)Ω1-u^,𝝉𝐧𝒮-σ^,v𝒮=(f,v)Ω1,
(2.5)(α¯u2-𝜷u2,w)Ω2+(γu2,w)Ω2+σ^,wΓ=(f,w)Ω2

for any (𝒗,w)V×HD1(Ω2) with

𝒗=(v,𝝉)andV:=H1(𝒯1)×𝐇(div,𝒯1).

This formulation can be used to define the combined DPG-FEM discretization, but requires that 𝒯1 is compatible across Γ with the finite element mesh in Ω2. We therefore replace the continuity constraint u^|Γ=u2|Γ by a variational coupling on Γ that is similar to a DG-bilinear form involving jumps and fluxes across element boundaries. To this end, we abbreviate

b(𝒖,𝒗):=(u1,div𝒯𝝉+𝜷α¯-T𝝉+γv)Ω1+(𝝈,𝒯v+α¯-T𝝉)Ω1-u^,𝝉𝐧𝒮-σ^,v𝒮,
(2.6)c2(u2,w2):=(α¯u2-𝜷u2,w2)Ω2+(γu2,w2)Ω2,
L1(𝒗):=(f,v)Ω1,L2(w2):=(f,w2)Ω2,

and define the coupling bilinear form

(2.7)d(𝒖,u2;𝒘,w2):=σ^,w2Γ+χ^,u^-u2Γ+12𝜷𝐧Ω1(u^-u2),w^+w2Γ
for(𝒖,u2),(𝒘,w2)Uwith𝒖=(u1,𝝈,u^,σ^),𝒘=(w1,𝝌,w^,χ^).

We will also need the bilinear form for Ω1 that corresponds to c2(,):

(2.8)c1(u1,w1):=(α¯u1-𝜷u1,w1)Ω1+(γu1,w1)Ω1(u1,w1H1(Ω1)).

The construction of d(;) fulfills two objectives. First, for a function (𝒖,u2)U with 𝒖=(u1,𝝈,u^,σ^) and continuity u^=u2 on Γ, there holds d(𝒖,u2;𝒘,w2)=σ^,w2Γ for any (𝒘,w2)U. Therefore, this bilinear form simply replaces the duality σ^,wΓ in (2.1). Second, by selecting (𝒘,w2)=(𝒖,u2)=(u1,𝝈,u^,σ^,u2)U, the bilinear form reduces to

d(𝒖,u2;𝒖,u2)=σ^,u^Γ+12𝜷𝐧Ω1u^,u^Γ-12𝜷𝐧Ω1u2,u2Γ.

The first term on the right-hand side is needed for consistency. The last two terms are the ones that generate coercivity of the bilinear forms c1(,) and c2(,); cf. Lemma 3.1 below.

The final combined ultra-weak primal formulation of (2.1) then reads

(2.9)

(𝒖,u2)=(u1,𝝈,u^,σ^,u2)U such that
(2.9a)b(𝒖,𝒗)=L1(𝒗)for all 𝒗V,
(2.9b)c2(u2,w2)+d(𝒖,u2;𝒘,w2)=L2(w2)for all (𝒘,w2)U.

Note that, although in (2.9b) we test with functions (𝒘,w2) from U, only w2 and the restrictions to Γ of the w^ and χ^-components of 𝒘 are being used there. Therefore, (2.9a) represents problem (2.1) restricted to Ω1 without interface condition on Γ, whereas (2.9b) represents the problem on Ω2 and relates the Cauchy data u^|Γ and σ^|Γ.

For reference, we explicitly specify the strong form of (2.9a):

(2.10)𝒖:=(u1,𝝈,u^,σ^)U1 such that B𝒖=L1.

By following [12], one can show that (2.9) is equivalent to (2.1), so that, in particular, (2.9) has a unique solution; see also Remark 2.3 below. However, since we will use different strategies for solving (2.9a) and (2.9b), we need a slightly different representation. To this end, we define the trial-to-test operator Θ:U1V by

Θ𝒖,𝒗V=b(𝒖,𝒗)for all 𝒗V.

Here, ,V denotes the canonical inner product in V. Note that Θ=-1B with the Riesz operator :VV. Since B is defined on U1 without boundary condition along Γ, it has a non-trivial kernel, and so does Θ. Still, Θ:U1V is surjective. Therefore, denoting by Θκ:=κΘ the scaled trial-to-test operator (for κ>0 to be chosen) yields the following equivalent formulation: For given κ>0, find (𝒖,u2)U such that

(2.11)a(𝒖,u2;𝒘,w2)=L(𝒘,w2)for all (𝒘,w2)U,

with

a(𝒖,u2;𝒘,w2):=b(𝒖,Θκ𝒘)+c2(u2,w2)+d(𝒖,u2;𝒘,w2)

and

L(𝒘,w2):=L1(Θκ𝒘)+L2(w2).

One of our main results is the following theorem.

Theorem 2.2.

Let κ>0 be sufficiently large. Then the variational formulation (2.11) is well posed, and equivalent to problem (2.1) in the following sense: If uH01(Ω) solves (2.1), then (𝐮,u2)=(u1,𝛔,u^,σ^,u2), with ui:=u|Ωi (i=1,2) and 𝛔,u^,σ^ defined by (2.2), satisfies (𝐮,u2)U and solves (2.11).

Vice versa, if (𝐮,u2)=(u1,𝛔,u^,σ^,u2)U solves (2.11), then u defined by u|Ωi:=ui (i=1,2) satisfies uH01(Ω) and solves (2.1).

Furthermore, the bilinear form a(,) is U-coercive, i.e.,

(2.12)a(𝒖,u2;𝒖,u2)(𝒖,u2)U2for all (𝒖,u2)U.

Proof.

By the assumptions on Ω, f, α¯, 𝜷, and γ, problem (2.1) is uniquely solvable. Furthermore, by the derivation of (2.11), if uH01(Ω) solves (2.1), then (𝒖,u2) as specified in the statement solves (2.11). This can be seen by integrating by parts and noting that d(𝒖,u2;𝒘,w2)=σ^,w2Γ since u^|Γ=u2|Γ; cf. (2.7).

The coercivity of a(,) will be shown in Section 3.1 under the assumption that κ>0 is large enough. It is also straightforward to show that this bilinear form is bounded on U×U, as is the linear functional L on U. The Lax–Milgram lemma proves the well-posedness of (2.11) for κ>0 being large enough. ∎

Remark 2.3.

Theorem 2.2 implies the well-posedness of formulation (2.9), and the equivalence of (2.9) with problem (2.1), in the same form as indicated in the theorem. It is also equivalent to formulation (2.11) when κ>0 is large enough. Indeed, the solution u of (2.1) defines a solution (𝒖,u2) of (2.9). Furthermore, any solution (𝒖,u2)U of (2.9) also solves (2.11), which, for large κ, is well posed by Theorem 2.2. The well-posedness of (2.9) and its equivalence with (2.1) and (2.11) thus follow (under the condition of κ being large). Since (2.1) and (2.9) are independent of κ, they are well posed and equivalent without further restriction.

As previously mentioned, the continuity constraint u^|Γ=u2|Γ can also be imposed strongly. In this case, the solution space is

U0:={(u1,𝝈,u^,σ^;u2)U;u^|Γ=u2|Γ}

and the coupling bilinear form reduces to

d0(𝒖,w2):=d(𝒖,u2;𝒘,w2)=σ^,w2Γfor all (𝒖,u2)=(u1,𝝈,u^,σ^,u2),(𝒘,w2)U0.

The variational formulation becomes the following: For given κ>0, find (𝒖,u2)U0 such that

(2.13)a0(𝒖,u2;𝒘,w2)=L(𝒘,w2)for all (𝒘,w2)U0,

with

(2.14)a0(𝒖,u2;𝒘,w2):=b(𝒖,Θκ𝒘)+c2(u2,w2)+d0(𝒖;w2)

and

L(𝒘,w2):=L1(Θκ𝒘)+L2(w2).

Analogously to Theorem 2.2, one obtains the well-posedness of (2.13) and coercivity of a0(,).

Corollary 2.4.

Let κ>0 be sufficiently large. Then the variational formulation (2.13) is well posed, and equivalent to problem (2.1) in the following sense: If uH01(Ω) solves (2.1), then (𝐮,u2)=(u1,𝛔,u^,σ^,u2), with ui:=u|Ωi (i=1,2) and 𝛔,u^,σ^ defined by (2.2), satisfies (𝐮,u2)U0 and solves (2.13).

Vice versa, if (𝐮,u2)=(u1,𝛔,u^,σ^,u2)U0 solves (2.13), then u defined by u|Ωi:=ui (i=1,2) satisfies uH01(Ω) and solves (2.1).

Furthermore, for sufficiently large κ>0, the bilinear form a0(,) is U0-coercive, i.e.,

a(𝒖,u2;𝒖,u2)(𝒖,u2)U2for all (𝒖,u2)U0.

2.2 Combined DPG-FEM Discretization

The coupled DPG-FEM method consists in solving (2.11) within finite-dimensional subspaces UhpU. The indices h and p indicate that this can be piecewise polynomial, conforming spaces of certain polynomial degrees. Specifically, the components of Uhp that belong to the unknowns u1,𝝈,u^,σ^ will be piecewise polynomial with respect to the mesh 𝒯1 and its skeleton 𝒮. On the other hand, the component of Uhp that approximates u2 is piecewise polynomial with respect to a mesh 𝒯2 in Ω2. In the current form, we do not need compatibility of the meshes 𝒯1, 𝒯2 along Γ. The discrete scheme then reads: For given κ>0, find (𝒖hp,u2,hp)Uhp such that

(2.15)a(𝒖hp,u2,hp;𝒘,w2)=L(𝒘,w2)for all (𝒘,w2)Uhp.

Note that this formulation includes the use of optimal test functions for the discretization in Ω1; cf. (2.9a) and the corresponding terms in (2.11) with trial-to-test operator Θκ. On the other hand, the part of the problem that belongs to Ω2 is solved by standard finite elements; cf. the corresponding relation (2.9b).

Our second main result is the following theorem.

Theorem 2.5.

If κ>0 is sufficiently large, then the scheme (2.15) is uniquely solvable and converges quasi-optimally, i.e.,

𝒖-𝒖hpU1+u2-u2,hpH1(Ω2)inf{𝒖-𝒘U1+u2-w2H1(Ω2);(𝒘,w2)Uhp}.

Proof.

The statement is a direct implication of the U-coercivity of a(,) for large κ by Theorem 2.2, the Lax–Milgram lemma, and Cea’s estimate. ∎

Remark 2.6.

We note that also the discrete scheme can be changed to impose strongly the continuity of the approximations of u^ and u2 across Γ. This only requires compatibility of the meshes 𝒯1 and 𝒯2 along the interface, conforming subspaces UhpU0, and replacing the bilinear form a(;) in (2.15) by the bilinear form a0(;); cf. (2.14). The quasi-optimal error estimate from Theorem 2.5 then holds analogously.

3 Technical Details and Proof of Coercivity

We start with recalling the H01(Ω)-coercivity of the full differential operator A. This transforms into the following properties of the bilinear forms c2, c1; cf. (2.6), (2.8).

Lemma 3.1.

The bilinear forms c1(,) and c2(,) satisfy

ci(ui,ui)+12𝜷𝐧Ωiui,uiΓuiH1(Ωi)2

for all uiHD1(Ωi) (i=1,2).

Proof.

Noting that

(𝜷u,u)Ωi=-12((div𝜷)u,u)Ωi+12𝜷𝐧Ωiu,uΩifor all uH1(Ωi),i=1,2,

we obtain

ci(ui,ui)=(α¯ui-𝜷ui,ui)Ωi+(γui,ui)Ωi
=(α¯ui,ui)Ωi+((12div𝜷+γ)ui,ui)Ωi-12𝜷𝐧Ωiui,uiΓ

for uiHD1(Ωi) (i=1,2). The coercivity property then follows with the positivity of the symmetric part of α¯ and by using either the Poincaré–Friedrichs inequality and 12div𝜷+γ0 in Ωi (if meas(Γi)0) or Assumption 2.1, i.e., 12div𝜷+γβi>0 in Ωi (i=1,2). ∎

We continue with some properties of the operator B, cf. (2.10), when restricted to the space incorporating homogeneous Dirichlet boundary conditions on the whole of Ω1, that is,

(3.1)B:U1,0:=L2(Ω1)×𝐋2(Ω1)×H001/2(𝒮)×H-1/2(𝒮)V.

Lemma 3.2.

The operator B:U1,0V is an isomorphism with

B(U1,0,V)𝑎𝑛𝑑B-1(V,U1,0)

bounded independently of the mesh T1.

Proof.

This is a particular case of the different variational formulations studied in [5, Example 3.7]. More generally, in [5], Carstensen, Demkowicz and Gopalakrishnan proved that by “breaking” a continuous variational formulation of a well-posed problem (by introducing broken test spaces) and using canonical trace norms, this does not alter the well-posedness of the formulation. ∎

Let us introduce the trace space H001/2(Γ):=H01(Ω)|Γ with canonical norm. To simplify the presentation of some technical details we will need the following trace operator:

trΓ:U1H001/2(Γ),trΓ(u,𝝈,u^,σ^):=u^|Γ.

The boundedness of this operator is immediate, and analogous to the case of the Laplacian on a single domain considered in [14, Lemma 4].

Lemma 3.3.

The operator trΓ is bounded with bound independent of T1.

In the following, we identify the kernel of B when acting on the full space U1. Let us recall that A is the operator of our problem (2.1). For given φH001/2(Γ), we define its A-harmonic extension (φ):=(u1,𝝈,u^,σ^)U1 by

(3.2){u1HD1(Ω1)such thatAu1=0 in Ω1,u1=φ on Γ,𝝈=α¯u1-𝜷u1,u^=u1 on 𝒮,σ^=𝝈𝐧T on T for all T𝒯1.

Lemma 3.4.

The operator E:H001/2(Γ)U1 is bounded with bound independent of T1. Moreover, E is a right-inverse of trΓ, and the image of E is the kernel of B, that is, kerB=EH001/2(Γ).

Proof.

These statements can be proved analogously to the case of the Laplacian; cf. [14, Lemmas 11 and 12]. ∎

Of course, we also need continuity of the bilinear forms b(,), c2(,), and d(,). This is straightforward to show and has already been used in the initial part of the proof of Theorem 2.2. We only give the statement.

Lemma 3.5.

The bilinear forms b:U1×VR, c2:H1(Ω2)×H1(Ω2)R, and d:U×UR are bounded with bounds independent of the mesh T1.

3.1 Proof of Coercivity, Statement (2.12) in Theorem 2.2

We are now ready to prove the U-coercivity of the bilinear form a(,); cf. (2.11). We adapt the procedure from [14] to our situation.

Let (𝒖,u2)=(u1,𝝈,u^,σ^,u2)U be given. We start with the simple estimate

(𝒖,u2)U𝒖U1+u2H1(Ω2)
(3.3)𝒖-trΓ(𝒖)U1+trΓ(𝒖)U1+u2H1(Ω2).

By Lemma 3.4, the u^-component of 𝒖-trΓ(𝒖) has zero trace on Ω1, i.e., 𝒖-trΓ(𝒖)U1,0; cf. (3.1). By combining Lemmas 3.2 and 3.4, this gives

(3.4)𝒖-trΓ(𝒖)U1𝒖-trΓ(𝒖)U1,0B𝒖V=b(𝒖,Θ𝒖)1/2.

The last identity is due to the well-known relations of the trial-to-test operator Θ,

B𝒖V=sup𝒗Vb(𝒖,𝒗)𝒗V=b(𝒖,Θ𝒖)Θ𝒖V,Θ𝒖V=-1B𝒖V=B𝒖V.

According to Lemma 3.4, the operator is bounded,

(3.5)trΓ(𝒖)U1u^H001/2(Γ).

A combination of (3.1), (3.4), and (3.5) then gives

(3.6)(𝒖,u2)U2b(𝒖,Θ𝒖)+u^H001/2(Γ)2+u2H1(Ω2)2.

We continue by considering 𝒖e:=(u1e,𝝈e,u^e,σ^e):=trΓ(𝒖)=u^|Γ. In particular, there holds

u^H001/2(Γ)u1eH1(Ω1).

Noting that, cf. (2.7),

d(𝒖e,u2;𝒖e,u2)=σ^e,u2Γ+σ^e,u^-u2Γ+12𝜷𝐧Ω1(u^-u2),u^+u2Γ
(3.7)=σ^e,u^Γ+12𝜷𝐧Ω1u^,u^Γ-12𝜷𝐧Ω1u2,u2Γ,

an application of Lemma 3.1 gives

u^H001/2(Γ)2+u2H1(Ω2)2c1(u1e,u1e)+c2(u2,u2)+12𝜷𝐧Ω1u1e,u1eΓ+12𝜷𝐧Ω2u2,u2Γ
=c1(u1e,u1e)+c2(u2,u2)+d(𝒖e,u2;𝒖e,u2)-σ^e,u^Γ.

Relation (3.1) can also be written as

d(𝒖e,u2;𝒖e,u2)=d(𝒖,u2;𝒖,u2)+σ^e-σ^,u^Γ,

so that the previous bound becomes

u^H001/2(Γ)2+u2H1(Ω2)2c1(u1e,u1e)+c2(u2,u2)+d(𝒖,u2;𝒖,u2)-σ^,u^Γ.

Now, by recalling the definitions of c1(,) (see (2.8)) and the extension operator (see (3.2)), integration by parts yields the expected relation c1(u1e,u1e)=σ^e,u^Γ. Therefore, continuing the estimate, we obtain

(3.8)u^H001/2(Γ)2+u2H1(Ω2)2c2(u2,u2)+d(𝒖,u2;𝒖,u2)+σ^e-σ^,u^Γ.

The last term in (3.8) can be bounded by duality, the continuity of H-1/2(𝒮)(σ^e-σ^)(σ^e-σ^)|ΓH-1/2(Γ), and relation (3.4). This gives

σ^e-σ^,u^Γ𝒖-trΓ(𝒖)U1u^H001/2(Γ)b(𝒖,Θ𝒖)1/2(𝒖,u2)U.

Combining this bound with (3.6) and (3.8) and applying Young’s inequality, we find that

(𝒖,u2)U2κb(𝒖,Θ𝒖)+c2(u2,u2)+d(𝒖,u2;𝒖,u2)
=b(𝒖,Θκ𝒖)+c2(u2,u2)+d(𝒖,u2;𝒖,u2)

for a sufficiently large constant κ>0. This proves the stated coercivity of a(,).

4 Numerical Experiments

In this section, we report on two numerical experiments. In both of them we choose d=2 and, starting from a manufactured solution, we compute the right-hand side function f. The solution of the second experiment does not satisfy the homogeneous Dirichlet boundary condition. In this case, we use a standard approach and extend the inhomogeneous Dirichlet datum into the domain and then shift the resulting terms to the right-hand side. As discrete trial space we use

Uhp:=P0(𝒯1)×[P0(𝒯1)]2×SD1(𝒮)×P0(𝒮)×SD1(𝒯2),

where 𝒯1 and 𝒮 are a mesh and its skeleton in Ω1 and 𝒯2 is a mesh in Ω2. Throughout, we use meshes 𝒯1 and 𝒯2 which are compatible on the interface Γ (although this is not necessary in our analysis). In the definition of Uhp, the term Pk(𝒯1) denotes the space of 𝒯1-piecewise polynomials of degree k, the term P0(𝒮) denotes the space of piecewise constant functions on 𝒮, and the term SD1(𝒮)HD1/2(𝒮) denotes the space of piecewise affine and continuous functions on 𝒮 which vanish on Γ1. The space SD1(𝒯2)HD1(Ω2) is the space of piecewise affine, globally continuous functions on 𝒯2 which vanish on Γ2. The trial-to-test operator Θκ=κ-1B with :VV being the Riesz operator is approximated using the discrete Riesz operator hp:VhpVhp with a finite dimensional space VhpV, which we choose to be

Vhp:=P2(𝒯1)×[P2(𝒯1)]2.

The resulting method is called practical DPG method, and was analyzed in [16]. In that work, it was shown that the additional discretization error of using Vhp instead of V does not degrade the convergence order. Throughout, we use κ=1 and do not encounter difficulties with this choice. Note that if (u1,𝝈,u^,σ^,u2) denotes the exact solution of (2.11) and the discrete solution of (2.15) is (u1,hp,𝝈hp,u^hp,σ^hp,u2,hp)Uhp, then by definition of the norm H1/2(𝒮) there holds

u^-u^hpH1/2(𝒮)u-I𝒯1u^hpH1(Ω1)=:err(u^),

where I𝒯1u^hpSD1(𝒯1) is the nodal interpolant of u^hp with (I𝒯1u^hp)|𝒮=u^hp. Likewise,

σ^-σ^hpH-1/2(𝒮)𝝈-𝐈𝒯1σ^hp𝐇(div,Ω1)=:err(σ^),

where 𝐈𝒯1σ^hp𝒯0(𝒯1) is the lowest-order Raviart–Thomas interpolant of σ^hp, i.e.,

(𝐧T𝐈𝒯1σ^hp)|T=σ^hp|T

for any T𝒯1. Furthermore, we plot the errors

err(u1):=u1-u1,hpL2(Ω1),
err(𝝈):=𝝈-𝝈hp𝐋2(Ω1),
err(u2):=u2-u2,hpH1(Ω2),

as well as the so-called energy error of the DPG part

err(𝒖):=sup𝒗Vb(𝒖-𝒖hp,𝒗)𝒗V=Θκ(𝒖-𝒖hp)V;

cf. (3.4). In both experiments, we use a sequence of meshes resulting from uniform mesh refinements. The quasi-optimality result of Theorem 2.5 and well-known approximation results then show that

𝒖-𝒖hpU1+u2-u2,hpH1(Ω2)=𝒪(h)=𝒪(N-1/2).

Here, N denotes the overall number of degrees of freedom of Uhp. Hence, err()=𝒪(N-1/2) for all of the errors defined above.

4.1 Experiment 1

We choose Ω1:=(0,1)×(0,1) and Ω2:=(1,2)×(0,1) and use the exact solution

u(x,y):=x(2-x)y(1-y).

The remaining parameters in equation (2.1) are chosen as α¯=id, 𝜷=(xy,1), and γ=1-sin(πx). In Figure 2, we plot the errors versus the degrees of freedom on a double logarithmic scale. As expected, all the errors behave like 𝒪(N-1/2), which is plotted in black without markers. In Figure 3, we plot the error u^-u2 on the coupling boundary Γ for the case with mesh width 1/32.

Figure 2 Error plots for Experiment 1. The black line without markers denotes 𝒪⁢(N-1/2){\mathcal{O}(N^{-1/2})}, and N is the total number of degrees of freedom.
Figure 2

Error plots for Experiment 1. The black line without markers denotes 𝒪(N-1/2), and N is the total number of degrees of freedom.

Figure 3 The jump u^-u2{\widehat{u}-u_{2}} on the interface Γ in Experiment 1 with mesh width 1/32{1/32}.
Figure 3

The jump u^-u2 on the interface Γ in Experiment 1 with mesh width 1/32.

4.2 Experiment 2

Figure 4 The exact solution from Experiment 2.
Figure 4

The exact solution from Experiment 2.

Figure 5 Error plots for Experiment 2. The black line without markers denotes
𝒪⁢(N-1/2){\mathcal{O}(N^{-1/2})}, and N is the total number of degrees of freedom.
Figure 5

Error plots for Experiment 2. The black line without markers denotes 𝒪(N-1/2), and N is the total number of degrees of freedom.

Figure 6 The jump u^-u2{\widehat{u}-u_{2}} on the interface Γ in Experiment 2 with mesh width 1/32{1/32}.
Figure 6

The jump u^-u2 on the interface Γ in Experiment 2 with mesh width 1/32.

We choose Ω1:=(0.2,0.7)×(0.2,1.2) and Ω2:=(0.7,1.2)×(0.2,1.2) and use the exact solution

u(x,y):=arctan(1-|(x,y)|ε).

The remaining parameters in equation 2.1 are chosen as α¯=εid, 𝜷=exp(x)(siny,cosy), γ=0, and ε=0.05. The exact solution u has a curved layer of moderate width inside Ω; see Figure 4. In Figure 5, we plot the errors versus the degrees of freedom on a double logarithmic scale. Again, as expected, we obtain the convergence order 𝒪(N-1/2). In Figure 6, we plot the error u^-u2 on the coupling boundary Γ, again for mesh width 1/32. Note that the layer of the exact solution cuts through Γ, and this is reflected in Figure 6.

Award Identifier / Grant number: 1150056

Award Identifier / Grant number: 1170672

Award Identifier / Grant number: ACT1118

Award Identifier / Grant number: PFB-03

Funding statement: We gratefully acknowledge financial support by CONICYT through FONDECYT projects 1150056 and 1170672, Anillo ACT1118 (ANANUM), and BASAL project CMM, Universidad de Chile, Chile.

References

[1] J. W. Barrett and K. W. Morton, Approximate symmetrization and Petrov–Galerkin methods for diffusion-convection problems, Comput. Methods Appl. Mech. Engrg. 45 (1984), no. 1–3, 97–122. 10.1016/0045-7825(84)90152-XSearch in Google Scholar

[2] D. Broersen and R. Stevenson, A robust Petrov–Galerkin discretisation of convection-diffusion equations, Comput. Math. Appl. 68 (2014), no. 11, 1605–1618. 10.1016/j.camwa.2014.06.019Search in Google Scholar

[3] D. Broersen and R. P. Stevenson, A Petrov–Galerkin discretization with optimal test space of a mild-weak formulation of convection-diffusion equations in mixed form, IMA J. Numer. Anal. 35 (2015), no. 1, 39–73. 10.1093/imanum/dru003Search in Google Scholar

[4] C. Carstensen, L. Demkowicz and J. Gopalakrishnan, A posteriori error control for DPG methods, SIAM J. Numer. Anal. 52 (2014), no. 3, 1335–1353. 10.1137/130924913Search in Google Scholar

[5] C. Carstensen, L. Demkowicz and J. Gopalakrishnan, Breaking spaces and forms for the DPG method and applications including Maxwell equations, Comput. Math. Appl. 72 (2016), no. 3, 494–522. 10.1016/j.camwa.2016.05.004Search in Google Scholar

[6] P. Causin and R. Sacco, A discontinuous Petrov–Galerkin method with Lagrangian multipliers for second order elliptic problems, SIAM J. Numer. Anal. 43 (2005), no. 1, 280–302. 10.1137/S0036142903427871Search in Google Scholar

[7] J. Chan, N. Heuer, T. Bui-Thanh and L. Demkowicz, A robust DPG method for convection-dominated diffusion problems II: Adjoint boundary conditions and mesh-dependent test norms, Comput. Math. Appl. 67 (2014), no. 4, 771–795. 10.1016/j.camwa.2013.06.010Search in Google Scholar

[8] L. Demkowicz and J. Gopalakrishnan, A class of discontinuous Petrov–Galerkin methods. Part I: The transport equation, Comput. Methods Appl. Mech. Engrg. 199 (2010), no. 23–24, 1558–1572. 10.1016/j.cma.2010.01.003Search in Google Scholar

[9] L. Demkowicz and J. Gopalakrishnan, A class of discontinuous Petrov–Galerkin methods. II. Optimal test functions, Numer. Methods Partial Differential Equations 27 (2011), no. 1, 70–105. 10.1002/num.20640Search in Google Scholar

[10] L. Demkowicz, J. Gopalakrishnan and A. H. Niemi, A class of discontinuous Petrov-Galerkin methods. Part III: Adaptivity, Appl. Numer. Math. 62 (2012), no. 4, 396–427. 10.1016/j.apnum.2011.09.002Search in Google Scholar

[11] L. Demkowicz and N. Heuer, Robust DPG method for convection-dominated diffusion problems, SIAM J. Numer. Anal. 51 (2013), no. 5, 2514–2537. 10.1137/120862065Search in Google Scholar

[12] F. Fuentes, B. Keith, L. Demkowicz and P. Le Tallec, Coupled variational formulations of linear elasticity and the DPG methodology, J. Comput. Phys. 348 (2017), 715–731. 10.1016/j.jcp.2017.07.051Search in Google Scholar

[13] T. Führer and N. Heuer, Robust coupling of DPG and BEM for a singularly perturbed transmission problem, Comput. Math. Appl. (2016), 10.1016/j.camwa.2016.09.016. 10.1016/j.camwa.2016.09.016Search in Google Scholar

[14] T. Führer, N. Heuer and M. Karkulik, On the coupling of DPG and BEM, Math. Comp. 86 (2017), no. 307, 2261–2284. 10.1090/mcom/3170Search in Google Scholar

[15] T. Führer, N. Heuer and E. P. Stephan, On the DPG method for Signorini problems, preprint (2016), https://arxiv.org/abs/1609.00765; to appear in IMA J. Numer. Anal., DOI 10.1093/imanum/drx048. 10.1093/imanum/drx048Search in Google Scholar

[16] J. Gopalakrishnan and W. Qiu, An analysis of the practical DPG method, Math. Comp. 83 (2014), no. 286, 537–552. 10.1090/S0025-5718-2013-02721-4Search in Google Scholar

[17] N. Heuer and M. Karkulik, A robust DPG method for singularly perturbed reaction-diffusion problems, SIAM J. Numer. Anal. 55 (2017), no. 3, 1218–1242. 10.1137/15M1041304Search in Google Scholar

[18] N. V. Roberts, Camellia: A software framework for discontinuous Petrov–Galerkin methods, Comput. Math. Appl. 68 (2014), no. 11, 1581–1604. 10.1016/j.camwa.2014.08.010Search in Google Scholar

Received: 2017-04-25
Revised: 2017-08-31
Accepted: 2017-09-02
Published Online: 2017-09-30
Published in Print: 2018-10-01

© 2018 Walter de Gruyter GmbH, Berlin/Boston

Downloaded on 13.12.2024 from https://www.degruyter.com/document/doi/10.1515/cmam-2017-0041/html
Scroll to top button