A weak Galerkin finite element scheme with boundary continuity for second-order elliptic problems

# A weak Galerkin finite element scheme with boundary continuity for second-order elliptic problems

Qilong Zhai Department of Mathematics, Jilin University, Changchun China    Xiu Ye Department of Mathematics, University of Arkansas at Little Rock, Little Rock, AR 72204 (xxye@ualr.edu). The research of Ye was supported in part by National Science Foundation Grant DMS-1115097.    Ruishu Wang Department of Mathematics, Jilin University, Changchun China    Ran Zhang Department of Mathematics, Jilin University, Changchun, China (zhangran@mail.jlu.edu.cn). The research of Zhang was supported in part by China Natural National Science Foundation(11271157, 11371171, 11471141), and by the Program for New Century Excellent Talents in University of Ministry of Education of China.
###### Abstract

A new weak Galerkin (WG) finite element method for solving the second-order elliptic problems on polygonal meshes by using polynomials of boundary continuity is introduced and analyzed. The WG method is utilizing weak functions and their weak derivatives which can be approximated by polynomials in different combination of polynomial spaces. Different combination gives rise to different weak Galerkin finite element methods, which makes WG methods highly flexible and efficient in practical computation. This paper explores the possibility of certain combination of polynomial spaces that minimize the degree of freedom in the numerical scheme, yet without losing the accuracy of the numerical approximation. Error estimates of optimal order are established for the corresponding WG approximations in both a discrete norm and the standard norm. In addition, the paper also presents some numerical experiments to demonstrate the power of the WG method. The numerical results show a great promise of the robustness, reliability, flexibility and accuracy of the WG method.

Key words. weak Galerkin finite element methods, weak gradient, second-order elliptic equation, polygonal meshes.

AMS subject classifications. Primary, 65N30, 65N15, 65N12, 74N20; Secondary, 35B45, 35J50, 35J35

## 1 Introduction

Nowadays, finite element methods (FEMs) are widely used in almost every field of engineering and industrial analysis. Since R. Courant  formulated the essence of what is now called a finite element in 1943, this method has been getting more and more attractive with the development of computers and is now recognized as one of the most versatile and powerful methods for approximating the solutions of boundary-value problems, especially for problems over complicated domains. Among the different finite elements methods, the conforming finite element method  with continuous, piecewise polynomial approximating spaces, has long been employed to approximate solutions for partial differential equations. Within the past few decades, however, a number of researchers have investigated Galerkin methods based on fully discontinuous approximating spaces, such as the discontinuous Galerkin (DG) methods [5, 1, 2, 3], the Hybridized discontinuous Galerkin (HDG) methods [6, 8], the weak Galerkin (WG) methods [19, 20], etc.

The WG method was first introduced in 2012  for the second order elliptic problem and further developed with other applications, such as the Stokes, Helmholtz, Maxwell, biharmonic [13, 12, 9, 10, 17, 16, 19, 20, 23], etc. Its central idea is that the shape function in the interior of each element is simply polynomial. The weak finite element function could be totally discontinuous across elements. The continuity is compensated by the stabilizer through a suitable boundary integral defined on the boundary of elements. That is, we know more information on the shape function by sacrificing the continuity.

Comparing with the conforming FEMs, there are a lot of advantages for the discontinuous methods, for instance, high-order accuracy, multiphysics capability, the finite element partition can be of polygon or polyhedral type, the weak finite element space is easy to construct with any approximation requirement and suit for any given stability requirement. All of this, however, comes at a price: most notably through an increase in the total degrees of freedom (d.o.f.) as a direct result of the decoupling of the elements. For linear elements, this yields a doubling in the total number of degrees of freedom compared to the continuous FEM. For WG finite element methods, in order to reduce the d.o.f., people have made some efforts: in [14, 22], the possibility of optimal combination of polynomial spaces was explored to minimize the number of unknowns in the numerical scheme, yet without compromising the accuracy of the numerical approximation; in , some hybridization of finite element methods has been introduced by utilizing the Lagrange multiplier. The distinctive feature of the methods in this framework is that the only globally coupled degrees of freedom are those of an approximation of the solution defined only on the boundaries of the elements, then the global unknowns are the numerical traces of the field variables. Thus, one can reduce the number of the globally coupled degrees of freedom of WG methods.

The goal of this paper is to explore the possibility of certain combination of the polynomial spaces to reduce the number of unknowns without compromising the rate of convergence, which utilizes the boundary continuity for the WG finite element spaces. This idea is motivated by Chen  that was presented in ICIAM 2015. Combining with the Schur complement technique, we further eliminate the interior unknowns and produce a much reduced system of linear equations involving only the unknowns representing the interface variables. In fact, if we take the triangle mesh, the d.o.f. of the new scheme is the same with that of conforming FEMs.

Next, we introduce this WG methods. For the sake of simplicity and easy presentation of the main ideas, we restrict ourselves to the following model problem

 −∇⋅(a∇u) = f in Ω, \hb@xt@.01(1.1) u = 0 on ∂Ω, \hb@xt@.01(1.2)

where is an open bounded polygonal domain in . We assume that are given, sufficiently smooth functions, and is a symmetric matrix-valued function. Suppose there exists a positive number such that

 ξtaξ≥λξtξ,∀ξ∈R2,

where is a column vector and is the transpose of .

The paper is organized as follows. In Section 2, we present some standard notations in Sobolev spaces and preliminaries. A weakly-defined differential operator is also introduced. The weak Galerkin finite element scheme is developed and some properties for the error analysis are discussed in Section 3. In Section 4, we shall derive an error equation for the WG approximations. Optimal-order error estimates of and for the WG finite element approximations are also derived in this Section. The equivalence of WG formulation and its Schur complement formulation is proved in Section 5. In Section 6, numerical experiments are conducted. Finally, we present some technical estimates for quantities related to the local projections into various finite element spaces and some approximation properties which are useful in the convergence analysis in ¡°Appendix¡±.

## 2 Notations and preliminaries

In this section, we shall introduce some notations used in this paper.

We use the standard Soblev space notations. For an open set , and stand for the norm and inner-product, namely. We shall drop the subscripts when and .

Let be a partition of the domain consisting of polygons in two dimension or polyhedra in three dimension satisfying a set of conditions specified in . Denote by the set of all edges or at faces in . For every element , we denote by its diameter and mesh size for .

For a given partition , denotes a piecewise polynomial on whose degree is no more than on each . Similarly, denotes a piecewise polynomial on whose degree is no more than on each .

Now, we define the weak finite element space as follows:

 Vh={(v0,vb):v0|T∈Pk(T),vb|e∈Pk(e),vb is continuous on Eh,vb=0 on ∂Ω}.

It should be noticed that is single-valued on each edge, and is continuous on , which means that share the same value on each node of .

Similar to the definition in , we can define the following weak gradient operator on .

###### Definition 2.1

For any , define the discrete weak gradient satisfying

 (∇wvh,q)T=−(v0,∇⋅q)T+⟨vb,q⋅n⟩∂T,∀q∈[Pk−1(T)]2, \hb@xt@.01(2.1)

where n is the outward unit normal vector along , stands for the -inner product in , and is the inner product in

## 3 A weak Galerkin finite element scheme

In this section, we shall propose a WG scheme for the second-order elliptic problem, and verify the wellposedness of the scheme.

For any , define the following bilinear forms

 s(vh,wh) = ∑T∈Thh−1T⟨v0−vb,w0−wb⟩∂T, aw(vh,wh) = (a∇wvh,∇wwh)+s(vh,wh).
###### Weak Galerkin Algorithm 1

The weak Galerkin numerical solution of problem (LABEL:problem-eq1) -(LABEL:problem-eq2) can be obtained by seeking such that

 aw(uh,vh)=(f,v0)∀vh∈Vh. \hb@xt@.01(3.1)

The following semi-norm can be inducted from directly

 |||vh|||2=aw(vh,vh),∀vh∈Vh.

We claim that defines a norm on indeed.

Notice that when , we have on each and on each , and it follows that ,

 0 = (∇wvh,∇v0)T = −(v0,∇⋅∇v0)T+⟨vb,∇v0⋅n⟩∂T = (∇v0,∇v0)T−⟨v0−vb,∇v0⋅n⟩∂T = ∥∇v0∥2T,

which implies that is constant on . With the condition that on each and on , we can conclude that , and thus defines a norm on . From this property we can arrive the wellposedness of WG scheme (LABEL:wg-scheme) directly.

###### Theorem 3.1

The WG scheme (LABEL:wg-scheme) has a unique solution.

Next, we shall discuss some properties for the error analysis.

Define the following projections operators:

 Q0:L2(T)→Pk(T),∀T∈Th, Qb:L2(e)→Pk(e),∀e∈Eh, Qh:[L2(T)]2→[Pk−1(T)]2,∀T∈Th.

It is known that on a one-dimension edge , for any smooth function we can get a polynomial interpolation of degree with different points. Suppose the interpolation points include the two endpoints of the edge, and we have the interpolation operator .

Combining and , we can define

 ~Qh={Q0,Ib}:H1(Ω)→Vh.

In the rest of this paper, we denote the local projection of . Obviously, is symmetric and bounded. With these projections, it arrives at the following community property.

###### Lemma 3.2

For any , ,

 (∇w~Qhw,q) = (Qh∇w,q)+∑T∈Th⟨Ibw−w,q⋅n⟩∂T.

Proof. For any , it follows the definition (LABEL:wgradient) and the integration by parts that

 (∇w~Qhw,q)T = −(Q0w,∇⋅q)T+⟨Ibw,q⋅n⟩∂T = −(Q0w,∇⋅q)T+⟨w,q⋅n⟩∂T+⟨Ibw−w,q⋅n⟩∂T = (∇w,q)T+⟨Ibw−w,q⋅n⟩∂T = (Qh∇w,q)T+⟨Ibw−w,q⋅n⟩∂T.

Summing over all and the proof is completed.

## 4 Error analysis

In this section, we shall present the and errors of the weak Galerkin numerical scheme (LABEL:wg-scheme) with optimal orders.

Denote the error , and we shall show the error equation for .

###### Lemma 4.1

Suppose is the solution of (LABEL:problem-eq1)- (LABEL:problem-eq2), and is the WG numerical solution of (LABEL:wg-scheme). Then for any ,

 aw(~Qhu−uh,vh)=l(u,vh),

where

 l(u,vh) = (a∇w~Qhu−Qh(a∇u),∇wvh) −∑T∈Th⟨(Qh(a∇u)−a∇u)⋅n,v0−vb⟩∂T+s(~Qhu,vh).

Proof. For any , it follows the integration by parts that

 aw(uh,vh) = (f,v0) \hb@xt@.01(4.1) = (a∇u,∇v0)−∑T∈Th⟨a∇u⋅n,v0⟩∂T = (a∇u,∇v0)−∑T∈Th⟨a∇u⋅n,v0−vb⟩∂T.

From the definition (LABEL:wgradient), we can obtain

 (a∇u,∇v0) \hb@xt@.01(4.2) = (Qh(a∇u),∇v0) = (Qh(a∇u),∇wvh)+∑T∈Th⟨Qh(a∇u)⋅n,v0−vb⟩∂T

Substituting (LABEL:err-est12) into (LABEL:err-est3), we get

 aw(uh,vh) = (Qh(a∇u),∇wvh)+∑T∈Th⟨(Qh(a∇u)−a∇u)⋅n,v0−vb⟩∂T. \hb@xt@.01(4.3)

Notice that

 aw(~Qhu,vh) = (a∇w~Qhu,∇wvh)+s(~Qhu,vh). \hb@xt@.01(4.4)

Subtracting (LABEL:err-est13) from (LABEL:err-est14), we have

 aw(~Qhu−uh,vh) = (a∇w~Qhu−Qh(a∇u),∇wvh) −∑T∈Th⟨(Qh(a∇u)−a∇u)⋅n,v0−vb⟩∂T+s(~Qhu,vh),

which completes the proof.

With the error equation and the technique tools in Appendix, we can arrive at the estimate for the error.

###### Theorem 4.2

Suppose is the solution of (LABEL:problem-eq1)- (LABEL:problem-eq2), and is the WG numerical solution of (LABEL:wg-scheme). The following estimate holds true,

 |||~Qhu−uh|||≤Chk∥u∥k+1.

Proof. Take in Lemma LABEL:error-eqn, then from Lemma LABEL:err-est4 we can obtain

 |||~Qhu−uh|||2 = l(u,~Qhu−uh) ≤ Chk|||~Qhu−uh|||,

which completes the proof.

Next, we shall show the optimal error order using the well-known Nistche’s argument. Consider the following dual problem: find satisfying

 −∇⋅(a∇φ)=e0. \hb@xt@.01(4.5)

Assume the dual problem (LABEL:dual-eq) has regularity, i.e.

 ∥φ∥2≤C∥e0∥. \hb@xt@.01(4.6)
###### Theorem 4.3

Suppose is the solution of (LABEL:problem-eq1)- (LABEL:problem-eq2), is the WG numerical solution of (LABEL:wg-scheme), and the dual problem (LABEL:dual-eq) satisfies the regularity assumption (LABEL:regularity-assump). The following estimate holds true,

 ∥Q0u−u0∥≤Chk+1∥u∥k+1.

Proof. Testing equation (LABEL:dual-eq) by and we can obtain

 ∥e0∥2 \hb@xt@.01(4.7) = (a∇φ,∇e0)−∑T∈Th⟨∇φ⋅n,e0⟩∂T = (Qh(a∇φ),∇e0)−∑T∈Th⟨∇φ⋅n,e0−eb⟩∂T = (Qh(a∇φ),∇weh)+∑T∈Th⟨(Qh(a∇φ)−a∇φ)⋅n,e0−eb⟩∂T.

Similar to the derivative of formula (LABEL:err-est14), we have

 aw(~Qhφ,eh) = (a∇w~Qhφ,∇weh)+s(~Qhφ,eh). \hb@xt@.01(4.8)

Subtracting (LABEL:err-est6) from (LABEL:err-est5) yields

 ∥e0∥2=aw(~Qhφ,eh)−l(φ,eh).

By substituting in Lemma LABEL:error-eqn we get

 aw(eh,~Qhφ)=l(u,~Qhφ),

which implies

 ∥e0∥2=l(u,~Qhφ)−l(φ,eh). \hb@xt@.01(4.9)

From Lemma LABEL:err-est4 and (LABEL:regularity-assump), we derive that

 l(φ,eh) ≤ Ch|||eh|||∥φ∥2 \hb@xt@.01(4.10) ≤ Chk+1∥u∥k+1∥e0∥.

Now we just need to estimate the following three terms of . As to the first term, we split it into two formulas and estimate them separately as follows

 (a∇w~Qhu−Qh(a∇u),∇w~Qhφ) = (∇w~Qhu−∇u,a∇w~Qhφ) = (∇w~Qhu−∇u,a∇w~Qhφ−Qh(a∇w~Qhφ)) +(∇w~Qhu−∇u,Qh(a∇w~Qhφ)) = (∇w~Qhu−∇u,a∇w~Qhφ−Qh(a∇w~Qhφ)) +(∇w~Qhu−Qh∇u,Qh(a∇w~Qhφ)).

For the first formula, from Lemma LABEL:Qh-dQ and Lemma LABEL:Qh-H2 we can obtain

 |(∇w~Qhu−∇u,a∇w~Qhφ−Qh(a∇w~Qhφ))| ≤ ∥∇w~Qhu−∇u∥∥a∇w~Qhφ−Qh(a∇w~Qhφ)∥ ≤ Chk∥u∥k+1∥∇w~Qhφ∥1 ≤ Chk+1∥u∥k+1∥φ∥2.

As to the second formula, it follows Lemmas LABEL:commu-prop and LABEL:Qh-dQ that

 (∇w~Qhu−Qh∇u,Qh(a∇w~Qhφ)) = ∑T∈Th⟨Ibu−u,Qh(a∇w~Qhφ)⋅n⟩∂T = ∑T∈Th⟨Ibu−u,Qh(a∇w~Qhφ−a∇φ)⋅n⟩∂T +∑T∈Th⟨Ibu−u,(Qh(a∇φ)−a∇φ)⋅n⟩∂T ≤ Chk+1∥u∥k+1∥φ∥2.

For the second term, from Lemma LABEL:err-est1 and the trace inequality we can obtain

 ∑T∈Th⟨(Qh(a∇u)−a∇u)⋅n,Q0φ−Ibφ⟩∂T ≤ C⎛⎝∑T∈ThhT∥Qh(a∇u)−a∇u∥2∂T⎞⎠12⎛⎝∑T∈Thh−1T∥Q0φ−φ∥2∂T⎞⎠12 +C⎛⎝∑T∈Th∥Qh(a∇u)−a∇u∥2T⎞⎠12⎛⎝∑T∈Thh−1T∥φ−Ibφ∥2∂T⎞⎠12 ≤ Chk+1∥u∥k+1∥φ∥2.

Similarly, for the third term we have

 s(~Qhu,~Qhφ) = ∑T∈Thh−1T⟨Q0u−Ibu,Q0φ−Ibφ⟩∂T = ∑T∈Thh−1T⟨Q0u−u,Q0φ−φ⟩∂T+∑T∈Thh−1T⟨u−Ibu,Q0φ−φ⟩∂T +∑T∈Thh−1T⟨Q0u−u,φ−Ibφ⟩∂T+∑T∈Thh−1T⟨u−Ibu,φ−Ibφ⟩∂T ≤ Chk+1∥u∥k+1∥φ∥2,

 l(u,~Qhφ) ≤ Chk+1∥u∥k+1∥φ∥2 \hb@xt@.01(4.11) ≤ Chk+1∥u∥k+1∥e0∥.

By substituting (LABEL:err-est10) and (LABEL:err-est11) into (LABEL:err-est7), the proof is completed.

## 5 Schur Complement

In this section, we shall show the technique in practical which can eliminate the degree of freedoms of , and then only the degree of freedoms of are involved in the total linear system.

Let be the solution of the WG scheme (LABEL:wg-scheme), then satisfies the following equation

 aw(uh,vh)=(f,v0),∀vh∈Vh.

This equation can be rewritten as the following form

 aw(uh,vh) = (f,v0),∀vh={v0,0}∈Vh, \hb@xt@.01(5.1) aw(uh,vh) = 0,∀vh={0,vb}∈Vh. \hb@xt@.01(5.2)

For given on , one can solve equation (LABEL:schur-1) for the , and note

 u0=D(ub,f). \hb@xt@.01(5.3)

It should be noticed that can be calculated elementwisely. We can solve the local problem

 aw(uh,v0) = (f,v0),∀v0∈Pk(T), \hb@xt@.01(5.4)

and get

 u0=DT(ub,f). \hb@xt@.01(5.5)

Then we substitute (LABEL:schur-3) into (LABEL:schur-2) and get the following equation

 aw({D(ub,f),ub},vh)=0,∀vh={0,vb}∈Vh. \hb@xt@.01(5.6)

Thus we can conclude a complement for the WG scheme as follows.

Step 1. Solve equation (LABEL:schur-6) on each element and get .

Step 2. Solve the global system (LABEL:schur-2) for .

Step 3. Recover by on each element.

When we take the uniform triangle mesh on the unit square and set , the degree of freedom of different schemes are listed in Table LABEL:table-3.

stands for the degree of freedom in the scheme in , stands for scheme (LABEL:wg-scheme), stands for the scheme in  with Schur complement, stands for scheme (LABEL:wg-scheme) with Schur complement, and stands for the conforming Galerkin method.

## 6 Numerical experiment

In this section, we shall present some numerical results to show the efficiency and accuracy of the WG scheme.

Example 1 Consider the problem (LABEL:problem-eq1) - (LABEL:problem-eq2) on the unit square . The analytic solution is set to be , the right side and the Dirichlet boundary condition is computed to match the exact solution. In this numerical experiment, the uniform triangle mesh is employed. We can conclude from Table LABEL:table-1 that the convergence rate for and errors are of and , respectively. This is coincide with the theoretical analysis in this paper.

Example 2 Consider the problem (LABEL:problem-eq1)- (LABEL:problem-eq2) on the unit square . The analytic solution is set to be , the right-hand side side and the Dirichlet boundary condition is computed to match the exact solution. In this numerical experiment, the uniform rectangle mesh is employed. We can conclude from Table LABEL:table-4 that the convergence rate for and errors are of and , respectively. This is coincide with the theoretical results.

## 7 Appendix

In this section, we shall give some technique tools which are applied in the error estimate.

The trace inequality and the inverse inequality for the weak Galerkin method are proved in .

###### Lemma 7.1

Assume that the partition satisfies the assumptions (A1), (A2), and (A3) as specified in . Then, there exists a constant such that for any and edge/face , we have for any ,

 ∥θ∥pe≤Ch−1T(∥θ∥pT+hpT∥∇θ∥pT). \hb@xt@.01(7.1)
###### Lemma 7.2

Assume that the partition satisfies the assumptions (A1), (A2), (A3) and (A4) as specified in . Then, there exists a constant such that for any , we have for any ,

 ∥∇φ∥T≤C(n)h−1T∥φ∥T. \hb@xt@.01(7.2)

The following estimate of interpolation error plays an essential role in the error estimate.

###### Lemma 7.3

Assume that the partition satisfies the assumptions (A1), (A2), (A3) and (A4) as specified in . For any , the following estimate holds true

 ∑T∈Th∥w−Ibw∥2∂T≤Ch2k+1∥w∥2k+1.

Proof. For any edge , according to assumption (A3), there exists a triangle , whose base is identical to and height is proportional to . Then there exists a interpolation operator onto which contains all the interpolation points of . It follows the trace inequality (LABEL:trace-ineq) and the property of interpolation operator that for any ,

 ∥w−Ibw∥2e = ∥w−IPew∥2e ≤ C(h−1Pe∥w−IPew∥2Pe+hPe∥∇(w−IPew)∥2Pe) ≤ Ch2k+1Pe∥w∥2k+1,Pe ≤ Ch2k+1T∥w∥2k+1,T.

Summing over all and the proof is completed.

With Lemma LABEL:err-est1, we can describe the difference between and in the following lemma.

###### Lemma 7.4

Assume that the partition satisfies the assumptions (A1), (A2), (A3) and (A4) as specified in . For any , the following estimate holds true

 ∥∇w~Qhw−Qh∇w∥≤Chk∥w∥k+1.

Proof. It follows Lemma LABEL:commu-prop that for any ,

 (∇w~Qhw−Qh∇w,q) = ∑T∈Th⟨Ibw−w,q⋅n⟩∂T ≤ C⎛⎝∑T∈Th∥Ibw−w∥2∂T⎞⎠12⎛⎝∑T∈Th∥q∥2∂T⎞⎠12.

From the trace inequality and the inverse inequality, we can obtain

 ∑T∈Th∥q∥2∂T ≤ C∑T∈Thh−1T∥q∥2T+C∑T∈ThhT∥∇q∥2T ≤ C∑T∈Thh−1T∥q∥2T.

Taking and applying Lemma LABEL:err-est1 yield

 ∥∇w~Qhw−Qh∇w∥2≤Chk∥w∥k+1∥∇w~Qhw−Qh∇w∥,

which completes the proof.

###### Lemma 7.5

Assume that the partition satisfies the assumptions (A1), (A2), (A3) and (A4) as specified in . For any , the following estimate holds true

 ∥∇w~Qhw∥1≤C∥w∥2.

Proof. When , from definition LABEL:wgradient we know and , so that

 ∥∇w~Qhw−∇w∥1=0,∀w∈P1(T).

Then, from Bramble-Hilbert Theorem it follows that

 ∥∇w~Qhw−∇w∥1≤Ch∥w∥2,∀w∈H2(T).

From the triangle inequality we can obtain

 ∥∇w~Qhw∥1 ≤ ∥∇w~Qhw−∇w∥1+∥∇w∥1 ≤ C∥w∥2,

which completes the proof.

In order to give the error estimate, we need to estimate the remainder of the error equation in Lemma LABEL:error-eqn.

###### Lemma 7.6

Suppose and , then the following estimates hold true

 (a∇w~Qhw−Qh(a∇w),∇wvh) ≤ Chk∥w∥k+1|||vh|||, ∑T∈Th⟨(Qh(a∇w)−a∇w)⋅n,v0−vb⟩∂T ≤ Chk∥w∥k+1|||vh|||, s(~Qhw,vh) ≤ Chk∥w∥k+1|||vh|||.

Moreover, we have

 l(w,vh)≤Chk∥w∥k+1|||vh|||.

Proof. For the first inequality, from the triangle inequality and Lemma LABEL:err-est1 it follows that

 (a∇w~Qhw−Qh(a∇w),∇wvh) ≤ ∥a∇w~Qhw−Qh(a∇w)∥∥∇wvh∥ ≤ ∥a∇w~Qhw−a(Qh∇w)∥|||vh|||+∥a(Qh∇w)−a∇w∥|||vh||| +∥a∇w−Qh(a∇w)∥|||v