# A Reduced Local Discontinuous Galerkin Method for Nearly Incompressible Linear Elasticity

Authors: Xuehai Huang
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

A reduced local discontinuous Galerkin (RLDG) method for nearly incompressible
linear elasticity is proposed in this paper, which is locking-free. RLDG method can be
formally regarded as a special case of LDG method withC11=0. However, RLDG method is actually not covered by LDG method, where C11 must be chosen to be positive to ensure the stability of LDG method. RLDG method can also be considered as
the localization of some symmetric nonconforming mixed finite element method. The
implementation of RLDG method is discussed. By introducing a lifting operator as
LDG method, RLDG method can be rewritten as primal formulation with unknown
displacement only. Next, we obtain that the convergence rates of the approximation
to stress tensor in energy norm and displacement in L2-norm are O(hk) and O(hk+1), respectively, which are both uniform with respect to λ. Moreover, we obtain a H(div)-conforming displacement by projecting the displacement and corresponding numerical
trace of RLDG method into the Raviart-Thomas element space. And then we analyze
the error estimates of this postprocessed displacement in H(div)-seminorm and L2-norm, which are also uniform with respect to λ. Finally, some numerical results are
shown to demonstrate the theoretical results.


---

## Body

## 1. Introduction

Volume locking phenomenon will happen when solving nearly incompressible linear elasticity problem by standard displacement method using continuous lower order finite elements (cf. [1–3]); that is, the convergence rates for these lower order conforming finite element methods will deteriorate as the Lamé constant λ tends to +∞. For this reason, robust numerical methods with respect to λ are necessary to deal with volume locking phenomenon, and these methods are also called locking-free methods. As for standard displacement method, Babuška and Suri showed detailedly in [3] that it is locking-free whenever the polynomial order p of the local shape space is no less than 4 and locked for p<4 on triangular meshes, and locking phenomenon always takes place for any p≥1 on quadrilateral meshes. Besides higher order h-version conforming finite element methods, p-version finite element method in [4], nonconforming finite element methods in [5–8], and mixed finite element methods in [9–11] are also suggested to overcome the volume locking phenomenon.On the other hand, discontinuous Galerkin (DG) methods which possess the flexibility of constructing feasible local shape function spaces and the advantage of capturing non-smooth or oscillatory solutions effectively have been applied to solve nearly incompressible linear elasticity. Hansbo and Larson analyzed symmetric interior penalty discontinuous Galerkin (SIPDG) methods in [12, 13], whereas a dimensionless penalty parameter appeared in the stabilized term cannot be quantified in theory and can only be determined by numerical experience. A locking-free nonsymmetric interior penalty discontinuous Galerkin (NIPDG) method was developed by Wihler in [14]. And Wihler also designed an adaptive algorithm for NIPDG method in [15], whose robust a posteriori error analysis was given by Bridgeman and Wihler in [15, 16]. A local discontinuous Galerkin (LDG) method based on a first-order system with three mathematical equations was presented and discussed by Cockburn et al. in [17], and a post-processed H(div)-conforming displacement was achieved by projecting the displacement and corresponding trace of LDG method into the Raviart-Thomas element space. Using symmetric normal-normal continuous element for stress and tangential continuous Nédélec element for displacement, Pechstein and Scho¨berl proposed a locking-free mixed finite element method in [18], which is in fact some kind of discontinuous Galerkin method. Through introducing a DG derivative, a family of discontinuous Galerkin mixed methods were designed by Shen and Lew in [19] which were proved to be locking-free as well. They also observed that locking might appear if Dirichlet boundary conditions were imposed strongly rather than weakly.In this paper, we intend to propose a reduced local discontinuous Galerkin (RLDG) method for nearly incompressible linear elasticity. RLDG method can be formally regarded as a special case of LDG method proposed in [20] with C11=0. However, RLDG method is actually not covered by LDG method, where C11 must be chosen to be positive to ensure the stability of LDG method, whereas RLDG method drops this requirement by choosing large enough finite element space for stress to match finite element space for displacement. On the other hand, RLDG method can also be considered as the localization of the symmetric nonconforming mixed finite element method given in [21]. Based on this observation, the unique solvability of RLDG method is proved using the similar argument as the method in [21]. Since there is not any continuity across the interior faces for the finite element space for stress, Schur complement is easy to obtain, and thus implementation of RLDG method is not difficult. On the other hand, RLDG method can be rewritten as primal formulation with unknown displacement only through introducing a lifting operator as LDG method in [20]. According to the error equations, an interpolation operator designed in [21] and BDM interpolation operator, we obtain that the convergence rate of the approximation to stress tensor and displacement in energy norm is O(hk). And by duality argument, the convergence rate of the approximation to displacement in L2-norm is shown to be O(hk+1). These two convergence rates are both uniform with respect to λ; hence our RLDG method is locking-free. Moreover, we obtain a H(div)-conforming displacement by projecting the displacement and corresponding numerical trace of RLDG method into the Raviart-Thomas element space, just as the postprocessing technique in [17]. And then we analyze the error estimates of this post-processed displacement in H(div)-seminorm and L2-norm, which are also uniform with respect to λ. Finally, some numerical results are shown to demonstrate the theoretical results.The rest of this paper is organized as follows. The reduced local discontinuous Galerkin method for nearly incompressible linear elasticity is presented in Section2. In Section 3, we give a robust error analysis for RLDG method. A postprocessing of the displacement of RLDG method is discussed in Section 4. And in Section 5, some numerical results are included to confirm our theoretical convergence orders. In the final section, we give some concluding statements.


## 2. The RLDG Method for Linear Elasticity

Assume thatΩ⊂ℝd(d=2,3) is a bounded polygon or polyhedron. Let σ=(σij)d×d be the stress, let u=(u1,…,ud)t be the displacement, and let f=(f1,…,fd)t be the body force. Denote by ε(u):=(εij(u))d×d the linearized strain tensor with εij(u):=(∂ui/∂xj+∂uj/∂xi)/2, tr the trace operator, and ∇· the divergence operator. Consider linear elasticity in the stress-displacement formulation:
(1)𝒜σ-ε(u)=0inΩ,-∇·σ=finΩ,u=0on∂Ω,
where 𝒜 is the compliance tensor of the fourth order defined by
(2)𝒜σ=12μ(σ-λdλ+2μ(trσ)δ).
Here, tr is the trace operation of a second-order tensor (a matrix), δ:=(δij)d×d the usual Kronecker delta, and the positive constants λ∈(0,+∞) and μ∈[μ0,μ1] are the Lamé coefficients of the elastic material under consideration with 0<μ0<μ1<+∞.To define a reduced LDG method for solving problem (1), we first introduce some notations frequently used later on. Denote the space of all symmetric d×d tensor by 𝕊. Given a bounded domain G⊂ℝd and a nonnegative integer m, let Hm(G) be the usual Sobolev space of functions on G, and let Hm(G,𝕏) be the usual Sobolev space of functions taking values in the finite-dimensional vector space 𝕏 for 𝕏 being 𝕊 or ℝd. The corresponding norm and seminorm are denoted, respectively, by ∥·∥m,G and |·|m,G. If G is Ω, we abbreviate them by ∥·∥m and |·|m, respectively. Let H0m(G) be the closure of C0∞(G) with respect to the norm ∥·∥m,G. We also denote by H(div,G,𝕊) the Sobolev space consisting of all L2(G,𝕊) functions whose divergence are square integrable. Let
(3)H0(div,Ω):={v∈L2(Ω,ℝd);∇·v∈L2(Ω),v·n|∂Ω=0}.Pm(G) stands for the set of all polynomials in G with the total degree no more than m, and Pm(G,𝕏) denotes the tensor or vector version of Pm(G) for 𝕏 being 𝕊 or ℝd, respectively.Let{𝒯h}h>0 be a regular family of triangulations of Ω (cf. [22, 23]). For each K∈𝒯h, we denote by hK the diameter of K and by ρK the diameter of the biggest ball included in K. Let ℰh be the union of all faces of the triangulation 𝒯h and let ℰhi be the union of all interior faces of the triangulation 𝒯h. For a face e∈ℰh, he is set to be the diameter of e. The triangulations we consider can have hanging nodes, but they have to be regular. Moreover, we restrict the ratio of the sizes of neighbor element domains. To formally state this property, we need to introduce the set 〈K,K′〉 defined as follows:
(4)〈K,K′〉:={∅,ifmeas(∂K∩∂K′)=0,interiorof∂K∩∂K′,otherwise.
Thus we assume that there exists a positive constant C2<1 such that, for each element K∈𝒯h,
(5)C2≤hK′hK≤1C2∀K′satisfying〈K,K′〉≠∅.
This assumption forbids the situation where the mesh is indefinitely refined in only one of two adjacent subdomains. Denote by ∇h and εh the elementwise counterparts of ∇ and ε related to 𝒯h, respectively.LetK+ and K- be two adjacent elements of 𝒯h. Let x be an arbitrary point of the set e′=〈K+,K-〉, and let n+ and n- be the corresponding outward unit normals at that point. For a vector-valued function v smooth inside each element K±, let us denote by v± the trace of v on e′ from the interior of K±. Then we define averages and jumps at x∈e′ as follows:
(6){v}=12(v++v-),[v]=v+·n++v-·n-,{τ}=12(τ++τ-),[τ]=τ+·n++τ-·n-.
If x is on a face e lying on the boundary ∂Ω, the above terms are defined by
(7){v}=v,[v]=v·n,{τ}=τ,[τ]=τ·n,
where n is the unit outward normal vector on ∂Ω. We define a matrix valued jump ⟦·⟧ of a vector v as follows:
(8)⟦v⟧=12(v+⊗n++n+⊗v++v-⊗n-+n-⊗v-),ifx∈e∈ℰhi,⟦v⟧=12(v⊗n+n⊗v),ifx∈e∈ℰh∩∂Ω,
where v⊗n denote the matrix whose (i,j)th component is vinj for two vectors v and n.Based on the triangulation𝒯h, let
(9)Σ:={τ∈L2(Ω,𝕊);τ|K∈H(div,K,𝕊)∀K∈𝒯h},V:={v∈L2(Ω,ℝd);v|K∈H1(K,ℝd)∀K∈𝒯h}.
The corresponding finite element spaces are given by
(10)Σh:={τ∈L2(Ω,𝕊);τ|K∈SK∀K∈𝒯h},Vh:={v∈L2(Ω,ℝd);v|K∈Pk(K,ℝd)∀K∈𝒯h},
where, for each K∈𝒯h, SK is a finite-dimensional symmetric tensor space of polynomials in K containing Pk+1(K,𝕊), with k≥1. Let
(11)Wh≔{∫e‍[τ]τ∈L2(Ω,𝕊);τ|K∈Pk+1(K,𝕊)∀K∈𝒯h,and∫e‍[τ]·vds=0∀v∈Pk(e,ℝd)foreache∈ℰhi},
which was used to construct a family of symmetric nonconforming mixed finite elements for linear elasticity in [21]. It is obvious that Wh⊂Σh.With former notations, we propose a reduced local discontinuous Galerkin (RLDG) method for nearly incompressible linear elasticity problem (1). Find (σh,uh)∈Σh×Vh such that
(12)a(σh,τ)+b(τ,uh)=0∀τ∈Σh,(13)b(σh,v)=-∫Ω‍f·vdx∀v∈Vh,
where
(14)a(σ,τ):=∫Ω‍𝒜σ:τdx,b(τ,v):=∑K∈𝒯h‍∫K‍(∇·τ)·vdx-∫ℰhi‍[τ]·{v}ds.
For comparison, let us recall the LDG method proposed in [20] in the following. Find (σh,uh)∈Σh×Vh such that
(15)a(σh,τ)+b(τ,uh)=0∀τ∈Σh,b(σh,v)-c(uh,v)=-∫Ω‍f·vdx∀v∈Vh,
where
(16)c(u,v)=∫ℰh‍C11⟦u⟧:⟦v⟧ds
with C11>0. Formally, RLDG method (12)-(13) can be viewed as a special case of LDG method (15) with C11=0. These two methods are actually different. LDG method requires C11 to be positive to ensure the stability, whereas RLDG method drops this requirement by choosing large enough finite element space Σh to match Vh in the sense that Σh⊃Wh.On the other hand, according to the definition ofWh, the following holds:
(17)b(τ,v)=∑K∈𝒯h‍∫K‍(∇·τ)·vdx∀τ∈Wh,v∈Vh.
With (17), RLDG method (12)-(13) will be reduced to the symmetric nonconforming mixed finite element method designed in [21] by taking Σh=Wh. Thus RLDG method can be also regarded as the discretization of the symmetric nonconforming mixed finite element method when SK=Pk+1(K,𝕊), for there is not any continuity of finite element space Σh across interior faces. The merit of RLDG method over the symmetric nonconforming mixed finite element method is that we do not need to solve a painful saddle problem in RLDG method. By the way, the symmetric nonconforming mixed finite element method is solved by the hybridized technique in [21].Theorem 1.
Reduced local discontinuous Galerkin method (12)-(13) is uniquely solvable.Proof.
SinceΣh and Vh are both finite-dimensional, it is enough to show that if f=0, then σh=0, uh=0. Taking τ=σh in (12) and v=uh in (13) and subtracting (13) from (12), it follows that a(σh,σh)=0. Thus σh=0 from the definition of a(·,·). Then (12) is reduced to
(18)b(τ,uh)=0∀τ∈Σh.
Since Wh⊂Σh and (17), thus
(19)∑K∈𝒯h‍∫K‍(∇·τ)·uhdx=0∀τ∈Wh.
Then we can get uh=0 by adopting the same argument used in the proof of Corollary 2.4 in [21], which ends the proof.Next, let us consider the implementation of RLDG method (12)-(13). Let {φi}i=1N1 and {ϕj}j=1N2 be the bases of Σh and Vh, respectively. Write
(20)σh=∑i=1N1‍σiφi,uh=∑j=1N2‍ujϕj.
And denote
(21)Aij:=a(φi,φj),Bij:=b(φj,ϕi),Fi:=∫Ω‍fϕidx.
Then RLDG method (12)-(13) has the matrix representation
(22)(ABTB0)(YX)=(0-F)
with A=(Aij)N1×N1, B=(Bij)N2×N1, F=(Fi)N2×1, X=(ui)N2×1, and Y=(σi)N1×1. Since Σh is a piecewise polynomial space without any continuity across faces in ℰhi, we can rearrange the basis {φi}i=1N1 such that A is a block diagonal matrix associated with triangulation 𝒯h. Thus it is trivial to compute the inverse of matrix A. Then the last symmetric indefinite linear system can be reduced to
(23)BA-1BTX=F.
The new linear system is much easier to solve compared to original system, for new linear system has only unknown X and BA-1BT is symmetric and positive definite. After obtaining X by solving the new linear system, we can get Y from
(24)Y=-A-1BTX.
Therefore, RLDG method (12)-(13) can be implemented efficiently.In the end of this section, we will intend to drive the primal formulation of RLDG method (12)-(13) using the argument taken in [20]. For this, assume that
(25)𝒜Σh⊂Σh
in this paragraph. Define a global lifting operator r:L2(ℰh,𝕊)→Σh by (cf. [20])
(26)∫Ω‍r(ϕ):τdx=-∫ℰh‍ϕ:{τ}ds∀τ∈Σh,ϕ∈L2(ℰh,𝕊).
It holds the following inequality for lifting operator (cf. [20]):
(27)∥r(ϕ)∥02≲∑e∈ℰh‍he-1∥ϕ∥0,e2∀ϕ∈L2(ℰh,𝕊).
From (12), by integration by parts and definition of lifting operator r we have
(28)∫Ω‍𝒜σh:τdx=∑K∈𝒯h‍∫K‍ε(uh):τdx-∫ℰh‍{τ}:⟦uh⟧ds=∑K∈𝒯h‍∫K‍ε(uh):τdx+∫Ω‍r(⟦uh⟧):τdx∀τ∈Σh.
Thus we have from (25)
(29)σh=𝒜-1(εh(uh)+r(⟦uh⟧)).
Substituting (29) into (13) and using integration by parts and definition of lifting operator r, we gain the primal formulation of RLDG method (12)-(13) as follows. Find uh∈Vh such that
(30)∫Ω‍𝒜-1(εh(uh)+r(⟦uh⟧)):(εh(v)+r(⟦v⟧))dx=∫Ω‍f·vdx∀v∈Vh.


## 3. A Priori Error Analysis

In this section, we provide an a priori error analysis for RLDG method (12)-(13). For a function v∈L2(Ω) with v|K∈Hm(K) for all K∈𝒯h, let ∥v∥m,h and |v|m,h be the usual broken Hm-type norm and seminorm of v:
(31)∥v∥m,h:=(∑K∈𝒯h‍∥v∥m,K2)1/2,|v|m,h:=(∑K∈𝒯h‍|v|m,K2)1/2.
If v is a vector-valued or tensor-valued function, the corresponding ∥·∥m,h and |·|m,h are defined in the similar manners. For a vector or tensor v, its length |v| is (v·v)1/2 or (v:v)1/2. Here the symbol : denotes the double dot product operation of tensors. Throughout this paper, we use the notation “≲⋯” to mean that “≤C⋯,” where C is a generic positive constant independent of local element sizes and λ, which may take different values at different appearances.According to the definition of𝒜, we have, for any second-order tensor τ,
(32)𝒜-1τ=2μτ+λ(trτ)δ.
Given a bounded domain G⊂ℝd, we define two energy-type norms for L2(G,𝕊) in the following. For any τ∈L2(G,𝕊), let
(33)∥τ∥𝒜,G2:=∫G‍𝒜τ:τdx,∥τ∥𝒜-1,G2:=∫G‍𝒜-1τ:τdx.
By the definition of 𝒜, it is obvious that
(34)∥τ∥𝒜,G≤12μ∥τ∥0,G∀τ∈L2(G,𝕊).
Then we denote by (·,·)𝒜,G the Hilbert inner product corresponding to energy norm ∥·∥𝒜,G. If G is Ω, we abbreviate ∥τ∥𝒜,Ω, ∥τ∥𝒜-1,Ω, and (·,·)𝒜,Ω by ∥τ∥𝒜, ∥τ∥𝒜-1, and (·,·)𝒜, respectively. For any τ1,τ2∈L2(G,𝕊), we obtain from Cauchy-Schwarz inequality
(35)∫Ω‍τ1:τ2dx=(τ1,𝒜-1τ2)𝒜≤∥τ1∥𝒜∥𝒜-1τ2∥𝒜=∥τ1∥𝒜∥τ2∥𝒜-1.
For later uses, we introduce two interpolation operators. Let Πh be the linear operator from H(div,Ω,𝕊)∩Lp(Ω,𝕊) for any p>2 to Wh defined by (3.12) in [21], and let Qh be the L2 projection operator from V onto the finite element space Vh. From (3.5) in [21], it holds for all τ∈H(div,Ω,𝕊)∩Lp(Ω,𝕊)(36)∇h·(Πhτ)=Qh(∇·τ).Lemma 2.
For allτ∈H(div,Ω,𝕊)∩Lp(Ω,𝕊) and v∈Vh, it follows
(37)b(τ-Πhτ,v)=0.Proof.
Since[τ]|e=0 for any e∈ℰhi and Πhτ∈Wh, we have by (17)
(38)b(τ-Πhτ,v)=∫Ω‍(∇·τ)·vdx-∫Ω‍(∇h·(Πhτ))·vdx.
According to (36) and the definition of Qh,
(39)b(τ-Πhτ,v)=∫Ω‍(∇·τ)·vdx-∫Ω‍(Qh(∇·τ))·vdx=0.
This ends the proof.Next, let us study the Galerkin orthogonality of RLDG method (12)-(13). Using integration by parts, we can get from problem (1) that
(40)a(σ,τ)+b(τ,u)=0∀τ∈Σ,b(σ,v)=-∫Ω‍f·vdx∀v∈V.
Subtracting (12)-(13) from (40), we obtain the following error equations:
(41)a(σ-σh,τ)+b(τ,u-uh)=0∀τ∈Σh,(42)b(σ-σh,v)=0∀v∈Vh.To derive the a priori error estimates of RLDG method (12)-(13), we require as in [17] some regularity result for the solution (σ,u) of problem (1). We say the solution (σ,u) of problem (1) is (s+2)-regular, where s is a nonnegative integer whenever σ∈Hs+1(Ω,𝕊),u∈Hs+2(Ω,ℝd) and there holds (cf. [4, 5, 17, 24])
(43)∥σ∥s+1+∥u∥s+2+λ∥∇·u∥s+1≲∥f∥s
for f∈Hs(Ω,ℝd). We also assume in this paper that
(44)tr(Σh)δ⊂Σh.
Moreover, introduce the usual BDM element space (cf. [25, 26])
(45)VhBDM:={v∈Vh;[v]|e=0∀e∈ℰh}.
Let Ih:(H01(Ω))d→VhBDM be the BDM nodal interpolation operator. For simplicity, we still write Πh and Ih for Πh|K and Ih|K. Using Lemma 2, error equations (41)-(42) can be rewritten as
(46)a(σ-σh,τ)+b(τ,Ihu-uh)=b(τ,Ihu-u)∀τ∈Σh,(47)b(Πhσ-σh,v)=0∀v∈Vh.
The error estimates for interpolation operators Πh and Ih are summarized in the following lemma (cf. [21, 26]).Lemma 3.
For allv∈Hm+1(K,ℝd), τ∈Hm(Ω,𝕊) with m a positive integer and all K∈𝒯h, one has the estimates
(48)∥v-Ihv∥0,K+hK|v-Ihv|1,K+hK1/2∥v-Ihv∥0,∂K≲hKmin{m,k}+1∥v∥m+1,K,∥∇·(v-Ihv)∥0,K≲hKmin{m,k}∥∇·v∥m,K,∥τ-Πhτ∥0,K+hK1/2|τ-Πhτ|0,∂K+hK|τ-Πhτ|1,K≲hKmin{m,k+2}∥τ∥m,ωK,
where ωK:=⋃‍{K-′∈𝒯h;K-′∩K-≠∅}.Lemma 4.
For allτ∈Σh and v∈H0(div,Ω)∩V, it follows
(49)b(τ,v)≤∥τ∥𝒜(∥ε(v)∥𝒜-1+2μ∥r(⟦v⟧)∥0).Proof .
Applying integration by parts,(50)b(τ,v)=-∫Ω‍τ:ε(v)dx+∫ℰh‍{τ}:⟦v⟧vds=-∫Ω‍τ:(ε(v)+r(⟦v⟧))dx.
Then using (35) and triangle inequality,
(51)b(τ,v)≤∥τ∥𝒜(∥ε(v)∥𝒜-1+∥r(⟦v⟧)∥𝒜-1).
On the other hand, according to (44) and the definitions of trace operator tr and lifting operator r,
(52)∥tr(r(⟦v⟧))∥02=∫Ω‍tr(r(⟦v⟧))δ:r(⟦v⟧)dx=-∫ℰh‍{tr(r(⟦v⟧))δ}:⟦v⟧ds.
Since v∈H0(div,Ω) implies [v]=0 on ℰh,
(53)∥tr(r(⟦v⟧))∥02=-∫ℰh‍tr({r(⟦v⟧)})δ:⟦v⟧ds=-∫ℰh‍tr({r(⟦v⟧)})[v]ds=0,
which means tr(r(⟦v⟧))=0 in Ω. Thus by the definition of 𝒜-1,
(54)∥r(⟦v⟧)∥𝒜-12=2μ∥r(⟦v⟧)∥02+λ∥tr(r(⟦v⟧))δ∥02=2μ∥r(⟦v⟧)∥02.
Therefore, we can finish the proof by combining last inequality and (51).Lemma 5.
Assume that the solution(σ,u) of problem (1) is (s+2)-regular, where s is a nonnegative integer. For all τ∈Σh, it follows
(55)b(τ,Ihu-u)≲hmin{k,s+1}∥τ∥𝒜∥f∥s.Proof.
It follows from the definition of𝒜-1, (27), Lemma 3, and regularity (43) that
(56)∥ε(u-Ihu)∥𝒜-12+2μ∥r(⟦u-Ihu⟧)∥02=2μ(∥ε(u-Ihu)∥02+∥r(⟦u-Ihu⟧)∥02)+λ∥∇·(u-Ihu)∥02≲2μ(∥ε(u-Ihu)∥02+∑e∈ℰh‍he-1∥⟦u-Ihu⟧∥0,e2)+λ∥∇·(u-Ihu)∥02≲h2min{k,s+1}(∥u∥s+22+λ∥∇·u∥s+12)≲h2min{k,s+1}∥f∥s2.
Noting the fact that Ihu-u∈H0(div,Ω)∩V, we have from Lemma 4(57)b(τ,Ihu-u)≤∥τ∥𝒜(2μ∥ε(u-Ihu)∥𝒜-1+2μ∥r(⟦u-Ihu⟧)∥0)≲hmin{k,s+1}∥τ∥𝒜∥f∥s,
which completes the proof.With previous preparation, we are ready to give the a priori error estimates.Theorem 6.
Assume that the solution(σ,u) of problem (1) is (s+2)-regular, where s is a nonnegative integer. Let (σh,uh) be the solution of RLDG method (12)-(13). Then
(58)∥σ-σh∥𝒜≲hmin{k,s+1}∥f∥s.Proof .
Choosingτ=Πhσ-σh in (46) and v=Ihu-uh in (47) and subtracting (47) from (46), we get
(59)a(σ-σh,Πhσ-σh)=b(Πhσ-σh,Ihu-u).
Then we have
(60)∥Πhσ-σh∥𝒜2=a(Πhσ-σh,Πhσ-σh)=a(Πhσ-σ,Πhσ-σh)+a(σ-σh,Πhσ-σh)=a(Πhσ-σ,Πhσ-σh)+b(Πhσ-σh,Ihu-u).
Applying Cauchy-Schwarz inequality and Lemma 5, it follows
(61)∥Πhσ-σh∥𝒜2≲∥Πhσ-σ∥𝒜∥Πhσ-σh∥𝒜+hmin{k,s+1}∥Πhσ-σh∥𝒜∥f∥s,
that is,
(62)∥Πhσ-σh∥𝒜≲∥σ-Πhσ∥𝒜+hmin{k,s+1}∥f∥s.
By triangle inequality and (34), we get
(63)∥σ-σh∥𝒜≤∥σ-Πhσ∥𝒜+∥Πhσ-σh∥𝒜≲∥σ-Πhσ∥0+hmin{k,s+1}∥f∥s.
The proof is finished from Lemma 3.Under the assumption (25), we can get the following error estimate for uh in energy norm from Theorem 6 and (29).Corollary 7.
Assume that (25) is true and the solution (σ,u) of problem (1) is (s+2)-regular, where s is a nonnegative integer. Let (σh,uh) be the solution of RLDG method (12)-(13). Then
(64)∥εh(u-uh)+r(⟦u-uh⟧)∥𝒜-1=∥σ-σh∥𝒜≲hmin{k,s+1}∥f∥s.Using the usual duality argument, we can additionally derive an optimal order error estimate betweenuh and u in L2(Ω,ℝd)-norm.Theorem 8.
LetΩ be a convex bounded polygonal domain. Assume that the solution (σ,u) of problem (1) is (s+2)-regular, where s is a nonnegative integer. Let (σh,uh) be the solution of RLDG method (12)-(13). Then
(65)∥u-uh∥0≲hmin{k+1,s+2}∥f∥s.Proof .
Let(σ~,u~) be the solution of the auxiliary problem it follows:
(66)𝒜σ~-ε(u~)=0inΩ,-∇·σ~=u-uhinΩ,u~=0on∂Ω.
For Ω is a convex bounded polygon or polyhedron, the solution (σ~,u~) of problem (66) has 2-regularity estimate (cf. [24]):
(67)∥σ~∥1+∥u~∥2+λ∥∇·u~∥1≲∥u-uh∥0.
Multiplying u-uh on both sides of the second equation of (66) and integrating on Ω, we have
(68)∥u-uh∥02=-b(σ~,u-uh)=-b(σ~-Πhσ~,u-uh)-b(Πhσ~,u-uh).
By Lemma 2, we get
(69)b(σ~-Πhσ~,u-uh)=b(σ~-Πhσ~,u-Ihu).
Choosing τ=Πhσ~ in error equation (41), it holds
(70)-b(Πhσ~,u-uh)=a(σ-σh,Πhσ~)=a(σ-σh,Πhσ~-σ~)+a(σ-σh,σ~).
Then it follows from the last three equalities that
(71)∥u-uh∥02=-b(σ~-Πhσ~,u-Ihu)+a(σ-σh,Πhσ~-σ~)+a(σ-σh,σ~).
For the first term of the right hand side in (71), applying Cauchy-Schwarz inequality, Lemma 3, and regularity (67), it holds
(72)-b(σ~-Πhσ~,u-Ihu)=-∑K∈𝒯h‍∫K‍(∇·(σ~-Πhσ~))·(u-Ihu)dx+∫ℰhi‍[σ~-Πhσ~]·{u-Ihu}ds≲hmin{k+1,s+2}∥σ~∥1∥u∥s+2≲hmin{k+1,s+2}∥f∥s∥u-uh∥0.
For the second term of the right hand side in (71), we get from Cauchy-Schwarz inequality, (34), Theorem 6, Lemma 3, and regularity (67)
(73)a(σ-σh,Πhσ~-σ~)≤∥σ-σh∥𝒜∥Πhσ~-σ~∥𝒜≤∥σ-σh∥𝒜∥Πhσ~-σ~∥0≲hmin{k+1,s+2}∥f∥s∥σ~∥1≲hmin{k+1,s+2}∥f∥s∥u-uh∥0.
For the third term of the right hand side in (71), by the first equation of (66), integration by parts, and using error equation (42), we get
(74)a(σ-σh,σ~)=∫Ω‍(σ-σh):ε(u~)dx=-b(σ-σh,u~)=-b(σ-σh,u~-Ihu~)=-b(σ-Πhσ,u~-Ihu~)-b(Πhσ-σh,u~-Ihu~).
According to Lemma 5 with regularity (67), triangle inequality, (34), Theorem 6, and Lemma 3, we have
(75)-b(Πhσ-σh,u~-Ihu~)≲h∥Πhσ-σh∥𝒜∥u-uh∥0≲h(∥σh-σ∥𝒜+∥σ-Πhσ∥0)∥u-uh∥0≲hmin{k+1,s+2}(∥f∥s+∥σ∥s+1)∥u-uh∥0≲hmin{k+1,s+2}∥f∥s∥u-uh∥0.
On the other hand, using Cauchy-Schwarz inequality, Lemma 3, and regularity (67),
(76)-b(σ-Πhσ,u~-Ihu~)=-∑K∈𝒯h∫K(∇·(σ-Πhσ))·(u~-Ihu~)dx+∫ℰhi[σ-Πhσ]·{u~-Ihu~}ds≲hmin{k+1,s+2}∥σ∥s+1∥u~∥2≲hmin{k+1,s+2}∥f∥s∥u-uh∥0.
Thus we know from the last two inequalities that
(77)a(σ-σh,σ~)≲hmin{k+1,s+2}∥f∥s∥u-uh∥0.
Combining (71)–(77) gives
(78)∥u-uh∥02≲hmin{k+1,s+2}∥f∥s∥u-uh∥0.
Finally, we can conclude the result by cancelling ∥u-uh∥0 on both sides.


## 4. Postprocessing of RLDG Solution

In this section, we will focus on reconstructing an approximate displacement in the Raviart-Thomas element space via the ideas in [17]. Assume that triangulation 𝒯h is conforming in this section; that is, there is no hanging nodes in 𝒯h. The Raviart-Thomas element space is defined by (cf. [26, 27])
(79)RTh≔{v∈H0(div,Ω);v|K∈Pk(K,ℝd)+xPk(K)∀K∈𝒯h}.
With this space, the postprocessing can be described as follows. Find u-h∈RTh such that for any K∈𝒯h,
(80)∫K‍(u-h-uh)·vdx=0∀v∈Pk-1(K,ℝd),∫e‍(u-h-{uh})·nKvds=0∀v∈Pk(e,ℝd),∀facee⊂∂K∖∂Ω,
where nK is the unit outward normal on ∂K, pointing at outside of K. Since (80) are defined locally, the function u-h can be obtained very conveniently and in parallel. To show its approximation property, we first establish a relation between the divergence of the post-processed displacement u-h and the RLDG displacement uh. For this, let
(81)Vh:={v∈L2(Ω);v|K∈Pk(K)∀K∈𝒯h}.
And Qh denote the L2-projection operator from L2(Ω) onto the finite element space Vh. We first recall two results which have been proved in [28]
(82)∇·u-h=∇h·uh+Qh(tr(r(⟦uh⟧))),(83)∥u-h-uh∥0,K≲hK1/2∥tr(⟦uh⟧)∥0,∂K=hK1/2∥[uh]∥0,∂K∀K∈𝒯h.Now we can give the error estimate foru-h as follows.Theorem 9.
Let(σh,uh) be the solution of RLDG method (12)-(13) and let u-h be the solution of (80). Under the same assumption as in Corollary 7, one has
(84)λ∥∇·(u-u-h)∥0≲hmin{k,s+1}∥f∥s.
Moreover, if Ω is a convex bounded polygon or polyhedron,
(85)∥u-u-h∥0≲hmin{k+1,s+2}∥f∥s.Proof .
From (82), and the definition of Qh, it follows that
(86)λ∥∇·(u-u-h)∥0=λ∥∇h·(u-uh)-Qh(tr(r(⟦uh⟧)))∥0=λ∥∇h·(u-uh)+Qh(tr(r(⟦u-uh⟧)))∥0=λ∥∇·u-Qh(∇·u)+Qh(∇h·(u-uh)+tr(r(⟦u-uh⟧)))∥0.
Then by triangle inequality, and the boundedness of Qh (cf. [22, 23]), it follows that
(87)λ∥∇·(u-u-h)∥0=λ∥∇·u-Qh(∇·u)∥0+λ∥Qh(∇h·(u-uh)+tr(r(⟦u-uh⟧)))∥0=λ∥∇·u-Qh(∇·u)∥0+λ∥∇h·(u-uh)+tr(r(⟦u-uh⟧))∥0≤λ∥∇·u-Qh(∇·u)∥0+∥εh(u-uh)+r(⟦u-uh⟧)∥𝒜-1.
On the other hand, it is well known that L2-projection operator Qh has the following estimate (cf. [22, 23]):
(88)∥∇·u-Qh(∇·u)∥0≲hmin{k+1,s+1}∥∇·u∥s+1.
Hence, we obtain from the last two inequalities and Corollary 7 that
(89)λ∥∇·(u-u-h)∥0≲λhmin{k+1,s+1}∥∇·u∥s+1+hmin{k,s+1}∥f∥s≲hmin{k+1,s+1}∥f∥s+hmin{k,s+1}∥f∥s≲hmin{k,s+1}∥f∥s.
This completes the proof of (84).
As for (85), noting the fact that [Ihu]=0 for any e∈ℰh, we get from (83) and inverse inequality
(90)∥u-h-uh∥0,K≲hK1/2∥[Ihu-uh]∥0,∂K≲∥Ihu-uh∥0,ωK.
Summing over all elements in 𝒯h, it gives
(91)∥u-h-uh∥0≲∥Ihu-uh∥0.
Then using triangle inequality,
(92)∥u-u-h∥0≤∥u-uh∥0+∥u-h-uh∥0≲∥u-uh∥0+∥Ihu-uh∥0≤∥u-uh∥0+∥u-Ihu∥0.
Therefore, we can obtain (85) from Theorem 8 and Lemma 3.


## 5. Numerical Results

In this section, some numerical results are shown to validate the theoretical convergence rates of RLDG method. LetΩ=(-1,1)×(-1,1), μ=1, and
(93)f(x1,x2)=(-8(x1+x2)((3x1x2-2)(x12+x22)+5(x1x2-1)2-2x12x22)-8(x1-x2)((3x1x2+2)(x12+x22)-5(x1x2+1)2+2x12x22)).
It is not difficult to check that the exact solution of (1) for this load is
(94)u(x1,x2)=(-4x2(1-x22)(1-x12)24x1(1-x12)(1-x22)2)+12+λ(-4x1(1-x12)(1-x22)2-4x2(1-x22)(1-x12)2),σ(x1,x2)=2ε(u)+λ(∇·u)δ.We make a uniform triangulation𝒯h on Ω. Set SK=Pk+1(K,𝕊) for all K∈𝒯h and k=1.First, let us consider the convergence of RLDG method for fixed Lamé coefficientλ. Errors between (σh,uh) and (σ,u) for λ=1 are listed in Table 1, from which we can see that the convergence rates of ∥u-uh∥0 and ∥σ-σh∥𝒜 are O(h2) and O(h), respectively. These numerical results agree with theoretical results in Theorems 6 and 8. We also present the numerical convergence rates of postprocessing solution u-h in Table 2 for λ=1. As Theorem 9 indicates that the convergence rates of ∥u-u-h∥0 and λ∥∇·(u-u-h)∥0 are O(h2) and O(h), respectively.  To illustrate the robustness of RLDG method with respect to λ, we give numerical results for λ=109 in Table 3, which is very close to the incompressible case. It can be observed from Table 3 that the convergence rates of ∥u-uh∥0 and ∥σ-σh∥𝒜 being again O(h2) and O(h) do not deteriorate even for such large λ; hence our RLDG method is locking-free.Table 1
Errors between (σh, uh) and (σ, u) for different h when λ=1.
h
∥
u
-
u
h
∥
0
Order
∥
σ
-
σ
h
∥
𝒜
Order
1
1.0305
E
+
00
—
5.0397
E
+
00
—
2
-
1
5.1927
E
-
01
0.99
3.8102
E
+
00
0.40
2
-
2
1.8366
E
-
01
1.50
2.3137
E
+
00
0.72
2
-
3
5.3527
E
-
02
1.78
1.2586
E
+
00
0.88
2
-
4
1.4183
E
-
02
1.92
6.5126
E
-
01
0.95
2
-
5
3.6238
E
-
03
1.97
3.3031
E
-
01
0.98
2
-
6
9.1385
E
-
04
1.99
1.6620
E
-
01
0.99Table 2
Errors betweenu-h and u for different h when λ=1.
h
∥
u
-
u
-
h
∥
0
Order
λ
∥
∇
·
(
u
-
u
-
h
)
∥
0
Order
1
1.0762
E
+
00
—
1.4544
E
+
00
—
2
-
1
5.3145
E
-
01
1.02
8.9662
E
-
01
0.40
2
-
2
1.8647
E
-
01
1.51
6.0571
E
-
01
0.72
2
-
3
5.4186
E
-
02
1.78
3.5225
E
-
01
0.88
2
-
4
1.4342
E
-
02
1.92
1.8690
E
-
01
0.95
2
-
5
3.6625
E
-
03
1.97
9.5603
E
-
02
0.98
2
-
6
9.2341
E
-
04
1.99
4.8248
E
-
02
0.99Table 3
Errors between(σh, uh) and (σ, u) for different h when λ=109.
h
∥
u
-
u
h
∥
0
Order
∥
σ
-
σ
h
∥
𝒜
Order
1
1.1624
E
+
00
—
5.2687
E
+
00
—
2
-
1
7.1547
E
-
01
0.70
4.2548
E
+
00
0.31
2
-
2
2.7852
E
-
01
1.36
2.6741
E
+
00
0.67
2
-
3
8.4169
E
-
02
1.73
1.4674
E
+
00
0.87
2
-
4
2.2530
E
-
02
1.90
7.5886
E
-
01
0.95
2
-
5
5.7907
E
-
03
1.96
3.8401
E
-
01
0.98Then, we fixh=2-4 and let λ vary from 1 to 109. The corresponding errors ∥u-uh∥0, ∥σ-σh∥𝒜, ∥u-u-h∥0, and λ∥∇·(u-u-h)∥0 are shown in Table 4. From Table 4, we note that the errors ∥u-uh∥0, ∥σ-σh∥𝒜, and ∥u-u-h∥0 are insensitive to the choice of λ. And λ∥∇·(u-u-h)∥0 becomes smaller as λ grows larger; thus it is of course bounded. In fact, the quantity λ∥∇·(u-u-h)∥0 is insensitive to the choice of λ by an in-depth observation. In a word, all the errors in Table 4 are uniform with respect to λ, which coincide with the theoretical results.Table 4
Numerical results for differentλ when h=2-4.
λ
∥
u
-
u
h
∥
0
∥
σ
-
σ
h
∥
𝒜
∥
u
-
u
-
h
∥
0
λ
∥
∇
·
(
u
-
u
-
h
)
∥
0
1
1.4183
E
-
02
6.5126
E
-
01
1.4342
E
-
02
1.8690
E
-
01
10
1.8877
E
-
02
7.0622
E
-
01
1.8957
E
-
02
2.0032
E
-
01
1
0
2
2.1407
E
-
02
7.3973
E
-
01
2.1447
E
-
02
9.0172
E
-
02
1
0
3
2.2021
E
-
02
7.4962
E
-
01
2.2048
E
-
02
5.0191
E
-
02
1
0
4
2.2380
E
-
02
7.5621
E
-
01
2.2390
E
-
02
3.5219
E
-
02
1
0
5
2.2504
E
-
02
7.5852
E
-
01
2.2505
E
-
02
1.4447
E
-
02
1
0
6
2.2520
E
-
02
7.5883
E
-
01
2.2521
E
-
02
4.7216
E
-
03
1
0
7
2.2522
E
-
02
7.5886
E
-
01
2.2522
E
-
02
1.4982
E
-
03
1
0
8
2.2521
E
-
02
7.5886
E
-
01
2.2521
E
-
02
4.7392
E
-
04
1
0
9
2.2530
E
-
02
7.5886
E
-
01
2.2530
E
-
02
1.4987
E
-
04In the end, let us compare RLDG method and LDG method. Setλ=1 again. Errors ∥u-uh∥0 and ∥σ-σh∥𝒜 of RLDG and LDG methods with various choices of C11 are shown in Tables 5 and 6, from which we can see that errors ∥u-uh∥0 and ∥σ-σh∥𝒜 increase as C11 becomes larger. To be more precise, errors ∥u-uh∥0 and ∥σ-σh∥𝒜 grow slowly with respect to C11 when C11|e≲he-1 for all e∈ℰh; however, they become much larger when C11 reaches O(h-2) especially for finer mesh. We use PCG method with diagonal preconditioner to solve the resulting linear systems, and the stop criterion requires that the relative residual is less than 10-10. The numbers of diagonal-PCG iterations for RLDG and LDG methods are given in Table 7, which show similar variation as those in Tables 5 and 6. Therefore, we can conclude that RLDG method performs better than LDG method, especially for large C11. The performance of LDG method is almost robust with respect to C11 when C11≲h-1 and deteriorates when C11=O(h-2). In addition, if stress tensor and displacement are both approximated by piecewise linear polynomials, RLDG method will be ill-posed, whereas LDG method is still well-posed whose numerical results can be found in [20].Table 5
Error∥u-uh∥0 of RLDG and LDG methods with various choices of C11.
h
RLDG
LDG
C
11
|
e
=
h
e
C
11
|
e
=
1
C
11
|
e
=
h
e
-
1
C
11
|
e
=
h
e
-
2
1
1.0305
E
+
00
1.0694
E
+
00
1.0676
E
+
00
1.0663
E
+
00
1.0654
E
+
00
2
-
1
5.1927
E
-
01
5.2359
E
-
01
5.2737
E
-
01
5.3452
E
-
01
5.4783
E
-
01
2
-
2
1.8366
E
-
01
1.8401
E
-
01
1.8498
E
-
01
1.8856
E
-
01
2.0068
E
-
01
2
-
3
5.3527
E
-
02
5.3551
E
-
02
5.3707
E
-
02
5.4837
E
-
02
6.1566
E
-
02
2
-
4
1.4183
E
-
02
1.4185
E
-
02
1.4206
E
-
02
1.4515
E
-
02
1.7520
E
-
02
2
-
5
3.6238
E
-
03
3.6239
E
-
03
3.6267
E
-
03
3.7068
E
-
03
4.8648
E
-
03
2
-
6
9.1385
E
-
04
9.1386
E
-
04
9.1421
E
-
04
9.3461
E
-
04
1.3268
E
-
03Table 6
Error∥σ-σh∥𝒜 of RLDG and LDG methods with various choices of C11.
h
RLDG
LDG
C
11
|
e
=
h
e
C
11
|
e
=
1
C
11
|
e
=
h
e
-
1
C
11
|
e
=
h
e
-
2
1
5.0397
E
+
00
5.0404
E
+
00
5.0384
E
+
00
5.0372
E
+
00
5.0364
E
+
00
2
-
1
3.8102
E
+
00
3.8105
E
+
00
3.8112
E
+
00
3.8138
E
+
00
3.8224
E
+
00
2
-
2
2.3137
E
+
00
2.3137
E
+
00
2.3138
E
+
00
2.3155
E
+
00
2.3331
E
+
00
2
-
3
1.2586
E
+
00
1.2586
E
+
00
1.2587
E
+
00
1.2595
E
+
00
1.2853
E
+
00
2
-
4
6.5126
E
-
01
6.5126
E
-
01
6.5126
E
-
01
6.5168
E
-
01
6.8287
E
-
01
2
-
5
3.3031
E
-
01
3.3031
E
-
01
3.3031
E
-
01
3.3052
E
-
01
3.6152
E
-
01
2
-
6
1.6620
E
-
01
1.6620
E
-
01
1.6620
E
-
01
1.6630
E
-
01
1.9153
E
-
01Table 7
The numbers of diagonal-PCG iterations for RLDG and LDG methods with various choices ofC11.
h
RLDG
LDG
C
11
|
e
=
h
e
C
11
|
e
=
1
C
11
|
e
=
h
e
-
1
C
11
|
e
=
h
e
-
2
1
27
27
26
25
25
2
-
1
64
66
66
65
67
2
-
2
128
129
130
131
139
2
-
3
247
247
247
254
295
2
-
4
473
473
477
502
670
2
-
5
924
924
924
992
1589
2
-
6
1848
1848
1848
1905
3844


## 6. Conclusions

In this paper, we develop a locking-free RLDG method for nearly incompressible linear elasticity. RLDG method can be formally regarded as a special case of LDG method withC11=0 proposed in [20] and can also be viewed as the localization of symmetric nonconforming mixed finite element method presented in [21]. Then we discuss the implementation of RLDG method and obtain robust error estimates in energy norm and L2-norm which are uniform with respect to λ. Error estimates of the postprocessing H(div)-nonconforming displacement in H(div)-norm and L2-norm are achieved, which are also uniform with respect to λ. According to numerical results, RLDG method performs a little better than LDG method for small C11 and much better than LDG method for large C11. The performance of LDG method is almost robust with respect to C11 when C11≲h-1 and deteriorates when C11=O(h-2).



---

*Source: 546408-2013-03-31.xml*