 # Convergence and supercloseness in a balanced norm of finite element methods on Bakhvalov-type meshes for reaction-diffusion problems

In convergence analysis of finite element methods for singularly perturbed reaction–diffusion problems, balanced norms have been successfully introduced to replace standard energy norms so that layers can be captured. In this article, we focus on the convergence analysis in a balanced norm on Bakhvalov-type rectangular meshes. In order to achieve our goal, a novel interpolation operator, which consists of a local weighted L^2 projection operator and the Lagrange interpolation operator, is introduced for a convergence analysis of optimal order in the balanced norm. The analysis also depends on the stabilities of the L^2 projection and the characteristics of Bakhvalov-type meshes. Furthermore, we obtain a supercloseness result in the balanced norm, which appears in the literature for the first time. This result depends on another novel interpolant, which consists of the local weighted L^2 projection operator, a vertices-edges-element operator and some corrections on the boundary.

## Authors

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

 −ε2Δu+bu=fin Ω:=(0,1)2,u=0on ∂Ω, (1.1)

where is a positive parameter. Assume that and are sufficiently smooth and

 b(x,y)≥2β2>0∀(x,y)∈¯Ω,

with a positive constant . Under these conditions on the data in problem (1.1), there exists a unique solution in for all . The solution to problem (1.1) typically exhibits boundary layers of width along all of in the singularly perturbed case of interest.

For singularly perturbed problems, it is popular to introduce layer-adapted meshes Linb:2010-Layer ; Roo1Sty2Tob3:2008-Robust ; Stynes:2005-Steady to fully resolve layers. Then uniform convergence with respect to singular perturbation parameter can be achieved for standard numerical methods. There are two kinds of layer-adapted meshes widely used in the literature, which are Bakhvalov-type mesh (B-mesh) and Shishkin-type mesh (S-mesh) (see Linb:2010-Layer ). There are a lot of research results on convergence theories of finite element methods on S-meshes; see Roo1Sty2Tob3:2008-Robust ; Franz:2008-Singularly ; Dur1Lom2Pri3:2013-Supercloseness ; Zhan1Liu2:2016-Analysis-CL ; Zhan1Liu2:2017-Supercloseness-EL ; Zhan1Styn2:2017-Supercloseness ; Liu1Sty2Zha3:2018-Supercloseness-CIP-EL and references therein.

Although B-meshes usually have better performances than S-meshes, there are very few articles on uniform convergence of finite element methods on the former. One main reason is that B-meshes have specific transition points between the fine and coarse parts, which are independent of mesh parameter. In the meantime, they bring great difficulties to convergence analysis. For example, Lagrange interpolant does not work well for B-meshes. Recently, Zhang and Liu Zhan1Liu2:2020-Optimal ; Zhan1Liu2:2021-Supercloseness proposed an variant of Lagrange interpolant for finite element methods on B-meshes in the case of convection-diffusion equations and succeeded to obtain a uniform convergence of optimal order.

For reaction-diffusion problems, they also proved optimal order of uniform convergence in the natural energy norm in Zhan1Liu2:2020-Convergence . However, energy norm is not strong enough to capture layers as the singular perturbation parameter tends to zero. Thus balanced norms, which are stronger than standard energy norms and characterize layers in an more appropriate way, were introduced in Lin1Styn2:2012-balanced for a mixed finite element method and Roos1Scho2:2015-Convergence for a finite element method. The authors Roos1Scho2:2015-Convergence introduced an

-projection to obtain desired estimations by

-stability of the -projection on Shishkin mesh. To improve estimations in Roos1Scho2:2015-Convergence , the authors Fran1Roos2:2019-Error introduced a new interpolation, which consists a local weighted projection defined on the uniform part of S-meshes. Unfortunately, unlike S-meshes, there is little development on convergence theories in balanced norms on B-meshes.

In this manuscript we analyze convergence theories in the balanced norm introduced in Roos1Scho2:2015-Convergence for th () order finite element method on Bakhvalov-type rectangular meshes. For this purpose, we propose a novel interpolant according to the structures of B-meshes and layer functions. This interpolant consists of a local weighted projection defined on a proper mesh subdomain and the Lagrange interpolant. To prove the convergence of optimal order in the balanced norm, we must take into account the scales of the meshes and different stabilities of the projection. The optimal order convergence is also supported by our numerical experiments. Furthermore, we propose another novel interpolant operator, which consists of the local weighted projection operator, a vertices-edges-element operator and some corrections on the boundary. By careful derivations, we obtain a supercloseness result, which appears in the literature for the first time. Here “supercloseness” means that the convergence order for the error between some interpolation of the solution and the numerical solution in some norm is greater than the order for in the same norm.

The rest of the paper is organized as follows. In Section 2 we present a priori information of the solution to (1.1), then introduce Bakhvalov-type meshes, finite element methods and some preliminary results. In Section 3 we give a new interpolant and prove uniform convergence of optimal order in the balanced norm. Supercloseness result is given in Section 4 by means of another novel interpolant. In Section 5, numerical results illustrate our theoretical results.

Let . In this article, we will write for the inner product in , , , and for the standard norms in , , and the standard seminorm in , respectively. If , the subscript will be omitted from the above norm designations. Throughout the paper, all constants and are independent of and ; the constants are generic while subscripted constants are fixed.

## 2 Finite element method on Bakhvalov-type mesh

### 2.1 Regularity results

To construct layer-adapted meshes and analyze uniform convergence, we need a priori information of the solution to (1.1), such as pointwise estimations of the derivatives of the solution, the locations and widths of layers.

For this aim, we give the following assumption on the solution to (1.1) according to Han1Kell2:1990-Differentiability ; Liu1Mad2Sty3etc:2009-two-scale ; Cla1Gra2Ori3:2005-parameter .

###### Assumption 1.

The solution of (1.1) can be decomposed as

 u=v0+4∑i=1wi+4∑i=1zi∀x∈¯Ω, (2.1)

where is the regular part, each is a boundary layer function and each is a corner layer function. For and , there exists a constant such that

 ∣∣∂mx∂nyv0(x,y)∣∣≤C(1+ε(k+1)−m−n)for 0≤m+n≤k+3, (2.2) ∣∣∂mx∂nyw1(x,y)∣∣≤C(1+ε(k+1)−m)ε−ne−βy/εfor 0≤m+n≤k+2, (2.3) ∣∣∂mx∂nyz1(x,y)∣∣≤Cε−(m+n)e−β(x+y)/εfor 0≤m+n≤k+2, (2.4)

and similarly for the remaining terms. Here denote by .

In the following analysis, we will denote by .

### 2.2 Bakhvalov-type meshes

Two Bakhvalov-type meshes will be discussed. Let be divisible by 4. The first Bakhvalov-type mesh is introduced in Roos-2006-Error and defined by

 xi=yi=ψ(i/N)=⎧⎪ ⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪ ⎪⎩−σεβln(1−4(1−ε)i/N)for i=0,…,N/4,d1(i/N−1/4)+d2(i/N−3/4)for i=N/4,…,3N/4,1+σεβln(1−4(1−ε)(1−i/N))for i=3N/4,…,N, (2.5)

where will be defined later and , are used to ensure the continuity of at and . The second Bakhvalov-type mesh is introduced in Kopteva:1999-the ; Kopt1Save2:2011-Pointwise and its mesh generating function is

 φ(t)=⎧⎪ ⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪ ⎪⎩−σεβln(1−4t)for t∈[0,ϑ],d3(t−ϑ)+d4(t−1+ϑ)for% t∈(ϑ,1−ϑ),1+σεβln(1−4(1−t))for t∈[1−ϑ,1] (2.6)

where will be specified later, with some positive constant independent of and , and are chosen so that is continuous at and . The original Bakhvalov mesh Bakhvalov:1969-Towards can be recovered from (2.6) by setting with .

###### Assumption 2.

Assume that in our analysis. In practice it is not a restriction.

For technical reasons, we also assume

 σ4eβ≤C1≤12max{σβ,14}. (2.7)

Assume . Under Assumption 2 and (2.7), we have and for meshes (2.5) and (2.6). The location of and conditions imposed on and will simplify our later analysis without changing the essential difficulties in our analysis.

The mesh points are or for . By drawing lines parallel to the axis through mesh points , we obtain a Bakhvalov-type rectangular mesh with equidistant cells in the coarse region and anisotropic cells in the layer region . The triangulation is denoted by . Denote by for the element and by for a generic rectangular element, which dimensions are written as and . Define for .

In the following lemma, we collect some important properties possessed by Bakhvalov-type meshes, which are important for convergence analysis. The reader is referred to Zhan1Liu2:2020-Optimal for the detailed proof.

###### Lemma 1.

Let Assumption 2 hold true. On Bakhvalov-type mesh (2.5) or (2.6), one has

 h0≤h1≤…≤hN/4−2, (2.8) σ4βε≤hN/4−2≤σβε, (2.9) Cε≤hN/4−1≤CN−1, (2.10) C4N−1≤hi≤C5N−1N/4≤i≤N/2, (2.11) xN/4−1=x3N/4+1≥CσεlnN,xN/4=x3N/4≥Cσεln(1/ε). (2.12)

Note for .

Let . Then one has

 hμie−βxi/ε≤CεμN−μfor 0≤i≤i∗ and 0≤μ≤σ. (2.13)

Similar bounds hold for the variable .

### 2.3 Finite element method

Now we present the finite element method for problem (1.1). First, the weak form of problem (1.1) is written as

 {find u∈V such that for all v∈Va(u,v):=ε2(∇u,∇v)+(bu,v)=(f,v), (2.14)

with . The natural energy norm derived from is

 ∥v∥ε:=(ε2|v|21+∥v∥2)1/2.

The bilinear form is coercive with respect to this energy norm, i.e.,

 a(v,v)≥min{2β2,1}∥v∥2ε∀v∈V. (2.15)

From the Lax-Milgram lemma, the weak formulation (2.14) has a unique solution.

Let denote th rectangular finite element. We introduce the finite element spaces

 VN:={v∈H1(Ω):v|τ∈Qk(τ)∀τ∈TN} and VN0:=VN∩H10(Ω).

Clearly, . When we replace the infinite dimensional space with the finite dimensional space , we get the th order finite element method

 {find uN∈VN0 such that for all vN0∈VNa(uN,vN)=(f,vN). (2.16)

Also, it is easy to very the coercivity

 a(vN,vN)≥C∥vN∥2ε∀vN∈VN0, (2.17)

and the Galerkin orthogonality

 a(u−uN,vN)=0∀vN∈VN0, (2.18)

where (2.14) and (2.16) have been used. Furthermore, we introduce the balanced norm in Roos1Scho2:2015-Convergence , which is defined by

 ∥v∥b:=(ε|v|21+∥v∥2)1/2. (2.19)

Clearly, the balanced norm is stronger than the energy norm in the case of . Furthermore, the former is better suited to capture of layers. For example, for a typical layer function , and are of order and of order , respectively. These orders imply that the balanced norm is more appropriate to capture layers than the energy norm when . The reader is also referred to Roos1Scho2:2015-Convergence ; Con1Fra2Lud3etc:2018-Finite for discussions on these two norms.

## 3 Uniform convergence

For convergence analysis in the balanced norm, we will present an interpolation operator, which components will be introduced at first.

Set and . Introduce a weighted -projection as follows: for , find such that

 (b(s−πs),vN)Ω∗0=0∀vN∈WN,

where . Of course, one has the -stability

 ∥πv∥Ω∗0≤C∥v∥Ω∗0. (3.1)

Denote by the Langrange interpolation operator from to . Furthermore, define by

 χ(sl,tm)={1(sl,tm)∈∂Ω∗0,0otherwise,

where are the interpolation points of the Lagrange interpolation.

Recall . Then the interpolation used in convergence analysis is defined by

 Pcu=P1v0+P2w, (3.2)

where

 P1v0:={πv0%in$Ω∗0$,I[(1−χ)v0+χπv0] in% Ω∖Ω∗0 (3.3)

and

 P2w:={0 in Ω∗0,I[(1−χ)w] in Ω∖Ω∗0. (3.4)

Clearly, .

The following lemma provides some pointwise bounds for errors between the Lagrange interpolant and the projection.

###### Lemma 2.

For introduced in Assumption 1, one has

 ∥Iv0−πv0∥∞,∂Ω∗0+∥Iv0−πv0∥∞,Ω∗0+∥v0−πv0∥∞,Ω∗0≤CN−(k+1), ∥∇(Iv0−πv0)∥Ω0≤CN−k,∥∇(Iv0−πv0)∥Ω∗0∖Ω0≤Cε−1/2N−(k+1).
###### Proof.

Standard interpolation theory and Lemma 1 imply . From the -stability of the -projection (Crou1Tome2:1987-stability, , Theorem 1), one has

 ∥Iv0−πv0∥∞,Ω∗0=∥π(Iv0−v0)∥∞,Ω∗0≤C∥Iv0−v0∥∞,Ω∗0≤CN−(k+1), ∥Iv0−πv0∥∞,∂Ω∗0≤∥Iv0−πv0∥∞,Ω∗0≤CN−(k+1), ∥v0−πv0∥∞,Ω∗0≤∥v0−Iv0∥∞,Ω∗0+∥Iv0−πv0∥∞,Ω∗0≤CN−(k+1).

Hölder inequalities, inverse inequalities, the -stability of the -projection (Crou1Tome2:1987-stability, , Theorem 1) and Lemma 1 yield

 ∥∇(Iv0−πv0)∥Ω∗0∖Ω0≤ meas(Ω∗0∖Ω0)1/2∥∇(Iv0−πv0)∥∞,Ω∗0∖Ω0 ≤ Ch−1/2N/4−1∥Iv0−πv0∥∞,Ω∗0∖Ω0 = Ch−1/2N/4−1∥π(Iv0−v0)∥∞,Ω∗0 ≤ Ch−1/2N/4−1∥Iv0−v0∥∞,Ω∗0 ≤ Cε−1/2N−(k+1).

Inverse inequalities and the -stability of the -projection (3.1) yield

 ∥∇(Iv0−πv0)∥Ω0≤CN∥Iv0−πv0∥Ω0≤CN∥Iv0−πv0∥Ω∗0 ≤ CN∥Iv0−v0∥Ω∗0≤CN−k.

###### Remark 1.

Define

 Vh={v∈C[xN/4−1,x3N/4+1];v|(sj,sj+1)∈Pk,j=0,…,N/2+1}

with for . In the same way as the proof of (Crou1Tome2:1987-stability, , Theorem 1), we could easily prove the -stability of the -projection , that is

 ∥πhv∥∞,(xN/4−1,x3N/4+1)≤C∥v∥∞,(xN/4−1,x3N/4+1)∀v∈L∞(xN/4−1,x3N/4+1).

Furthermore, from tensor product we could easily obtain

 ∥πv∥∞,Ω∗0≤C∥v∥∞,Ω∗0∀v∈L∞(Ω∗0).

The error is split as follows:

 u−uN=(u−Pcu)+(Pcu−uN)=:η+ξ.

From the coercivity (2.17) and the Galerkin orthogonality (2.18) we have

 C∥ξ∥2ε≤a(Pcu−u,ξ)≤ε1/2∥η∥b∥ξ∥ε+|(bη,ξ)|. (3.5)

In the following analysis, we give the estimations on each term in the right-hand side of (3.5).

###### Lemma 3.

Let . Under Assumptions 1 and 2 we have

 |(bη,ξ)|≤Cε1/2N−(k+1)ln1/2N∥ξ∥.
###### Proof.

Our arguments are based on the following splitting

 (bη,ξ)= (b(v0−P1v0),ξ)+(b(w−P2w),ξ) = (b(v0−πv0),ξ)Ω∗0+(b(v0−Iv0),ξ)Ω∖Ω∗0+(bI[χ(v0−πv0)],ξ)Ω∗∗0∖Ω∗0 +(bw,ξ)Ω∗0+(b(w−Iw),ξ)Ω∖Ω∗0+(bI(χw),ξ)Ω∗∗0∖Ω∗0 =:

From the definition of , we obtain

 I=0. (3.6)

From (Zhan1Liu2:2020-Convergence, , Lemma 4), one has

 |II|+|V|≤ C(∥v0−Iv0∥∞,Ω∖Ω∗0+∥w−Iw∥∞,Ω∖Ω∗0)∥ξ∥1,Ω∖Ω∗0 (3.7) ≤ CN−(k+1)meas1/2(Ω∖Ω∗0)∥ξ∥Ω∖Ω∗0 ≤ Cε1/2N−(k+1)ln1/2N∥ξ∥.

Lemmas 1 and 2, (2.3) and (2.4) yield

 |III|+|VI|≤ C(∥Iv0−πv0∥∞,∂Ω∗0+∥Iw∥∞,∂Ω∗0)∥ξ∥1,Ω∗∗0∖Ω∗0 (3.8) ≤ C(N−(k+1)+N−σ)h1/2N/4−2∥ξ∥Ω∗∗0∖Ω∗0 ≤ Cε1/2N−(k+1)∥ξ∥.

From (2.3) and (2.4), we get

 |IV|≤C∥w∥Ω∗0∥ξ∥Ω∗0≤Cε1/2N−σ∥ξ∥. (3.9)

Collecting (3.6)–(3.9), we are done. ∎

###### Lemma 4.

Let . Under Assumptions 1 and 2 we have

 ∥η∥b≤CN−k.
###### Proof.

From (3.2), (3.3) and (3.4) one has

 ∥η∥b≤ ∥v0−πv0∥b,Ω∗0+∥v0−Iv0∥b,Ω∖Ω∗0+∥I(χ(v0−πv0))∥b,Ω∗∗0∖Ω∗0 +∥w−Iw∥b,Ω∖Ω∗0+∥I(χw)∥b,Ω∗∗0∖Ω∗0 =: S1+S2+S3+S4+S5.

To analyze , we need the following bounds

 ∥∇(v0−πv0)∥Ω0≤ ∥∇(v0−πv0)∥Ω∗0∖Ω0≤ ∥∇(v0−Iv0)∥Ω∗0∖Ω0+∥∇(Iv0−πv0)∥Ω∗0∖Ω0 ≤ CN−k+Cε−1/2N−(k+1), ∥v0−πv0∥Ω∗0≤ ∥v0−Iv0∥Ω∗0+∥Iv0−πv0∥Ω∗0 ≤ CN−(k+1)+C∥Iv0−πv0∥∞,Ω∗0≤CN−(k+1),

which could be derived from standard interpolation theories and Lemma 2. Then we obtain

 |S1|≤ ε1/2∥∇(v0−πv0)∥Ω∗0+∥v0−πv0∥Ω∗0≤Cε1/2N−k+CN−(k+1). (3.10)

Similar to (Zhan1Liu2:2020-Convergence, , Lemma 4), we have

Imitating the proof of Lemma 5 in Zhan1Liu2:2020-Convergence and replacing (Zhan1Liu2:2020-Convergence, , (3.15)) by , we obtain

 |w1−Iw1|1,Ω∖Ω∗0≤Cε−1/2N−k,

and similarly have

 |w−Iw|1,Ω∖Ω∗0≤Cε−1/2N−k.

Thus we obtain

 S2+S4≤CN−k. (3.11)

Inverse inequalities, Lemmas 1 and 2 yield

 S3≤ ε1/2∥∇(I(χ(v0−πv0)))∥Ω∗∗0∖Ω∗0+∥I(χ(v0−πv0))∥Ω∗∗0∖Ω∗0 (3.12) ≤ C(ε1/2h−1N/4−2+1)∥I(χ(v0−πv0))∥Ω∗∗0∖Ω∗0 ≤ C(ε1/2h−1N/4−2+1)h1/2N/4−2∥Iv0−πv0∥∞,∂Ω∗0 ≤ CN−(k+1).

Similarly, one has

 S5≤CN−σ. (3.13)

Collecting (3.10)–(3.13), we are done. ∎

Similar to (Fran1Roos2:2019-Error, , Theorem 2.6), we obtain the following theorem from Lemmas 3 and 4.

###### Theorem 1.

Let . Let Assumptions 1 and 2 hold. Then for the exact solution to (1.1) and the numerical solution to (2.16) on Bakhvalov-type rectangular mesh (2.5) or (2.6), one has

 ∥u−uN∥b≤CN−k.

## 4 Supercloseness

In order to derive the supercloseness result in the balanced norm, we need another novel interpolant, which will be described in the following.

Instead of Langrange interpolant operator used in the previous section, we introduce an vertices-edges-element interpolation operator (see Lin1Yan2Zho3:1991-rectangle ; Styn1Tobi2:2008-Using ). This interpolant is used for superconvergence analysis of the diffusion term. First we define the interpolant operator on the reference element , whose vertices and edges are denoted by and respectively for . Let . The operator is determined by continuous linear functionals , which are defined by

 ^v→^v(^ai)i=1,\ldots,4, ^v→∫^ei^vqds∀q∈Pk−2(^ei)i=1,\ldots,4, ^v→∫^τ^vqdxdy∀q∈Qk−2(^τ).

From (Styn1Tobi2:2008-Using, , Lemma 3), the operator is uniquely determined. Then using the affine transformation to map from to an arbitrary , one obtains the corresponding interpolation operator . At last a continuous global interpolation operator is defined by setting

 (Av)|τ:=Aτ(v|τ)∀τ∈TN.

Besides, we denote by degree of freedom (DoF) of , which originates from the linear functional .

Recall and . The operator is defined by

 S1w1:={0 in Ω∖Ωw1,Aw1−B1w1 in Ωw1, (4.1)

where satisfies

 F(B1w1)={F(w1)if F is the DoF of VN attached to [0,1]×{y=yj∗+1},0otherwise.

If is the DoF of attached to , then it must be one of the following forms

 v→v(xi,yj∗+1)for i=0,…,N, v→j+1hj+1i∫xi+1xiv(s,yj∗+1)(s−xi)jdsfor i=0,…,N−1 and j=0,…,k−2.

In fact, is introduced for the continuity of at . Set . The operator is defined by

 T