# Stabilized Multiscale Nonconforming Finite Element Method for the Stationary Navier-Stokes Equations

Authors: Tong Zhang; Shunwei Xu; Jien Deng
Journal: Abstract and Applied Analysis
 (2012)

---

## Abstract

We consider a stabilized multiscale nonconforming finite element
method for the two-dimensional stationary incompressible Navier-Stokes problem.
This method is based on the enrichment of the standard polynomial space for the velocity
component with multiscale function and the nonconforming lowest equal-order
finite element pair. Stability and existence uniqueness of the numerical solution are
established, optimal-order error estimates are also presented. Finally, some numerical
results are presented to validate the performance of the proposed method.


---

## Body

## 1. Introduction

As the development of science and technology, finite element method has became an important and powerful tool for the complex fluid problems, such as for the Navier-Stokes equations. It is well known that the pressure and velocity pairs satisfy the discrete Inf-Sup condition [1] that plays the key role for simulating the Navier-Stokes equations. However, some unstable mixed finite element pairs which violate the so-called Inf-Sup condition are also popular, see [2–4]. In order to overcome this restriction, various of stabilized methods have been proposed, including the bubble condensation-based methods [5], pressure projection method (PPM) [6–8], the local Gauss integration method (LGIM) [9–11], multiscale method [12, 13], macroelement stabilized method [3, 14], and so on. Most of these stabilized methods necessarily need to introduce the stabilization parameters either explicitly or implicitly. In addition, some of these techniques are conditionally stable or are of suboptimal accuracy. Therefore, the development of mixed finite element methods free from stabilization parameters has become increasingly important.In 2005, Franca et al. gave a new multiscale method for the reaction-diffusion equation in [15]. The chief characteristic of their method is to use the Petrov-Galerkin approach to split the solution into two parts, and the trial function space is enriched with an unstable bubble-like function, which is the solution to a local problem. Later, Barrenechea and Valentin [13] considered the relationship between the enriched multiscale method and stabilized techniques for generalized Stokes problem based on the P1-P1 pair. By enriching the velocity space with an unusual bubble function, Araya et al. established the convergence for the Stokes problem in [16], their method is different from usual residual free bubble method in [5], in which one should choose local basis functions to enrich the standard finite element spaces by solving some local problem analytically. Furthermore, the method proposed in [15] can also be used to treat the unsteady reaction-diffusion problem (see [17]).Compared with conforming finite element method, the nonconforming finite element methods are more popular due to their simplicity and small support sets of basis functions. Crouzeix and Ravizrt in [18] used the nonconforming piecewise linear velocity and a piecewise constant pressure to solve the Stokes equations. In this paper, motivated by the ideas of [13, 15, 16], we will use the Petrov-Galerkin approach based on the nonconforming velocity space to handle with the steady Navier-Stokes equations. The main differences between [13, 15, 16, 19] and this work lie in the following: (i) the finite element spaces of velocity are different; it is nonconforming element in this paper; (ii) the treated problems are different; we consider the nonlinear problem; (iii) the finite element pairs are different; the NCP1-P1 pair is used in this paper.The outline of this work is arranged as follows. In the following section, the abstract functional setting for steady Navier-Stokes equations is recalled. Section3 is devoted to derive the general form of enriched multiscale method based on the NCP1-P1 pair. After providing the stability and existence uniqueness for the approximation solution, the optimal error estimates are established in Section 4. In Section 5, Some numerical results are presented to verify the established theoretical analysis. Finally, Some conclusions are made in Section 6.


## 2. Preliminaries

LetΩ be an open bounded domain of ℝ2 with Lipschitz continuous boundary ∂Ω and satisfy a further condition stated in (A1) below. The incompressible stationary Navier-Stokes equations with the homogeneous Dirichlet boundary condition are
(2.1)-νΔu+(u·∇)u+∇p=finΩ,divu=0,inΩ,u|∂Ω=0on∂Ω,
where u=(u1(x),u2(x))T represents the velocity, p=p(x) the pressure, f=f(x)∈L2(Ω)2 the prescribed body force, ν>0 the viscosity coefficient.In order to introduce the variational formulation for problem (2.1), we set
(2.2)X=H01(Ω)2,Y=L2(Ω)2,D(A)=H2(Ω)2∩X,M=L02(Ω)={q∈L2(Ω):∫Ω‍qdx=0}.The standard notations of Sobolev spaceWm,r(Ω) are used. To simplify, we use Hm(Ω) instead of Wm,r(Ω) as r=2 and ∥·∥m for ∥·∥m,2. The spaces L2(Ω)m(m=1,2) are endowed with the usual L2-scalar product (·,·) and L2-norm ∥·∥0. The spaces H01(Ω) and X are equipped with the scalar product (∇u,∇v) and the norm |u|1,Ω2, u, v∈H01(Ω)(orX).DefineAu=-Δu is the operator associated with the Navier-Stokes problem, it is positive self-adjoint operator from D(A) onto Y.Introducing the bilinear operator(2.3)B(u,v)=(u·∇)v+12(divu)v∀u,v∈X,
and defining a trilinear form on X×X×X as follows:
(2.4)b(u,v,w)=〈B(u,v),w〉X′×X=((u·∇)v,w)+12((divu)v,w)=12((u·∇)v,w)-12((u·∇)w,v).The variational formulation of problem (2.1) reads as: find (u,p)∈(X,M) such that for all (v,q)∈(X,M)(2.5)a(u,v)-d(v,p)+d(u,q)+b(u,u,v)=(f,v),
where
(2.6)a(u,v)=ν(∇u,∇v),d(v,q)=-(v,∇q)=(q,divv),B0((u,p);(v,q))=a(u,v)-d(v,p)+d(u,q).Clearly, the bilinear formsa(·,·) and d(·,·) are continuous on X×X and X×M, respectively. Moreover, d(·,·) also satisfies (see [20]):
(2.7)sup0≠v∈X|d(v,q)||v|1,Ω≥β∥q∥0,Ω,
where β is a positive constant depending only on Ω.It is easy to verify thatB0 satisfies the following important properties for all (u,p), (v,q)∈(X,M) (see [1]):
(2.8)B0((u,p);(u,p))=ν∥u∥1,Ω2,(2.9)|B0((u,p);(v,q))|≤C(∥u∥1,Ω+∥p∥0,Ω)(∥v∥1,Ω+∥q∥0,Ω),(2.10)β0(∥u∥1,Ω+∥p∥0,Ω)≤sup(v,q)∈(X,M)|B0((u,p);(v,q))|∥v∥1,Ω+∥q∥0,Ω,
where β0>0 is a constant. Here and below, the letter C (with or without subscript) denotes a generic positive constant, depending at most on the data ν, Ω and f. Furthermore, the following estimates about b(·,·,·) are hold [1, 20]:
(2.11)b(u,v,w)=-b(u,w,v),(2.12)|b(u,v,w)|≤12c0∥u∥0,Ω1/2|u|1,Ω1/2(|v|1,Ω∥w∥0,Ω1/2|w|1,Ω1/2+∥v∥0,Ω1/2|v|1,Ω1/2|w|1,Ω),
for all u, v, w∈X and
(2.13)|b(u,v,w)|+|b(v,u,w)|+|b(w,u,v)|≤C|u|1,Ω∥Av∥0,Ω∥w∥0,Ω,
for all u∈X,v∈D(A), w∈Y.As mentioned above, a further assumption aboutΩ is needed (see [1]).(
A
1
) Assume that Ω is regular so that the unique solution (v,q)∈(X,M) of the steady Stokes equations
(2.14)-νΔv+∇q=g,divv=0inΩ,v|∂Ω=0,
for a prescribed g∈Y exists and satisfies
(2.15)∥Av∥0,Ω+|q|1,Ω≤C∥g∥0,Ω.Under the assumption of(A1), if ∂Ω is of C2 or Ω is a two-dimensional convex polygon, it has been shown that (see [20])
(2.16)∥v∥0,Ω≤γ0|v|1,Ω,∀v∈X,|v|1,Ω≤γ0∥Av∥0,Ω,∀v∈D(A),
where γ0 is a positive constant only depending on Ω.The following existence and uniqueness results about problem (2.5) are classical (see [1, 20]).Theorem 2.1.
Assume thatν and f∈Y satisfy the following uniqueness condition:
(2.17)1-c0γ02ν2∥f∥0,Ω>0.
Then problem (2.5) admits a unique solution (u,p)∈(D(A),H1(Ω)∩M) with divu=0 such that
(2.18)|u|1,Ω≤γ0ν∥f∥0,Ω,∥Au∥0,Ω+|p|1,Ω≤C∥f∥0,Ω.


## 3. Enriched Nonconforming Finite Element Method

Let𝒯h be a regular triangulation of Ω into element {Kj}:Ω¯=∪K¯j, that is, |Kj|≃ChKj2, where |Kj| is the area of the element Kj and hKj is the diameter of Kj; the mesh parameter h is given by h=max{hKj:Kj∈𝒯h}. Denote the boundary segment and the interior boundary by γj=∂Ω∩∂Kj and γjk=γkj=∂Kj∩∂Kk, respectively. Let Γh and Γj be the sets of γjk and γj. The centers of γj and γjk are indicated by ξj and ξjk, respectively. The finite element spaces investigated in this paper are the following mixed finite element spaces:
(3.1)NCP1={v∈Y:v|Kj∈P1(Kj)2,v(ξjk)=v(ξkj),v(ξj)=0∀j,k,Kj∈𝒯h},P1={q∈H1(Ω):q|Kj∈P1(Kj),∀Kj∈𝒯h},
where P1(Kj) is the set of line polynomials on Kj, and noting that the nonconforming finite element space NCP1 is not a subspace of X. Defining the energy norm (3.2)∥v∥1,h=(∑Kj‍|v|1,Kj2)1/2,∀v∈NCP1.
The finite element spaces NCP1 and P1 satisfy the following approximation property (see [4, 21]): for (v,q)∈H2(Ω)×H1(Ω), there are two approximations vI∈NCP1 and qI∈P1 such that
(3.3)∥v-vI∥0,Ω+h(∥v-vI∥1,h+∥q-qI∥0,Ω)≤Ch2(∥Av∥0,Ω+|q|1,Ω),
and the compatibility conditions hold for all j and k:
(3.4)∫γjk‍[v]ds=0,∫Γj‍vds=0∀v∈NCP1,
where [v]=vγjk-vγkj denotes the jump of the function v across the boundary γjk.Set〈·,·〉j=(·,·)∂Kj and |·|m,j=|·|m,Kj. Then for all u, v∈H1(Kj)2, q∈L2(Ω), the discrete bilinear forms are
(3.5)ah(u,v)=∑Kj‍ν(∇u,∇v)Kj,dh(v,p)=∑Kj‍(divv,p)Kj.
For the nonconforming space NCP1, we define a local operator
(3.6)Πj:H1(Kj)2→NCP1(Kj),
satisfying
(3.7)∫∂Kj‍(v-Πjv)ds=0.
Then the local operator Πj satisfies (see [21])
(3.8)|v-Πjv|1,Kj≤Chi|v|i+1,Kj,v∈Hi+1(Kj),i=0,1,∥Πjv∥1,Kj≤C∥v∥1,Kj.
The global operator Πh:X→NCP1 is defined as Πhv|j=Πjv,v∈X.As noted, the choiceNCP1-P1 is an unstable pair that does not satisfy the discrete Inf-Sup condition. Therefore, we need to introduce the enrichment multiscale method to overcome this restriction.LetEh be a finite dimensional space, called multiscale space, such that
(3.9)Eh∈H1(𝒯h)2,Eh∩NCP1={0},whereH1(𝒯h)2={v∈Y:v|Kj∈H1(Kj)2}.
The discrete weak formulation of the Stokes equations is to find uh+ue∈NCP1⊕Eh and ph∈P1, such that
(3.10)ah(uh+ue,vh)-dh(vh,ph)+dh(uh+ue,qh)=(f,v)Ω,
for all vh∈NCP1⊕H01(𝒯h)2 and qh∈P1. Let ue|Kj=ueKj+ue∂Kj, we can solve it through the following local problem:
(3.11)-νΔueKj=f+νΔuh-∇phinKj,ueKj|∂Kj=0,-νΔue∂Kj=0inKj,ue∂Kj=geon∂Kj,-ν∂ssge=1he[ν∂nuh+phI·n]E,ge=0atthenodes,
where he denotes the length of the edge E∈∂Kj; n the normal outward vector on ∂Kj;∂s,∂n are the tangential and normal derivative operators, respectively; I is the ℝ2×2 identity matrix. Equation (3.11) is well posed, that is, ue can be expressed by uh, ph, and f on each element Kj. For convenience, we define two local operators ℳKj:L2(Kj)2→H01(Kj)2 and ℋKj:L2(∂Kj)2→H1(Kj)2 by
(3.12)ueKj=1νℳKj(f+νΔuh-∇ph),∀Kj∈𝒯h,ue∂K=1νℋKj([ν∂nuh+phI·n]E),∀Kj∈𝒯h,E∈Γh.
With Green formulation and (3.12), for all (uh,ph), (vh,qh)∈NCP1×P1, (3.10) can be rewritten as
(3.13)∑Kj‍[ν(∇uh,∇vh)Kj-(ph,∇·vh)Kj+(qh,∇·uh)Kj]+∑Kj‍1ν(ℳKj(-νΔuh+∇ph)-ℋKj([ν∂nuh+phI·n]E),νΔvh+∇qh)Kj+∑E∈Γh‍1ν(ℋKj([ν∂nuh+phI·n]E),ν∂nvh+qhI·n)E=∑Kj‍[(f,vh)Kj+1ν(ℳKj(f),νΔvh+∇qh)Kj].
With the help of (3.13), the enriched nonconforming finite element method for the stationary Navier-Stokes equations (2.1) is rewritten as follows: find (uh,ph)∈NCP1×P1 such that
(3.14)B((uh,ph);(vh,qh))+b(uh,uh,vh)=F(vh,qh)
for all (vh,qh)∈NCP1×P1, where
(3.15)B((uh,ph);(vh,qh))=Bh((uh,ph);(vh,qh))+∑Kj‍1ν(ℳKj(∇ph),∇qh)Kj+∑E∈Γh‍1ν(ℋKj([ν∂nuh]),[ν∂nvh])E≜Bh((uh,ph);(vh,qh))Ω+∑Kj‍τKj(∇ph,∇qh)Kj+∑E∈Γh‍τE([ν∂nuh],[ν∂nvh])E,F(vh,qh)=∑Kj‍[(f,vh)Kj+1ν(ℳKj(f),∇qh)Kj],Bh((uh,ph);(vh,qh))=ah(uh,vh)-dh(vh,ph)+dh(uh,qh).By applying the technique to one used in [16], we can obtain that (bKj,1)Kj/|Kj|⋍C~hKj2, (aKj,1)E/he⋍he/12, τKj⋍C~hKj2 and τE⋍he/(12ν). Moreover, if f is a piecewise constant, then we have ℳKj(f)=bKjf,
(3.16)(ℳKj(f),∇qh)Kj=(bKj,1)Kj|Kj|(f,∇qh)Kj≃C~hKj2(f,∇qh)Kj.
Define the mesh-dependent norms as follows:
(3.17)∥|u|∥h2=ν∥u∥1,h2+∑E∈Γh‍τE∥[ν∂nu]∥0,E2,∥q∥h2=∑Kj‍τKj|q|1,Kj2.Remark 3.1.
The assumption of piecewise constantf is made simply to analyze the problem (3.14), but this assumption does not affect the precision of this method, and (3.14) may be implemented as it is presented for a general function f∈L2(Ω)2. Here, we do not give the detail proof about this fact; readers can visit Appendix  B of the paper [16] for f∈H1(Ω)2.Remark 3.2.
Generally speaking, the following linear algebra equations can be obtained from the discrete system of original problem:(3.18)(A-DDT0)(UP)=(F0),
where the matrices A and D are deduced from the diffusion, convection, and incompressible terms; F is the variation of the source term. The norm of matrix A gets smaller as the convection increases; therefore, some unnecessary oscillations will be created. In order to eliminate these oscillations, we introduce the stabilized term, in this case, the coefficient matrix of discrete formulation transforms into
(3.19)(A-DDTG),
where G is derived from the stabilized term, that is, the term of (∇ph,∇qh). As the considered problem has strong convection, in order to obtain a good behavior of matrix (A-DDTG), we should choose a proper G. In this way, the singularly perturbed problem can be treated effectively. The reason that we treat the convection term not use enriched function technique is to simply the theoretical analysis and computation, and the discrete convection term has no influence about the stabilized term G(·,·).Lemma 3.3.
Let(vh,qh)∈NCP1×P1, then,
(3.20)B((vh,qh);(vh,qh))=∥|vh|∥h2+∥qh∥h2.Proof.
The results follow from the definition of (3.15) and the mesh-dependent norms in each K∈𝒯h.
Before establishing the stability of scheme (3.14), we introduce the local trace theorem (see [1]). There exists C>0, independent of h, such that
(3.21)∥u∥0,∂Kj2≤C(hKj-1∥u∥0,Kj2+hKj|u|1,Kj2),∀u∈H1(Kj).Theorem 3.4.
There exist two positive constantsC, β depending on ν, for all (uh,ph), (vh,qh)∈NCP1×P1 such that
(3.22)|B((uh,ph);(vh,qh))|≤C(∥uh∥1,h+∥ph∥0,Ω)(∥vh∥1,h+∥qh∥0,Ω),(3.23)sup0≠(vh,qh)∈(NCP1,P1)|B((uh,ph);(vh,qh))|∥vh∥1,h+∥qh∥0,Ω≥β(∥uh∥1,h+∥ph∥0,Ω).Proof.
It follows from(uh,ph), (vh,qh)∈NCP1×P1, inverse inequality, (3.15), and (3.21) that
(3.24)|B((uh,ph);(vh,qh))|≤ν∥uh∥1,h∥vh∥1,h+∥uh∥1,h∥qh∥0,Ω+∥ph∥0,Ω∥vh∥1,h+∑Kj‍τKj|ph|1,Kj|qh|1,Kj+∑E∈Γh‍τE∥[ν∂nuh]∥0,E∥[ν∂nvh]∥0,E≤ν∥uh∥1,h∥vh∥1,h+∥uh∥1,h∥qh∥0,Ω+∥ph∥0,Ω∥vh∥1,h+C1∥ph∥0,Ω∥qh∥0,Ω+C2∑E∈Γh‍τE(hKj-1/2∥∇uh∥0,Kj+hKj1/2|∇uh|1,Kj)(hKj-1/2∥∇vh∥0,Kj+hKj1/2|∇vh|1,Kj)≤C(ν)∥uh∥1,h∥vh∥1,h+∥uh∥1,h∥qh∥0,Ω+∥ph∥0,Ω∥vh∥1,h+C1∥ph∥0,Ω∥qh∥0,Ω≤C(∥uh∥1,h+∥ph∥0,Ω)(∥vh∥1,h+∥qh∥0,Ω),
that is, the continuity result (3.22) holds.
From the properties of the nonconforming finite element given in [18], for all ph∈L2(Kj), there exists a function w∈H1(Kj)2, such that ∥w∥1,h=∥ph∥0 and
(3.25)(∇·w,ph)Kj=∥ph∥0,Kj2,∥w∥1,h≤C∥ph∥0,Ω.
Using the Cauchy-Schwartz inequality and (3.25), we have
(3.26)|B((uh,ph);(-w,0))|=-ν(∇uh,∇w)+(ph,∇·w)-∑E∈ΓhτE([ν∂nuh],[ν∂nw])E≥-ν∥uh∥1,h∥w∥1,h+C0∥w∥1,h∥ph∥0,Ω-∑E∈ΓhτE∥[ν∂nuh]∥0,E∥[ν∂nw]∥0,E≥-(ν∥uh∥1,h2+∑E∈ΓhτE∥[ν∂nuh]0,E2∥)1/2(ν∥w∥1,h2+∑E∈ΓhτE∥[ν∂nw]∥0,E2)1/2+C0∥w∥1,h∥ph∥0,h.
Using (3.21) and inverse inequality, we obtain that
(3.27)τE∥[ν∂nw]∥0,E2≤he12ν(hK-1∥ν∇w·n∥0,K2+hK|ν∇w·n|1,K2)≤heν12hK|w|1,K2+CKνhe12hK|w|1,K2≤ν(1+CK)12|w|1,K2.
Combining (3.26) with (3.27) yields
(3.28)|B((uh,ph);(-w,0))|≥-Cν∥w∥1,h(ν|uh|1,h2+∑E∈ΓhτE∥[ν∂nuh]∥0,E2)1/2+C0∥w∥1,h∥ph∥0,Ω=-Cν∥w∥1,h∥|uh|∥h+C0∥ph∥0,Ω2≥-Cνγ1-1∥|uh|∥h2+(C0-γ1)∥ph∥0,Ω2,
where C=(1+C0)/12 with C0=maxk∈𝒯hCK, and γ1 is chosen small enough. Let
(3.29)(vh,qh)=(uh-δw,ph),δ>0.
Using (3.26) and Lemma 3.3 we have
(3.30)|B((uh,ph);(vh,qh))|=|B((uh,ph);(uh,ph))+δB((uh,ph);(-w,0))|≥∥|uh|∥h2+∥ph∥h2+δ(-Cνγ1-1∥|uh|∥h2+(C0-γ1)∥ph∥0,Ω2)≥(1-Cδνγ1-1)∥|uh|∥h2+∥ph∥h2+δ(C0-γ1)∥ph∥0,Ω2≥ν(1-Cδνγ1-1)∥uh∥1,h2+δ(C0-γ1)∥ph∥0,Ω2,
provided that 0<δ<γ1/(Cν) and 0<γ1<C0. Denote
(3.31)C(ν)≜min{ν(1-Cδνγ1-1),δ(C0-γ1)},C(δ)≜max{2,1+2δ2}.
Then we have
(3.32)∥vh∥1,h2+∥qh∥0,Ω2=∥uh-δw∥1,h2+∥ph∥0,Ω2≤2∥uh∥1,h2+(1+2δ2)∥ph∥0,Ω2≤C(δ)(∥uh∥1,h2+∥ph∥0,Ω2).
Taking β=C(ν)/(C(δ)), we obtain the desired result (3.23).Theorem 3.5.
Under the assumptions of Theorem2.1 and the following condition:
(3.33)thestronguniquenesscondition:1-c0γ0γ0+ν1/2ν2∥f∥0,Ω>0.
Problem (3.14) admits a unique solution (uh,ph)∈(NCP1,P1), and satisfying
(3.34)∥uh∥1,h≤γ0+ν1/2ν∥f∥0,Ω,∥ph∥0,Ω≤β-1[(γ0+Cν)∥f∥0,Ω+c0γ0ν-2(γ0+ν1/2)2∥f∥0,Ω2].Proof.
Let Hilbert spaceHh=(NCP1,P1) be with the scalar product and norm
(3.35)((v,q);(w,r))Hh=∑Kj(∇v,∇w)Kj+(q,r),
and Kh be a nonvoid, convex, and compact subset of Hh defined by
(3.36)Kh={c0γ0(γ0+ν1/2)2βν2(v,q)∈Hh:∥v∥1,h≤γ0+ν1/2ν∥f∥0,∥q∥0,Ω≤νγ0+Cβν∥f∥0,Ω+c0γ0(γ0+ν1/2)2βν2∥f∥0,Ω2}.
Defining a continuous mapping fromKh into Hh as follows: given (v¯,q¯)∈Kh,for all (w,r)∈Hh, find (v,q)=F(v¯,q¯) such that
(3.37)B((v,q);(w,r))+b(v¯,v,w)=(f,w)+∑KjτKj(f,∇r)Kj.
Taking (w,r)=(v,q), using (2.8)–(2.13) and inverse inequality yields
(3.38)ν∥v∥1,h2+∥q∥h2≤∥|v|∥h2+∥q∥h2≤γ0∥f∥0,Ω∥v∥1,h+∑KjτKj∥f∥0,Kj∥∇q∥0,Kj≤(γ022ν+h22)∥f∥0,Ω2+ν2∥v∥1,h2+12∑KjτKj∥∇q∥0,Kj2=(γ022ν+h22)∥f∥0,Ω2+ν2∥v∥1,h2+12∥q∥h2.
As a consequence, we have
(3.39)∥v∥1,h≤γ0+ν1/2ν∥f∥0,Ω.
Using again (2.17), (3.23), (3.37), and inverse inequality, we arrive at
(3.40)β(∥v∥1,h+∥q∥0,Ω)≤|f,w|+|∑KjτKj(f,∇r)Kj|∥w∥1,h+∥r∥0+c0γ0∥v¯∥1,h∥v∥1,h≤γ0∥f∥0,Ω+Chν∥f∥0,Ω+c0γ0ν-2(γ0+ν1/2)2∥f∥0,Ω2≤(γ0ν+Chν)∥f∥0,Ω+c0γ0ν-2(γ0+ν1/2)2∥f∥0,Ω2.
Hence, the two estimates imply (v,q)=F(v¯,q¯)∈Kh, thanks to the fixed point theorem, the mapping (v,q)=F(v¯,q¯)∈Kh has at least one fixed point (uh,ph)∈Kh; namely, (uh,ph)∈Kh is a numerical solution of problem (3.14).
Next, we shall prove that the problem (3.14) has a unique solution (uh,ph). In fact, if (vh,qh) also satisfies (3.14), then for all (w,r)∈(NCP1,P1) we have
(3.41)B((uh-vh,ph-qh);(w,r))=b(vh-uh,uh,w)+b(vh,vh-uh,w).
Taking (w,r)=(uh-vh,ph-qh) in (3.41) and using again (2.8)–(2.13), Lemma 3.3, it follows that
(3.42)ν∥uh-vh∥1,h2≤c0γ0∥uh∥1∥uh-vh∥1,h2≤c0γ0γ0+ν1/2ν∥f∥0,Ω∥uh-vh∥1,h2,
Which, together with the strong uniqueness condition
(3.43)ν-c0γ0γ0+ν1/2ν∥f∥0,Ω=ν(1-c0γ0γ0+ν1/2ν2∥f∥0,Ω)>0,
gives uh=vh. Using again (2.13), (3.23), and (3.41), we obtain β∥ph-qh∥0,Ω2≤0 which implies ph=qh.


## 4. Error Estimates

In order to derive the error estimates of the numerical solution(uh,ph), we introduce the Galerkin projection (Rh,Qh):(X,M)→(NCP1,P1) defined as follows: for all (vh,qh)∈(NCP1,P1)(4.1)B((Rh(v,q),Qh(v,q));(vh,qh))=B0((v,q);(vh,qh)).
Noting the Theorem 3.4, (Rh(v,q),Qh(v,q)) is well defined.By using a similar argument to the one used in [14, 22], we have the following lemma.Lemma 4.1.
Let(u,p)∈D(A)×(H1(Ω)∩M); under the assumptions of Theorems 3.4 and 3.5, the projection operator (Rh,Qh) satisfies
(4.2)∥u-Rh(u,p)∥0,Ω+h(∥u-Rh(u,p)∥1,h+∥p-Qh(u,p)∥0,Ω)≤Ch2(∥Au∥0,Ω+|p|1,Ω).Proof.
From(u,p)∈[H2(Ω)2∩X]×[H1(Ω)∩M], we have [ν∂nu]E=0. For all (vh,qh)∈NCP1×P1, using (4.1) yields
(4.3)B((Rh(u,p),Qh(u,p))(vh,qh))=B0((u,p);(vh,qh))=B((u,p);(vh,qh))-∑KjτKj(∇p,∇qh)Kj.
From the definition of (Rh(u,p),Qh(u,p)), (3.3), combining Theorem 3.4, (4.3), the triangular with inverse inequalities, we arrive at
(4.4)∥u-Rh(u,p)∥1,h+∥p-Qh(u,p)∥0,Ω≤∥u-uI∥1,h+∥p-pI∥0,Ω+∥uI-Rh(u,p)∥1,h+∥pI-Qh(u,p)∥0,Ω≤∥u-uI∥1,h+∥p-pI∥0,Ω+β-1sup(vh,qh)∈NCP1×P1|B((uI-Rh(u,p),pI-Qh(u,p));(vh,qh))|∥vh∥1,h+∥qh∥0,Ω≤β-1sup(vh,qh)∈NCP1×P1|B((uI-u,pI-p);(vh,qh))|+|∑KjτKj(∇p,∇qh)Kj|∥vh∥1,h+∥qh∥0,Ω+∥u-uI∥1,h+∥p-pI∥0,Ω.
It is easy to check that
(4.5)|B((uI-u,pI-p);(vh,qh))|≤C(∥u-uI∥1,h+∥p-pI∥0,Ω+h|p|1)(∥vh∥1,h+∥qh∥0,Ω).
Combining (4.4), (4.5), and inverse inequality yields
(4.6)∥u-Rh(u,p)∥1,h+∥p-Qh(u,p)∥0,Ω≤C(∥u-uI∥1,h+∥p-pI∥0,Ω+h|p|1)+β-1sup(vh,qh)∈NCP1×P1|∑KjτKj(∇p,∇qh)Kj|∥vh∥1,h+∥qh∥0,Ω≤C1h(∥Au∥0,Ω+|p|1,Ω)+C2hsup(vh,qh)∈NCP1×P1∑KjhKj∥∇p∥0,Kj∥∇qh∥0,Kj∥vh∥1,h+∥qh∥0,Ω≤Ch(∥Au∥0,Ω+∥p∥1,Ω).
In order to derive the estimate in the L2-norm, we consider the following dual problem with (e,η)=(u-Rh(u,p),p-Qh(u,p)):
(4.7)-ΔΦ+∇Ψ=einΩ,(4.8)divΦ=0inΩ,(4.9)Φ|∂Ω=0on∂Ω.
Based on the assumption of (A1), (4.7)–(4.9) have a unique solution and satisfy
(4.10)∥AΦ∥0,Ω+|Ψ|1,Ω≤C∥u-Rh(u,p)∥0,Ω.
Multiplying (4.7) and (4.8) by e and η, respectively, integrating over Ω, and using (4.3) with (vh,qh)=(ΦI,ΨI), we see that
(4.11)∥u-Rh(u,p)∥02=B0((u-Rh(u,p),p-Qh(u,p));(Φ,Ψ))-∑Kj〈∂Φ∂n,u-Rh(u,p)〉j+∑Kj〈(u-Rh(u,p))·n,Ψ〉j=B((u-Rh(u,p),p-Qh(u,p));(Φ,Ψ))-∑KjτKj(∇(p-Qh(u,p)),∇Ψ)Kj-∑Kj〈∂Φ∂n,u-Rh(u,p)〉j+∑Kj〈(u-Rh(u,p))·n,Ψ〉j=B((u-Rh(u,p),p-Qh(u,p));(Φ-ΦI,Ψ-ΨI))-∑KjτKj(∇(p-Qh(u,p)),∇Ψ)Kj+∑KjτKj(∇p,∇ΨI)Kj-∑Kj〈∂Φ∂n,u-Rh(u,p)〉j+∑Kj〈(u-Rh(u,p))·n,Ψ〉j,
where (ΦI,ΨI) is the finite element interpolation of (Φ,Ψ) in (NCP1,P1) and satisfies (3.3). For each E∈∂Kj, we define the mean value of u-Rh(u,p) and Ψ on E(4.12)u-Rh(u,p)¯=1he∫E‍(u-Rh(u,p))|Kjds;Ψ¯=1he∫E‍Ψ|Kjds.
Note that each interior edge appears twice in the sum of (4.11); u-Rh(u,p)¯ and Ψ¯ are constants. Then it follows from (4.11) that
(4.13)∥u-Rh(u,p)∥0,Ω2=B((u-Rh(u,p),p-Qh(u,p));(Φ-ΦI,Ψ-ΨI))+∑KjτKj(∇p,∇ΨI)Kj-∑KjτKj(∇(p-Qh(u,p)),∇Ψ)Kj+∑Kj∑E∈∂Kj((u-Rh(u,p))·n,Ψ-Ψ¯)E-∑Kj∑E∈∂Kj(∂Φ∂n,u-Rh(u,p)-u-Rh(u,p)¯)E.
Combining (3.3) with Lemma 4.1, we deduce that
(4.14)B((u-Rh(u,p),p-Qh(u,p));(Φ-ΦI,Ψ-ΨI))≤Ch2(∥Au∥0,Ω+|p|1,Ω)(∥AΦ∥0,Ω+|Ψ|1,Ω),(4.15)∑KjτKj(∇p,∇ΨI)Kj-∑KjτKj(∇(p-Qh(u,p)),∇Ψ)Kj≤∑KjτKj|p|1,Kj|ΨI|1,Kj+∑KjτKj|p-Qh(u,p)|1,Kj|Ψ|1,Kj.
With the help of (4.12), we have
(4.16)∫E‍[(u-Rh(u,p))-u-Rh(u,p)¯]ds=0.
Combining the definition of Πj, (4.16), and local trace theorem (3.21) with the standard argument for the nonconforming element (see [21]), we see that
(4.17)∑Kj‍∑E∈∂Kj(∂Φ∂n,u-Rh(u,p)-u-Rh(u,p)¯)E=∑Kj‍∑E∈∂Kj(∂Φ∂n-∂(ΠjΦ)∂n,u-Rh(u,p)-u-Rh(u,p)¯)E≤∑Kj‍∑E∈∂Kj‍∥∇(Φ-ΠjΦ)∥L2(E)∥u-Rh(u,p)-u-Rh(u,p)¯∥L2(E)≤Ch2(∥Au∥0,Ω+|p|1,Ω)∥AΦ∥0,Ω.
In a similar way, we have
(4.18)∑Kj‍∑E∈∂Kj((u-Rh(u,p))·n,Ψ-Ψ¯)E≤Ch2(∥Au∥0,Ω+|p|1)|Ψ|1,Ω.
By combining (4.13)–(4.15) with (4.17)-(4.18), we deduce that
(4.19)∥u-Rh(u,p)∥0,Ω≤Ch2(∥Au∥0,Ω+|p|1,Ω),
which, together with (4.6). We finish the proof.Theorem 4.2.
Assume that the conditions of Theorems3.4 and 3.5 are valid; let (u,p), (uh,ph) be the solutions of (2.1) and (3.14), respectively, then
(4.20)∥u-uh∥1,h+∥p-ph∥0,Ω≤Ch.Proof.
We get the following error equation by combining (2.1) with (3.14), for all (vh,qh)∈(NCP1,P1)(4.21)B((u-uh,p-ph);(vh,qh))+b(u-uh,u,vh)+b(uh,u-uh,vh)-∑Kj〈∂u∂n,vh〉j+∑Kj〈vh·n,p〉j=∑KjτKj(∇p,∇qh)Kj-∑KjτKj(f,∇qh)Kj.
With (4.3), (4.21) can be rewritten as
(4.22)B((eh,ηh);(vh,qh))+b(u-Rh(u,p)+eh,u,vh)+b(uh,u-Rh(u,p)+eh,vh)-∑Kj〈∂u∂n,vh〉j+∑Kj〈vh·n,p〉j=-∑KjτKj(f,∇qh)Kj,
where eh=Rh(u,p)-uh and ηh=Qh(u,p)-ph.
From Theorem3.5 and (4.22), we get that
(4.23)∥ηh∥0,Ω≤β-1sup0≠(vh,qh)∈(NCP1,P1)|B((eh,ηh);(vh,qh))|∥vh∥1,h+∥qh∥0,Ω≤β-1sup0≠(vh,qh)∈(NCP1,P1)1∥vh∥1,h+∥qh∥0,Ω×(|b(u-Rh(u,p)+eh,u,vh)|+|∑KjτKj(f,∇qh)Kj|+|b(uh,u-Rh(u,p)+eh,vh)|+|∑Kj〈∂u∂n,vh〉j|+|∑Kj〈vh·n,p〉j|).
Again, with (2.13), Theorem 2.1, inverse inequality, and Lemma 4.1, we have
(4.24)|b(u-Rh(u,p),u,vh)|+|b(uh,u-Rh(u,p),vh)|≤c0γ0(|u|1,Ω+∥uh∥1,h)∥u-Rh(u,p)∥1,h∥vh∥1,h≤Ch∥vh∥1,h,(4.25)|b(eh,u,vh)|+|b(uh,eh,vh)|≤c0γ0(|u|1,Ω+∥uh∥1,h)∥eh∥1,h∥vh∥1,h,(4.26)∑KjτKj|(f,∇qh)Kj|≤∑KjτKj∥f∥0,Kj∥∇qh∥0,Kj≤∑KjC~hKj2∥f∥0,Kj∥∇qh∥0,Kj≤C(ν)h∥f∥0,Ω∥qh∥0,Ω,
and using the similar arguments as for (4.17)-(4.18) yields
(4.27)|∑Kj〈∂u∂n,vh〉j|+|∑Kj〈vh·n,p〉j|≤Ch2(∥Au∥0,Ω+|p|1,Ω)∥vh∥1,h.
Combining (4.23)–(4.27) with Theorem 3.5, we arrive at
(4.28)∥ηh∥0,Ω≤Ch+2c0γ0γ0+ν1/2ν∥f∥0,Ω∥eh∥1,h.
Choosing (vh,qh)=(eh,ηh) in (4.22), we obtain that
(4.29)B((eh,ηh);(eh,ηh))+b(eh,u,eh)+∑KjτKj(f,∇ηh)Kj=-b(u-Rh(u,p),u,eh)-b(uh,u-Rh(u,p),eh)+∑Kj〈∂u∂n,vh〉j-∑Kj〈eh·n,p〉j.
Using (2.13), (2.17), Theorem 2.1, and Lemma 3.3, we get
(4.30)B((eh,ηh);(eh,ηh))-b(eh,u,eh)=∥|eh|∥h2+∥ηh∥h2-c0γ0|u|1∥eh∥1,h2≥ν∥eh∥1,h2-c0γ0|u|1∥eh∥1,h2≥ν(1-c0γ02ν-2∥f∥0)∥eh∥1,h2>0.
Combining (4.24)–(4.28) with (4.29) yields:
(4.31)∥eh∥1,h≤Ch.
From (4.28) and (4.31), we obtain that ∥ηh∥0,Ω≤Ch. Furthermore, we finish the proof by combining triangles inequality with Lemma 4.1, (4.28), and (4.31).Theorem 4.3.
Let(u,p) and (uh,ph) be the solutions of (2.1) and (3.14), respectively, then we have
(4.32)∥u-uh∥0,Ω≤Ch2(∥Au∥0,Ω+|p|1,Ω).Proof.
Using the duality argument for a linearized stationary Navier-Stokes problem; for some giveng∈Y and the solution (u,p) of (2.1), defining (Φ,Ψ)∈(X,M) by
(4.33)-νΔΦ+∇Ψ+B~(u,Φ)-B(u,Φ)=g,inΩ,(4.34)divΦ=0inΩ,(4.35)u|∂Ω=0on∂Ω,
where B~(u,Φ) is defined as 〈v,B~(u,Φ)〉X×X′=b(v,u,Φ), for all v∈X; multiplying (4.33) and (4.34) by v∈X and q∈M, respectively; integrating over Ω, from (2.8)–(2.11), it is easily to see that the bilinear form a(Φ,v)-d(v,Ψ)+d(Φ,q) is continuity and X×M coercive, by using the Lax-Milgram’s Lemma, (4.33)–(4.35) have a unique solution (Φ,Ψ).
Multiplying (4.33) and (4.34) by Φ and Ψ, respectively, using (2.13) and Theorem 2.1, we have
(4.36)ν|Φ|1,Ω2-b(Φ,u,Φ)≥ν|Φ|1,Ω2-c0γ0|u|1,Ω|Φ|1,Ω2≥ν(1-c0γ02ν-2∥f∥0,Ω)|Φ|1,Ω2>0.
On the other hand, estimating the right term yields
(4.37)(Φ,g)≤∥Φ∥0,Ω∥g∥0,Ω≤γ0|Φ|1,Ω∥g∥0,Ω.
By using (4.36) and (4.37), we arrive at
(4.38)|Φ|1,Ω≤C∥g∥0,Ω.
Setting Ψ=0 and taking the scalar product of (4.33) with AΦ in Y yields
(4.39)ν∥AΦ∥02+b(AΦ,u,Φ)-b(u,Φ,AΦ)=(g,AΦ).
Using the Gagliardo-Nirenberg inequality yields
(4.40)∥v∥L42≤C∥v∥0,Ω|v|1,Ω,∀v∈X,∥∇v∥L42≤C∥Av∥0,Ω|v|1,Ω,∀v∈H2(Ω)2∩X.
With the help of the Agmon’s inequality, we have
(4.41)∥v∥L∞2≤C∥v∥0,Ω∥Av∥0,Ω,∀v∈H2(Ω)2∩X.
Furthermore, the following estimates are hold:
(4.42)b(AΦ,u,Φ)≤C∥AΦ∥0,Ω(∥∇u∥L4∥Φ∥L4+∥u∥L∞∥∇Φ∥0,Ω)≤ν4∥AΦ∥0,Ω2+Cν-1∥Au∥0,Ω2∥∇Φ∥0,Ω2),b(u,Φ,AΦ)≤ν4∥AΦ∥0,Ω2+Cν-1∥Au∥0,Ω2∥∇Φ∥0,Ω2),(g,AΦ)≤ν4∥AΦ∥0,Ω2+Cν-1∥g∥0,Ω.
Combining above inequalities with (4.38) and (4.39), we arrive at
(4.43)∥AΦ∥0,Ω2≤cν-2(1+∥Au∥0,Ω2)∥g∥0,Ω2.
Applying the continuous Inf-Sup condition (2.8) yields
(4.44)|Ψ|1,Ω≤cν∥AΦ∥0,Ω+c∥Au∥0,Ω|Φ|1,Ω+c∥g∥0,Ω.
Combining (4.43)-(4.44) with (2.18), (4.38), we arrive at
(4.45)∥AΦ∥0,Ω+|Ψ|1,Ω≤C∥g∥0,Ω.
Taking g=Rh(u,p)-uh, multiplying (4.33) and (4.34) by Rh(u,p)-uh and Qh(u,p)-ph, respectively, using (4.1) yields
(4.46)∥Rh(u,p)-uh∥0,Ω2=B((Rh(u,p)-uh,Qh(u,p)-ph);(Φ,Ψ))+b(u,Rh(u,p)-uh,Φ)-∑KjτKj(∇(Qh(u,p)-ph),∇Ψ)Kj+∑Kj〈Ψ,(Rh(u,p)-uh)·n,〉j+b(Rh(u,p)-uh,u,Φ)-∑Kj〈∂Φ∂n,Rh(u,p)-uh〉j.
Setting (vh,qh)=(ΦI,ΨI) in (4.22), and using (4.46), we obtain that
(4.47)∥Rh(u,p)-uh∥0,Ω2=B((Rh(u,p)-uh,Qh(u,p)-ph);(Φ-ΦI,Ψ-ΨI))-∑KjτKj(∇(Qh(u,p)-ph),∇Ψ)Kj+b(u-uh,Rh(u,p)-uh,Φ)+∑Kj〈∂u∂n,vh〉j+b(u-Rh(u,p),u,ΦI)+b(Rh(u,p)-uh,u,Φ-ΦI)+b(uh,u-Rh(u,p),Φ-ΦI)+b(uh,Rh(u,p)-uh,Φ-ΦI)-b(uh,u-Rh(u,p),Φ)+∑Kj〈Ψ,(Rh(u,p)-uh)·n,〉j-∑KjτKj(f,∇ΨI)Kj-∑Kj〈∂Φ∂n,Rh(u,p)-uh〉j-∑Kj〈vh·n,p〉j.
We now estimate the right terms of (4.47). Thanks to (3.3); Theorems 2.1, 3.4-3.5, and 4.2; inverse inequality; Lemma 4.1; we know that
(4.48)B((Rh(u,p)-uh,Qh(u,p)-ph);(Φ-ΦI,Ψ-ΨI))≤C(∥Rh(u,p)-uh∥1,h+∥Qh(u,p)-ph∥0,Ω)(∥Φ-ΦI∥1,h+∥Ψ-ΨI∥0,Ω+h|Ψ|1,Ω)≤Ch2(∥Au∥0,Ω+|p|1,Ω)(∥AΦ∥0,Ω+|Ψ|1,Ω),b(Rh(u,p)-uh,u,Φ-ΦI)+b(uh,Rh(u,p)-uh,Φ-ΦI)≤2c0(|u|1,Ω+∥uh∥1,h)∥Rh(u,p)-uh∥1,h∥Φ-ΦI∥1,h≤Ch2(∥Au∥0,Ω+|p|1,Ω)∥AΦ∥0,Ω,∑KjτKj(∇(Qh(u,p)-ph),∇Ψ)Kj≤∑Kj‍C~hK2∥∇(Qh(u,p)-ph)∥0,Kj∥∇Ψ∥0,Kj≤Ch2(∥Au∥0,Ω+|p|1)|Ψ|1,Ω,b(u-Rh(u,p),u,ΦI)≤C∥u-Rh(u,p)∥0,Ω∥Au∥0,Ω|ΦI|1,h≤Ch2(∥Au∥0,Ω+|p|1,Ω)|Φ|1,Ω,b(u-uh,Rh(u,p)-uh,Φ)≤c0∥u-uh∥1,h∥Rh(u,p)-uh∥1,h|Φ|1,Ω≤Ch2(∥Au∥0,Ω+|p|1)|Φ|1,Ω,b(uh,u-Rh(u,p),Φ)≤C∥uh∥1,h∥u-Rh(u,p)∥0,Ω∥AΦ∥0,Ω,∑KjτKj(f,∇ΨI)Kj≤∑KjC~hKj2∥f∥0,Kj∥∇ΨI∥0,Kj≤Ch2∥f∥0,Ω|Ψ|1,Ω,b(uh,u-Rh(u,p),Φ-ΦI)≤c0∥uh∥1,h∥u-Rh(u,p)∥1,h∥Φ-ΦI∥1,h≤Ch2(∥Au∥0,Ω+|p|1,Ω)∥AΦ∥0,Ω.
Applying the argument used for (4.18)-(4.21), (4.28) gives
(4.49)∑Kj〈∂Φ∂n,Rh(u,p)-uh〉j+∑Kj〈Ψ,(Rh(u,p)-uh)·n,〉j≤Ch2(∥Au∥0,Ω+|p|1,Ω)(∥AΦ∥0,Ω+|Ψ|1,Ω),∑Kj〈∂u∂n,ΦI〉+∑Kj‍〈ΦI·n,p〉j≤Ch2(∥Au∥0,Ω+∥p∥1,Ω)∥ΦI∥1,h≤Ch2(∥Au∥0,Ω+∥p∥1,Ω)∥ΦI∥1,Ω.
Combining above inequalities, Theorems 2.1, 3.5, (4.46), and (4.51), we get that
(4.50)∥Rh(u,p)-uh∥0,Ω2≤Ch2(∥Au∥0,Ω+|p|1,Ω)(∥AΦ∥0,Ω+|Ψ|1,Ω)+C∥uh∥1,h∥u-Rh(u,p)∥0,Ω∥AΦ∥0,Ω≤Ch2(∥Au∥0,Ω+|p|1,Ω)∥Rh(u,p)-uh∥1,Ω+Cγ0+ν1/2ν∥f∥0,Ω∥u-Rh(u,p)∥0,Ω2.
Choosing the appropriate ν,Ω and f such that 1-C((γ0+ν1/2)/ν)∥f∥0>0, then we have
(4.51)∥Rh(u,p)-uh∥0,Ω≤Ch2(∥Au∥0,Ω+∥p∥1,Ω).
By applying the triangles inequality, Lemma 4.1 and (4.51), we finish the proof.Lemma 4.4.
Under the assumptions of Theorem4.2, the following estimate about u-uh in mesh-dependent norm
(4.52)∥|u-uh|∥h≤Ch,
holds, where u and uh are the solution of problem (2.1) and (3.14), respectively.Proof.
According to the definition of∥|·|∥h, uh∈Xh, with Theorems 2.1 and 4.2, inverse and local trace inequalities (3.21), we have
(4.53)∥|u-uh|∥h2=ν∥|u-uh|∥1,h2+∑E∈ΓhτE∥ν∂n(u-uh)∥0,E2=ν∥|u-uh|∥1,h2+∑E∈Γh‍he12ν∥ν∂n(u-uh)∥0,E2≤Ch2+C∑E∈Γh‍he[hKj-1∥∇(u-uh)·n∥0,Kj2+hKj|∇(u-uh)·n|1,Kj2]≤Ch2.


## 5.  Numerical Validations

In this section, we provide two numerical examples to illustrate the theoretical analysis of the method (3.14). In all experiments, we consider the domain Ω to be the square [0,1]×[0,1]. The mesh consists of triangular elements that are obtained by dividing Ω into subsquares of equal size and then drawing the diagonal in each sub-square, see Figure 1. The software Freefem++, developed by Hecht et al. [23], is used in our experiments.Figure 1
Uniform triangulation ofΩ into triangulars.


### 5.1. An Analytical Solution: Convergence Validation

For this test, our purpose is to verify the theoretical analysis which has been established in the previous section by setting the viscosity coefficientν=1 and f is given by the exact solution
(5.1)u1=10x2(x-1)2y(y-1)(2y-1),u2=-10x(x-1)(2x-1)y2(y-1)2,p=10(2x-1)(2y-1).In our numerical validation, The experimental rates of convergence with respect to the mesh sizeh are calculated by the formula (log(Ei/Ei+1))/(log(hi/hi+1)), where Ei and Ei+1 are the relative errors corresponding to the meshes of sizes hi and hi+1.In order to show the efficiency of the enriched multiscale method, we compare the numerical results obtained by using different methods, which are shown in Tables1, 2, 3, and 4. The compared methods include the pressure projection method (PPM) in [6], the local Gauss integration method (LGIM) in [7, 10, 19], and the standard Galerkin method (SGM) with MINI element (see [1]), respectively. From these tables, we can see that the stabilized multiscale method has good precision for pressure, and the precision of velocity, worse than other methods. Table 5 explains the CPU times that needed for solving the steady Navier-Stokes equations in different mesh sizes. From these data, we know that our method takes less time than other methods. Furthermore, from Tables 1–4, we can see that the numerical results reproduce the established theoretical analysis and show an 𝒪(h) order of convergence for |u-uh|1,Ω and ∥p-ph∥0,Ω, and an 𝒪(h2) convergence for ∥u-uh∥0,Ω.Table 1
Numerical results for enriched multiscale method withNCP1-P1 element.
1
/
h
∥
p
h
-
p
~
∥
0
/
∥
p
~
∥
0
∥
u
h
-
u
~
∥
0
/
∥
u
~
∥
0
|
u
h
-
u
~
|
1
/
|
u
~
|
1
p
L
2 rate
u
L
2 rate
u
H
1 rate
20
0.160667
0.15627
0.816235
25
0.126534
0.10132
0.665410
1.0703
1.9418
0.9155
30
0.105785
0.0708881
0.556179
0.9823
1.9591
0.9835
35
0.083427
0.0523941
0.478261
1.5403
1.9611
0.9791
40
0.0667358
0.0405399
0.420987
1.6717
1.9209
0.9552Table 2
Numerical results for PPM in [6] with NCP1-P1 element.
1
/
h
∥
p
h
-
p
∥
0
/
∥
p
~
∥
0
∥
u
h
-
u
∥
0
/
∥
u
∥
0
∥
u
h
-
u
∥
1
/
∥
u
∥
1
p
L
2 rate
u
L
2 rate
u
H
1 rate
20
0.208238
0.0182227
0.12718
1.8797
1.9629
0.9963
25
0.136834
0.0117122
0.101772
1.8818
1.9809
0.9988
30
0.0972068
0.0081486
0.0848213
1.8754
1.9898
0.9993
35
0.0727087
0.00599297
0.0727009
1.8837
1.9932
1.0002
40
0.056478
0.00459177
0.0636082
1.8918
1.9945
1.0006Table 3
Numerical results for SGM withP1b-P1 element.
1
/
h
∥
p
h
-
p
∥
0
/
∥
p
∥
0
∥
u
h
-
u
∥
0
/
∥
u
∥
0
∥
u
h
-
u
∥
1
/
∥
u
∥
1
p
L
2 rate
u
L
2 rate
u
H
1 rate
20
0.194855
0.0167771
0.120396
1.9352
2.0150
1.0040
25
0.127314
0.010703
0.0962298
1.9073
2.0144
1.0041
30
0.0903239
0.00741503
0.080144
1.8827
2.0130
1.0033
35
0.0675341
0.00543371
0.068793
1.8863
2.0168
0.9907
40
0.0527341
0.00416312
0.059691
1.8532
1.9947
1.0628Table 4
Numerical results for LGIM in [4, 7, 19] with NCP1-P1 element.
1
/
h
∥
p
h
-
p
∥
0
/
∥
p
∥
0
∥
u
h
-
u
∥
0
/
∥
u
∥
0
∥
u
h
-
u
∥
1
/
∥
u
∥
1
p
L
2 rate
u
L
2 rate
u
H
1 rate
20
0.263248
0.0186335
0.127546
1.7919
2.0088
1.0010
25
0.176286
0.0118944
0.102010
1.7970
2.0117
1.0012
30
0.128021
0.00824161
0.0849864
1.7547
2.0122
1.0014
35
0.0978302
0.00604517
0.0728229
1.7448
2.0106
1.0020
40
0.0776455
0.00462309
0.0637016
1.7305
2.0085
1.0022Table 5
CPU time for solving the steady Navier-Stokes equations with different mesh sizes.
1
/
h
20
25
30
35
40
Our method withNCP1-P1
1.6324
3.286
5.975
8.792
14.088
PPM in [6] with NCP1-P1
1.875
3.594
6.406
10.891
17.453
LGIM in [10] with NCP1-P1
1.844
3.547
6.391
10.593
17.047
SGM withP1b-P1
3.031
6.078
11.375
21.106
35.020


### 5.2. Lid-Driven Cavity Problem

In this test, we consider the incompressible lid-driven cavity flow problem defined on the unit square. Settingf=0 and the boundary condition u=0 on [{0}×(0,1)]∪[(0,1)×{0}]∪[{1}×(0,1)] and u=(1,0)T on (0,1)×{1}, see Figure 2. The mesh consists of triangular element and the mesh size h=1/50.Figure 2
Lid-driven cavity flow.Figure3 shows the pressure contours at different Reylond numbers, where the stopping criterion ∥uhn+1-uhn∥0,Ω/∥uhn+1∥0,Ω≤10-6 is employed, where uhn+1 is the approximation of uh at the n+1 Newton iterative. From Figure 3, we can see that the oscillations are absented for the pressure isovalues by the NCP1-P1 approximations, and compared with the results given in [24], we can see that our method has the effect to stabilized flow field (see Figure 4). In this sense, we say that the stabilized multiscale nonconforming finite element method is effective for the stationary Navier-Stokes problem.Pressure contours of the driven cavity by using the stabilized multiscale method withNCP1-P1 at different Reynolds numbers. (a) Re = 1, (b) Re = 100.
(a)
(b)The computed velocity profiles through the geometric center at Re = 100. (a) Horizontal velocity, (b) vertical velocity.
(a)
(b)


## 5.1. An Analytical Solution: Convergence Validation

For this test, our purpose is to verify the theoretical analysis which has been established in the previous section by setting the viscosity coefficientν=1 and f is given by the exact solution
(5.1)u1=10x2(x-1)2y(y-1)(2y-1),u2=-10x(x-1)(2x-1)y2(y-1)2,p=10(2x-1)(2y-1).In our numerical validation, The experimental rates of convergence with respect to the mesh sizeh are calculated by the formula (log(Ei/Ei+1))/(log(hi/hi+1)), where Ei and Ei+1 are the relative errors corresponding to the meshes of sizes hi and hi+1.In order to show the efficiency of the enriched multiscale method, we compare the numerical results obtained by using different methods, which are shown in Tables1, 2, 3, and 4. The compared methods include the pressure projection method (PPM) in [6], the local Gauss integration method (LGIM) in [7, 10, 19], and the standard Galerkin method (SGM) with MINI element (see [1]), respectively. From these tables, we can see that the stabilized multiscale method has good precision for pressure, and the precision of velocity, worse than other methods. Table 5 explains the CPU times that needed for solving the steady Navier-Stokes equations in different mesh sizes. From these data, we know that our method takes less time than other methods. Furthermore, from Tables 1–4, we can see that the numerical results reproduce the established theoretical analysis and show an 𝒪(h) order of convergence for |u-uh|1,Ω and ∥p-ph∥0,Ω, and an 𝒪(h2) convergence for ∥u-uh∥0,Ω.Table 1
Numerical results for enriched multiscale method withNCP1-P1 element.
1
/
h
∥
p
h
-
p
~
∥
0
/
∥
p
~
∥
0
∥
u
h
-
u
~
∥
0
/
∥
u
~
∥
0
|
u
h
-
u
~
|
1
/
|
u
~
|
1
p
L
2 rate
u
L
2 rate
u
H
1 rate
20
0.160667
0.15627
0.816235
25
0.126534
0.10132
0.665410
1.0703
1.9418
0.9155
30
0.105785
0.0708881
0.556179
0.9823
1.9591
0.9835
35
0.083427
0.0523941
0.478261
1.5403
1.9611
0.9791
40
0.0667358
0.0405399
0.420987
1.6717
1.9209
0.9552Table 2
Numerical results for PPM in [6] with NCP1-P1 element.
1
/
h
∥
p
h
-
p
∥
0
/
∥
p
~
∥
0
∥
u
h
-
u
∥
0
/
∥
u
∥
0
∥
u
h
-
u
∥
1
/
∥
u
∥
1
p
L
2 rate
u
L
2 rate
u
H
1 rate
20
0.208238
0.0182227
0.12718
1.8797
1.9629
0.9963
25
0.136834
0.0117122
0.101772
1.8818
1.9809
0.9988
30
0.0972068
0.0081486
0.0848213
1.8754
1.9898
0.9993
35
0.0727087
0.00599297
0.0727009
1.8837
1.9932
1.0002
40
0.056478
0.00459177
0.0636082
1.8918
1.9945
1.0006Table 3
Numerical results for SGM withP1b-P1 element.
1
/
h
∥
p
h
-
p
∥
0
/
∥
p
∥
0
∥
u
h
-
u
∥
0
/
∥
u
∥
0
∥
u
h
-
u
∥
1
/
∥
u
∥
1
p
L
2 rate
u
L
2 rate
u
H
1 rate
20
0.194855
0.0167771
0.120396
1.9352
2.0150
1.0040
25
0.127314
0.010703
0.0962298
1.9073
2.0144
1.0041
30
0.0903239
0.00741503
0.080144
1.8827
2.0130
1.0033
35
0.0675341
0.00543371
0.068793
1.8863
2.0168
0.9907
40
0.0527341
0.00416312
0.059691
1.8532
1.9947
1.0628Table 4
Numerical results for LGIM in [4, 7, 19] with NCP1-P1 element.
1
/
h
∥
p
h
-
p
∥
0
/
∥
p
∥
0
∥
u
h
-
u
∥
0
/
∥
u
∥
0
∥
u
h
-
u
∥
1
/
∥
u
∥
1
p
L
2 rate
u
L
2 rate
u
H
1 rate
20
0.263248
0.0186335
0.127546
1.7919
2.0088
1.0010
25
0.176286
0.0118944
0.102010
1.7970
2.0117
1.0012
30
0.128021
0.00824161
0.0849864
1.7547
2.0122
1.0014
35
0.0978302
0.00604517
0.0728229
1.7448
2.0106
1.0020
40
0.0776455
0.00462309
0.0637016
1.7305
2.0085
1.0022Table 5
CPU time for solving the steady Navier-Stokes equations with different mesh sizes.
1
/
h
20
25
30
35
40
Our method withNCP1-P1
1.6324
3.286
5.975
8.792
14.088
PPM in [6] with NCP1-P1
1.875
3.594
6.406
10.891
17.453
LGIM in [10] with NCP1-P1
1.844
3.547
6.391
10.593
17.047
SGM withP1b-P1
3.031
6.078
11.375
21.106
35.020


## 5.2. Lid-Driven Cavity Problem

In this test, we consider the incompressible lid-driven cavity flow problem defined on the unit square. Settingf=0 and the boundary condition u=0 on [{0}×(0,1)]∪[(0,1)×{0}]∪[{1}×(0,1)] and u=(1,0)T on (0,1)×{1}, see Figure 2. The mesh consists of triangular element and the mesh size h=1/50.Figure 2
Lid-driven cavity flow.Figure3 shows the pressure contours at different Reylond numbers, where the stopping criterion ∥uhn+1-uhn∥0,Ω/∥uhn+1∥0,Ω≤10-6 is employed, where uhn+1 is the approximation of uh at the n+1 Newton iterative. From Figure 3, we can see that the oscillations are absented for the pressure isovalues by the NCP1-P1 approximations, and compared with the results given in [24], we can see that our method has the effect to stabilized flow field (see Figure 4). In this sense, we say that the stabilized multiscale nonconforming finite element method is effective for the stationary Navier-Stokes problem.Pressure contours of the driven cavity by using the stabilized multiscale method withNCP1-P1 at different Reynolds numbers. (a) Re = 1, (b) Re = 100.
(a)
(b)The computed velocity profiles through the geometric center at Re = 100. (a) Horizontal velocity, (b) vertical velocity.
(a)
(b)


## 6. Conclusion

In this paper we have derived a theoretical analysis of enriched multiscale nonconforming finite element method for the steady incompressible Navier-Stokes equations. The analysis has extended the work in [16] from the linear problem to the nonlinear problem. The discretization uses nonconforming and conforming piecewise linear finite elements for velocity and pressure over triangles elements, respectively. Numerical tests show that this stabilized method is computationally efficient, and it can be performed locally at the element level with minimal additional cost; at the same time, our numerical results obtained are in good agreement with the established theoretical results.



---

*Source: 651808-2012-10-02.xml*