ABSTRACT

The Weak Galerkin (WG) finite element method for the unsteady Stokes equations in the primary velocity-pressure formulation is introduced in this paper. Optimal-order error estimates are established for the corresponding numerical approximation in an H¹ norm for the velocity, and L² norm for both the velocity and the pressure by use of the Stokes projection.

Keywords:

Weak Galerkin Finite Element Methods, Unsteady Stokes Equations, Stokes Projection

1. Introduction

The finite element method for the unsteady Stokes equations developed over the last several decades is based on the weak formulation by constructing a pair of finite element spaces satisfying the inf-sup condition of Babuska [1] and Brezzi [2] . Readers are referred to [3] [4] [5] [6] [7] for specific examples and details in the different finite element methods for the Stokes equations. The idea of weak Galerkin method was first introduced by the Professor Junping Wang in June 2011. Weak Galerkin refers to a general finite element technique for partial differential equations in which differential operators are approximated by weak forms as distributions for generalized functions. Thus, two of the key features in weak Galerkin methods are 1) the approximating functions are discontinuous, and 2) the usual derivatives are taken as distributions or approximations of distributions. The method was successfully applied to the second order elliptic equations [8] [9] , the Stokes equations [10] , Parabolic equations [11] , and Maxwell equations [12] . A posteriori error is effectively estimated, and proved the convergence of the WG finite element method in this paper.

2. Preliminaries

In this paper, we study the initial-boundary value problems of the Stokes.

$\{\begin{array}{l}{u}_t-\Delta u+\nabla p=f\left(x\mathrm{,}t\right)\mathrm{,\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>}t\in \left(\mathrm{0,}T\right)\hfill \\ \nabla \cdot u=\mathrm{0,\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>}\left(x\mathrm{,}t\right)\in \Omega \times \left(\mathrm{0,}T\right)\hfill \\ u=\mathrm{0,\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>}\left(x\mathrm{,}t\right)\in \partial \Omega \times \left(\mathrm{0,}T\right)\hfill \\ u\left(x\mathrm{,0}\right)=\psi \left(x\right)\mathrm{<mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>}x\in \Omega \mathrm{,}t=0\hfill \end{array}$ (1)

where $u={\left(u_1\mathrm{,}{u}_2\right)}^{\text{T}}$ is fluid velocity, p is pressure, $f={\left(f_1\mathrm{,}{f}_2\right)}^{\text{T}}$ is volumetric power density.

The solution of the Stokes equations forms an important aspect of both theoretical and computational fluid dynamics. A limited number of solutions of these non-linear partial differential equations mostly involving spatially one-dimensional problems are given in the literature. Solutions of practical interest have been obtained for cases where, with suitable approximations, the equations are reduced to linear partial differential equations.

Let W be a bounded domain in R². We introduce function spaces

$X={\left[H_0^1\left(\Omega \right)\right]}^2$ , $V=\left\{u\in X,div\mathrm{<mi>\; </mi>}u=0\right\}$ , $M=\left\{q\in L_2\left(\Omega \right);{\displaystyle {\int }_{\Omega }}q\text{d}x\text{d}y=0\right\}$ , then the unsteady Stokes problem would take the following form: seek $\left(u\mathrm{,}p\right)\in X\times M$ satisfying

$\{\begin{array}{l}\left(u_t,v\right)+\left(\nabla u,\nabla v\right)-\left(p,\nabla \cdot v\right)=\left(f,v\right),\mathrm{<mi>\; </mi>\; <mi>\; </mi>}\forall v\in X\hfill \\ \left(q,\nabla \cdot u\right)=0,\mathrm{<mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>}\forall q\in W\hfill \\ u\left(x,0\right)=\psi \left(x\right)\hfill \end{array}$ (2)

We use ${\Vert \cdot \Vert }_{s\mathrm{,}D}$ and ${\left|\cdot \right|}_{s\mathrm{,}D}$ to be denote the norm and Semi-norm in the Sobolev space $H^s\left(D\right)$ for any $s\ge 0$ , respectively. The inner product in $H^s\left(D\right)$ is denoted by ${\left(\cdot \mathrm{,}\cdot \right)}_{s\mathrm{,}D}$ . For example, for each $s\ge 0$ , the Semi-norm ${\left|\cdot \right|}_{s\mathrm{,}D}$ is given by

${\left|\phi \right|}_{s,D}={\left({\displaystyle \underset{\left|\gamma \right|=s}{\sum }}{\displaystyle {\int }_D}{\left|{\partial }^{\gamma }\phi \right|}^2\text{d}D\right)}^{\frac{1}{2}}$

and $\Vert \cdot \Vert$ is said to be the norm of $L_2$ .

For w is $\left[\mathrm{0,}T\right]$ to $H^s\left(D\right)$ , the definition is given by

${\Vert w\Vert }_{L^q\left(0,T;H^s\left(D\right)\right)}={\left({\displaystyle {\int }_0^T}{\Vert w\left(\cdot ,t\right)\Vert }_{s,D}^q\text{d}t\right)}^{\frac{1}{q}}$

for $1\le q\le \infty$ , we have

${\Vert w\Vert }_{L^{\infty }\left(0,T;H^s\left(D\right)\right)}=\underset{0\le t\le T}{sup}{\Vert w\left(\cdot \mathrm{,}t\right)\Vert }_{s\mathrm{,}D}$

The space $H^s\left(D\right)$ and the norm defined in the $H^s\left(D\right)$ defined as

$H\left(div\mathrm{,}\Omega \right)=\left\{q\mathrm{<mo>\; :\; </mo>}q\in \left[L^2\left(\Omega \right)\right]\mathrm{,}\nabla \cdot q\in L^2\left(\Omega \right)\right\}$

${\Vert q\Vert }_{H\left(div\mathrm{,}\Omega \right)}={\left({\Vert q\Vert }^2+{\Vert \nabla \cdot q\Vert }^2\right)}^{\frac{1}{2}}$

3. Weak Galerkin Finite Element Approximation Scheme

Let K be any polygonal or polyhedral domain with boundary $\partial K$ . A weak vector-valued function on the region K refers to a vector-valued function $v=\left\{v_0,v_b\right\}$ such that $v_0\in {\left[L^2\left(K\right)\right]}^d$ and $v_b\in {\left[H^{\frac{1}{2}}\left(\partial K\right)\right]}^d$ . The first component $v_0$ can be understood as the value of v in K, and the second component $v_b$ represents v on the boundary of K. Note that $v_b$ may not necessarily be related to the trace of $v_0$ on $\partial K$ should a trace be well-defined. Denote by $\nu \left(K\right)$ the space of weak functions on K;

$\nu \left(K\right)=\left\{v=\left\{v_0,v_b\right\}:v_0\in {\left[L^2\left(K\right)\right]}^d,v_b\in {\left[H^{\frac{1}{2}}\left(\partial K\right)\right]}^d\right\}$ (3)

Definition 1. For any $v\in \nu \left(K\right)$ , the weak gradient of v is defined as a linear functional ${\nabla }_{w}v$ in the dual space of ${\left[H\left(div\mathrm{,}K\right)\right]}^d$ , whose action on each $q\in {\left[H\left(div\mathrm{,}K\right)\right]}^d$ is given by

${\left({\nabla }_{w}v\mathrm{,}q\right)}_K=-{\left(v_0\mathrm{,}\nabla \cdot q\right)}_K+{\langle v_b\mathrm{,}q\cdot n\rangle }_{\partial K}$ (4)

where n is the outward normal direction to $\partial K$ , ${\left(v_0\mathrm{,}\nabla \cdot q\right)}_K={\displaystyle {\int }_K}{v}_0\left(\nabla \cdot q\right)\text{d}K$ is the action of $v_0$ on $\nabla \cdot q$ , and ${\langle v_b\mathrm{,}q\cdot n\rangle }_{\partial K}={\displaystyle {\int }_{\partial K}}{v}_{b}q\cdot n\text{d}s$ is the action of $q\cdot n$ on $v_b\in {\left[H^{\frac{1}{2}}\left(\partial K\right)\right]}^d$ .

The Sobolev space ${\left[H^1\left(K\right)\right]}^d$ can be embedded into the space $\nu \left(K\right)$ by an inclusion map $i\nu \mathrm{<mo>\; :\; </mo>}{\left[H^1\left(K\right)\right]}^d\to \nu \left(K\right)$ defined as follows

$i\nu \left(\varphi \right)=\left\{{\varphi |}_K\mathrm{,}{\varphi |}_{\partial K}\right\}\mathrm{,\; <mi>\; </mi>}\varphi \in {\left[H^1\left(K\right)\right]}^d$

With the help of the inclusion map $i\nu$ , the Sobolev space ${\left[H^1\left(K\right)\right]}^d$ can be viewed as a subspace of $\nu \left(K\right)$ by identifying each $\varphi \in {\left[H^1\left(K\right)\right]}^d$ with $i\nu \left(\varphi \right)$ .

Let $P_r\left(K\right)$ be the set of polynomials on K with degree no more than r.

Definition 2. The discrete weak gradient operator, denoted by ${\nabla }_{w\mathrm{,}r\mathrm{,}K}$ , is defined as the unique polynomial $\left({\nabla }_{w\mathrm{,}r\mathrm{,}K}v\right)\in {\left[P_r\left(K\right)\right]}^{d\times d}$ satisfying the following equation,

${\left({\nabla }_{w\mathrm{,}r\mathrm{,}K}v\mathrm{,}q\right)}_K=-{\left(v_0\mathrm{,}\nabla \cdot q\right)}_K+{\langle v_b\mathrm{,}q\cdot n\rangle }_{\partial K}\mathrm{,}$ (5)

for all $q\in {\left[P_r\left(K\right)\right]}^{d\times d}$ .

In what follows, we give the definition of weak divergence, first of all, we require weak function $v=\left\{v_0,v_b\right\}$ such that $v_0\in {\left[L^2\left(K\right)\right]}^d$ an $v_b\cdot n\in L^2\left(\partial K\right)$ Denote by $V\left(K\right)$ the space of weak vector-valued functions on K;

$V\left(K\right)=\left\{v=\left\{v_0,v_b\right\}:v_0\in {\left[L^2\left(K\right)\right]}^d,v_b\cdot n\in L^2\left(\partial K\right)\right\}$ (6)

Definition 3. For any $v\in V\left(K\right)$ , the weak divergence of v is defined as a linear functional ${\nabla }_w\cdot v$ in the dual space of $H^1\left(K\right)$ whose action on each $\phi \in H^1\left(K\right)$ is given by

${\left({\nabla }_w\cdot v\mathrm{,}\phi \right)}_K=-{\left(v_0\mathrm{,}\nabla \phi \right)}_K+{\langle v_b\cdot n\mathrm{,}\phi \rangle }_{\partial K}$ (7)

where n is the outward normal direction to $\partial K$ , ${\left(v_0\mathrm{,}\nabla \phi \right)}_K$ is the action of $v_0$

on $\nabla \phi$ , and ${\langle v_b\cdot n\mathrm{,}\phi \rangle }_{\partial K}$ is the action of $v_b\cdot n$ on $\phi \in H^{\frac{1}{2}}\left(\partial K\right)$ .

The Sobolev space ${\left[H^1\left(K\right)\right]}^d$ can be embedded into the space $V\left(K\right)$ by an inclusion map $iv\mathrm{<mo>\; :\; </mo>}{\left[H^1\left(K\right)\right]}^d\to V\left(K\right)$ defined as follows

$i\nu \left(\varphi \right)=\left\{{\varphi |}_K\mathrm{,}{\varphi |}_{\partial K}\right\}\mathrm{,\; <mi>\; </mi>}\varphi \in {\left[H^1\left(K\right)\right]}^d$

Definition 4. A discrete weak divergence operator, denoted by ${\nabla }_{w\mathrm{,}r\mathrm{,}K}$ , is defined as the unique polynomial $\left({\nabla }_{w\mathrm{,}r\mathrm{,}K}\cdot v\right)\in P_r\left(K\right)$ that satisfies the following equation.

${\left({\nabla }_{w\mathrm{,}r\mathrm{,}K}\cdot v\mathrm{,}\phi \right)}_K=-{\left(v_0\mathrm{,}\nabla \phi \right)}_K+{\langle v_b\cdot n\mathrm{,}\phi \rangle }_{\partial K}\mathrm{,}$ (8)

for all $\phi \in P_r\left(K\right)$ .

4. Weak Galerkin Finite Element Scheme

Let ${\mathcal{T}}_h$ be a partition of the domain W with mesh size h that consists of arbitrary polygons/polyhedra. In this paper, we assume that the partition ${\mathcal{T}}_h$ is WG shape regular-defined by a set of conditions as detailed in references. Denote by ${\epsilon }_h$ the set of all edges/flat faces in ${\mathcal{T}}_h$ , and let ${\epsilon }_h^0={\epsilon }_h\partial \Omega$ be the set of all interior edges/faces. For any integer $k\ge 1$ , we define a weak Galerkin finite element space for the velocity variable as follows,

$\text{}\text{}\text{}\text{}\text{}\text{}\text{}\text{}{V}_h=\left\{v=\left\{v_0,v_b\right\}:{\left\{v_0,v_b\right\}|}_T\in {\left[P_k\left(T\right)\right]}^d\times {\left[P_{k-1}\left(e\right)\right]}^d,e\subset \partial T\right\}$

We would like to emphasize that there is only a single value $v_b$ defined on each edge $e\in {\epsilon }_h$ . For the pressure variable, we have the following finite element space

$W_h=\left\{q:q\in L_0^2\left(\Omega \right),{q|}_T\in P_{k-1}\left(T\right)\right\}$

Denote by $V_h^0$ the subspace of $V_h$ consisting of discrete weak functions with vanishing boundary value;

$V_h^0=\left\{v=\left\{v_0,v_b\right\}\in V_h,v_b=0\text{on}\partial \Omega \right\}$

The discrete weak gradient ${\nabla }_{w\mathrm{,}k-1}$ and the discrete weak divergence $\left({\nabla }_{w\mathrm{,}k-1}\right)$ on the finite element space $V_h$ can be computed by using (5) and (8) on each element T, respectively. More precisely, they are given by

${\left({\nabla }_{w,k-1}v\right)|}_T={\nabla }_{w,k-1,T}\left({v|}_T\right),\forall v\in V_h$

${\left({\nabla }_{w\mathrm{,}k-1}\cdot v\right)|}_T={\nabla }_{w\mathrm{,}k-\mathrm{1,}T}\cdot \left({v|}_T\right)\mathrm{,}\forall v\in V_h$

For simplicity of notation, from now on we shall drop the subscript $k-1$ in the notation ${\nabla }_{w\mathrm{,}k-1}$ and $\left({\nabla }_{w\mathrm{,}k-1}\right)$ for the discrete weak gradient and the discrete weak divergence. The usual $L^2$ inner product can be written locally on each element as follows

$\left({\nabla }_{w}v\mathrm{,}{\nabla }_{w}w\right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\left({\nabla }_{w}v\mathrm{,}{\nabla }_{w}w\right)}_T$

$\left({\nabla }_w\cdot v\mathrm{,}q\right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\left({\nabla }_w\cdot v\mathrm{,}q\right)}_T$

Denote by $Q_0$ the L² projection operator from ${\left[L^2\left(T\right)\right]}^d$ onto ${\left[P_k\left(T\right)\right]}^d$ . For each edge/face $e\in {\epsilon }_h$ , denote by $Q_b$ the L² projection from ${\left[L^2\left(e\right)\right]}^d$ onto ${\left[P_{k-1}\left(e\right)\right]}^d$ . We shall combine $Q_0$ with $Q_b$ by writing $Q_h=\left\{Q_0{Q}_b\right\}$ .

We are now in a position to describe a weak Galerkin finite element scheme for the Stokes Equations (1). To this end, we first introduce three bilinear forms as follows

$s\left(v\mathrm{,}w\right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{h}_T^{-1}{\langle Q_b{v}_0-v_b\mathrm{,}{Q}_b{w}_0-w_b\rangle }_{{\partial }_T}$

$a\left(v\mathrm{,}w\right)=\left({\nabla }_{w}v\mathrm{,}{\nabla }_{w}w\right)+s\left(v\mathrm{,}w\right)$

$b\left(v\mathrm{,}q\right)=\left({\nabla }_w\cdot v\mathrm{,}q\right)$

WG Algorithm. Seek $\left(u_h=\left\{u_0,u_b\right\},p_h\right)\in V_h\times W_h$ satisfying

$\{\begin{array}{l}\left(u_{h\mathrm{,}t}\mathrm{,}v\right)+a\left(u_h\mathrm{,}v\right)-b\left(v\mathrm{,}{p}_h\right)=\left(f\mathrm{,}v\right)\mathrm{,}\forall v\in V_h^0\hfill \\ b\left(u_h\mathrm{,}q\right)=\mathrm{0,\; <mi>\; </mi>\; <mi>\; </mi>}\forall q\in W_h\hfill \\ u_h\left(x\mathrm{,0}\right)={\psi }_0\left(x\right)\hfill \end{array}$ (9)

In the following, the proof process of Lemma 1-6 refers to reference [10] [11] [12] .

Lemma 1. For any $v\in V_h^0$ , the following equation hold true,

${\Vert \left|v\right|\Vert }^2={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\Vert {\nabla }_{w}v\Vert }_T^2+{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{h}_T^{-1}{\Vert Q_b{v}_0-v_b\Vert }_{\partial T}^2$

Lemma 2. For any $v\mathrm{,}w\in V_h^0$ we have

$\left|a\left(v\mathrm{,}w\right)\right|\le \Vert \left|v\right|\Vert \Vert \left|w\right|\Vert$

$a\left(v\mathrm{,}v\right)={\Vert \left|v\right|\Vert }^2\mathrm{.}$

In addition to the projection $Q_h=\left\{Q_0,Q_b\right\}$ defined in the previous section, let $Q_h$ and $Q_h$ be two local L² projections onto $P_{k-1}\left(T\right)$ and ${\left[P_{k-1}\left(T\right)\right]}^{d\times d}$ , respectively.

Lemma 3. The projection operators $Q_h$ , $R_h$ , and $S_h$ satisfy the following commutative properties

${\nabla }_w\left(Q_{h}v\right)=R_h\left(\nabla v\right),\mathrm{<mi>\; </mi>}\forall v\in {\left[H^1\left(\Omega \right)\right]}^d$

${\nabla }_w\cdot \left(Q_{h}v\right)=S_h\left(\nabla \cdot v\right),\mathrm{<mi>\; </mi>}\forall v\in H\left(div\mathrm{<mi>\; </mi>}\Omega \right)$

Lemma 4. There exists a positive constant b independent of h such that

$\underset{v\in V_h^0}{sup}\frac{b\left(v\mathrm{,}\rho \right)}{\Vert \left|v\right|\Vert }\ge \beta \Vert \rho \Vert$

for all $\rho \in W_h$ .

Lemma 5. Poincare inequality of Weak gradient operator: If $v\in V_h^0$ , then exists a constant c satisfying

${\Vert v\Vert }^2\le c{\Vert {\nabla }_{w}v\Vert }^2\le c{\Vert \left|v\right|\Vert }^2$

First of all, we study the existence and uniqueness of the solution for (9). The space defined as follows

$T_h=\left\{v_h\in V_h\left(j\right)\mathrm{<mo>\; ;\; </mo>}\left(q\mathrm{,}{\nabla }_w\cdot v\right)=\mathrm{0,}\forall q\in W_h\right\}\mathrm{.}$

Then we need to seek $u_h\left(x\mathrm{,}t\right)\mathrm{<mo>\; :\; </mo>}\left(\mathrm{0,}T\right)\to T_h$ satisfying

$\{\begin{array}{l}\left(u_{h,t},v\right)+a\left(u_h,v\right)=\left(f,v\right),v\in T_h\hfill \\ u_h\left(x,0\right)={\psi }_0\left(x\right)\hfill \end{array}$ (10)

Let $u_h\left(x\mathrm{,}t\right)$ be the solution of (10) and which is unique, the linear bounded functional $l=l\left(u_h\right)$ on $V_h$ defined as follows.

$\langle l,v\rangle =\left(u_{h,t},v\right)+a\left(u_h,v\right)-\left(f,v\right)$ (11)

Then problem (9) is equivalent to seek $P_h\in W_h$ satisfying

$\left(P_h,{\nabla }_w\cdot v\right)=\langle l,v\rangle ,\forall v\in V_h$ (12)

Using LBB condition and Lax-Milgram Lemma, we know that the solution $P_h\in W_h$ of (12) is unique.

Combing (11) and (12), it is concluded that if initial approximation $u_h\left(x,0\right)={\psi }_0\left(x\right)\in T_h$ , the solution $\left(u_h\mathrm{,}{P}_h\right)\in V_h\times M_h$ of (9) is unique.

In what follows, we introduce Stokes projection, which is the important approximation of projection.

Lemma 6. First of all, we introduce Stokes projection of $\left(u\mathrm{,}p\right)\in X\times W$ , which is $\left(Q_h^{1}u\mathrm{,}{S}_h^{1}p\right)\in V_h\times W_h$ need satisfying

$\{\begin{array}{l}a\left(u_h-Q_h^{1}u,v\right)-b\left(v,p_h-S_h^{1}p\right)=-{\phi }_{u,p}\left(v\right),v\in V_h\\ b\left(u_h-Q_h^{1}u,q\right)=0,\mathrm{<mi>\; </mi>}q\in W_h\end{array}$ (13)

If let $f^{\ast }=-\Delta u+\Delta p$ , easy to know that $\left(Q_h^{1}u\mathrm{,}{S}_h^{1}p\right)\in V_h\times W_h$ satisfying

$\{\begin{array}{l}a\left(Q_h^{1}u\mathrm{,}v\right)-b\left(v\mathrm{,}{S}_h^{1}p\right)=\left(f^{\ast }\mathrm{,}{v}_0\right),v\in V_h\\ b\left(Q_h^{1}u\mathrm{,}q\right)=\mathrm{0,}q\in W_h\end{array}$ (14)

Then $\left(Q_h^{1}u\mathrm{,}{S}_h^{1}p\right)\in V_h\times W_h$ is the finite element approximation of $\left(u\mathrm{,}p\right)\in X\times W$ , so we have

$\{\begin{array}{l}\Vert u_h-Q_h^1\Vert +h\Vert \left|u_h-Q_h^1\right|\Vert \le C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k\right)\\ \Vert p_h-S_h^{1}p\Vert \le C{h}^k\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k\right)\\ \Vert u_h-Q_h^1\Vert +h\Vert \left|{\left(u_h-Q_h^1\right)}_t\right|\Vert \mathrm{<mi>\; </mi>}\le C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\right)\end{array}$ (15)

5. Error Equations

In what follows, we list Lemma 7 to prove the error estimation of approximate solution for Semi-discrete scheme.

We know that $\left(u\mathrm{,}p\right)\in X\times M$ and $\left(u_h\mathrm{,}{p}_h\right)\in X_h\times M_h$ be solution of (1) and Galerkin finite element solution of (9), respectively. The L² projection of u in the finite element space $V_h$ is given by $Q_{h}u=\left\{Q_{0}u,Q_{b}u\right\}$ . Similarly, the pressure p is projected into $W_h$ as $S_{h}p$ . Denote by $e_h$ and ${\epsilon }_h$ the corresponding error given by

$\{\begin{array}{l}{e}_h=\left\{e_0,e_b\right\}=\left\{Q_{0}u-u_0,Q_{b}u-u_b\right\}\\ {\epsilon }_h=S_{h}p-p_h\end{array}$ (16)

Lemma 7. Let $\left(w\mathrm{,}p\right)\in {\left[H^1\left(\Omega \right)\right]}^d\times L^2\left(\Omega \right)$ be sufficiently smooth and satisfy the following equation

$w_t-\Delta w+\nabla \rho =\eta$ (17)

in the domain W. Let $Q_{h}w=\left\{Q_{0}w,Q_{b}w\right\}$ and $S_h\rho$ be the L² projection of $\left(w\mathrm{,}p\right)$ into the finite element space $V_h\times W_h$ . Then, the following equation holds true

$\begin{array}{l}\left(Q_h{w}_t,v_0\right)+\left({\nabla }_w\left(Q_{h}w\right),{\nabla }_{w}v\right)-\left({\nabla }_w\cdot v,S_h\rho \right)\\ =\left(\eta ,v_0\right)+{\mathcal{l}}_w\left(v\right)-{\theta }_{\rho }\left(v\right),\end{array}$ (18)

for all $v\in V_h^0$ . Where ${\mathcal{l}}_w$ and ${\theta }_{\rho }$ are two linear functionals on $V_h^0$ defined by

${\mathcal{l}}_w\left(v\right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\nabla w\cdot n-R_h\left(\nabla w\right)\cdot n\rangle }_{\partial T}$

${\theta }_{\rho }\left(v\right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\left(\rho -S_h\rho \right)n\rangle }_{\partial T}$

Proof. Together Lemma 3, Equation (5) and integration by parts. we obtain

$\begin{array}{l}{\left({\nabla }_w\left(Q_{h}w\right),{\nabla }_{w}v\right)}_T\\ ={\left(R_h\left({\nabla }_w\right),{\nabla }_{w}v\right)}_T\\ =-{\left(v_0,\nabla \cdot R_h\left({\nabla }_w\right)\right)}_T+{\langle v_b,\nabla \cdot R_h\left({\nabla }_w\right)\cdot n\rangle }_{\partial T}\\ ={\left(\nabla v_0,R_h\left({\nabla }_w\right)\right)}_T-{\langle v_0-v_b,\nabla \cdot R_h\left({\nabla }_w\right)\cdot n\rangle }_{\partial T}\\ ={\left(\nabla w,\nabla v\right)}_T-{\langle v_0-v_b,\nabla \cdot R_h\left({\nabla }_w\right)\cdot n\rangle }_{\partial T}\end{array}$ (19)

Next, Combing Lemma 3 and Equation (8), the fact that ${\displaystyle {\sum }_{T\in {\mathcal{T}}_h}}{\langle v_b\mathrm{,}pn\rangle }_{\partial T}=0$ ,

then using integration by parts, we obtain

$\begin{array}{l}\left({\nabla }_w\cdot v\mathrm{,}{S}_h\rho \right)\\ =-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0\mathrm{,}\nabla \left(S_h\rho \right)\rangle }_T+{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_b\mathrm{,}\left(S_h\rho \right)n\rangle }_{\partial T}\\ ={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle \nabla \cdot v_0\mathrm{,}{S}_h\rho \rangle }_T-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\left(S_h\rho \right)n\rangle }_{\partial T}\\ ={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle \nabla \cdot v_0\mathrm{,}\rho \rangle }_T-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\left(S_h\rho \right)n\rangle }_{\partial T}\end{array}$

$\begin{array}{l}=-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0\mathrm{,}\nabla \rho \rangle }_T+{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0\mathrm{,}\rho n\rangle }_{\partial T}-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\left(S_h\rho \right)n\rangle }_{\partial T}\\ =-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0\mathrm{,}\nabla \rho \rangle }_T+{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\rho n\rangle }_{\partial T}-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\left(S_h\rho \right)n\rangle }_{\partial T}\\ =-\left(v_0\mathrm{,}\nabla \rho \right)+{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\left(\rho -S_h\rho \right)n\rangle }_{\partial T}\end{array}$

We can imply that

$\left(v_0\mathrm{,}\nabla p\right)=-\left({\nabla }_w\cdot v\mathrm{,}{S}_h\rho \right)+{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\left(\rho -S_h\rho \right)n\rangle }_{\partial T}$ (20)

Next, we test (17) by using $v_0$ in $v=\left\{v_0,v_b\right\}\in V_h^0$ to obtain, we can obtain

$\left(w_t\mathrm{,}{v}_0\right)-\left(\Delta w\mathrm{,}{v}_0\right)+\left(\nabla \rho \mathrm{,}{v}_0\right)=\left(\eta \mathrm{,}{v}_0\right)$ (21)

It follows from the usual integration by parts that

$-\left(\Delta w\mathrm{,}{v}_0\right)={\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\left(\nabla w\mathrm{,}\nabla v_0\right)}_T-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\nabla w\cdot n\rangle }_{\partial T}$

Where we have used the fact that ${\displaystyle {\sum }_{T\in {\mathcal{T}}_h}}{\langle v_b\mathrm{,}\nabla w\cdot n\rangle }_{\partial T}=0$ . using Equations (19) and (20), we have

$-\left(\Delta w\mathrm{,}{v}_0\right)=\left({\nabla }_w\left(Q_{h}w\right)\mathrm{,}{\nabla }_{w}v\right)-{\displaystyle \underset{T\in {\mathcal{T}}_h}{\sum }}{\langle v_0-v_b\mathrm{,}\nabla w\cdot n-R_h\left(\nabla w\right)\cdot n\rangle }_{\partial T}$ (22)

Substituting (20), (22) and $\left(Q_h{w}_t\mathrm{,}{v}_0\right)=\left(w_t\mathrm{,}{v}_0\right)$ into (21) yields

$\left(Q_h{w}_t\mathrm{,}{v}_0\right)+\left({\nabla }_w\left(Q_{h}w\right)\mathrm{,}{\nabla }_{w}v\right)-\left({\nabla }_w\cdot v\mathrm{,}{S}_h\rho \right)=\left(\eta \mathrm{,}{v}_0\right)+{\mathcal{l}}_w\left(v\right)-{\theta }_{\rho }(v)$

which completes the proof of the lemma.

In what follows, we give the derivation of the error equation of (9).

Lemma 8. Let $e_h$ and ${\epsilon }_h$ be the error of the weak Galerkin finite element solution arising from (9), as defined by (16). Then, we have

$\{\begin{array}{l}\left(e_{h,t},v\right)+a\left(e_h,v\right)-b\left(v,{\epsilon }_h\right)={\phi }_{u,p}\left(v\right)\hfill \\ b\left(e_h,q\right)=0,\hfill \end{array}$ (23)

for all $v\in V_h^0$ and $q\in W_h$ , where ${\phi }_{u\mathrm{,}p}\left(v\right)={\mathcal{l}}_u\left(v\right)-{\theta }_p\left(v\right)+s\left(Q_{h}u\mathrm{,}v\right)$ is a linear functional defined on $V_h^0$ .

Proof. Since $\left(u\mathrm{,}p\right)$ satisfies the Equation (17) with $\eta =f$ , then from Lemma 6 we have

$\left(Q_h{u}_t,v_0\right)+\left({\nabla }_w\left(Q_{h}u\right),{\nabla }_{w}v\right)-\left({\nabla }_w\cdot v,S_{h}p\right)=\left(f,v_0\right)+{\mathcal{l}}_u\left(v\right)-{\theta }_{\rho }(v)$

Adding $s\left(Q_{h}u\mathrm{,}v\right)$ to both side of the above equation give

$\left(Q_h{u}_t,v_0\right)+a\left(Q_{h}u,v\right)-b\left(v,S_{h}p\right)=\left(f,v_0\right)+{\mathcal{l}}_u\left(v\right)-{\theta }_{\rho }\left(v\right)+s\left(Q_{h}u,v\right)$ (24)

The difference of (24) and (9) yields the following equation,

$\left(e_{h,t},v_0\right)+a\left(e_h,v\right)-b\left(v,{\epsilon }_h\right)={\mathcal{l}}_u\left(v\right)-{\theta }_{\rho }\left(v\right)+s\left(Q_{h}u,v\right)$

for all $v\in V_h^0$ , where $e_h=\left\{e_0,e_b\right\}=\left\{Q_{0}u-u_0,Q_{b}u-u_b\right\}$ . This completes the derivation of (23).

As to (24), we test Equation (1) by $q\in W_h$ and use (9) to obtain

$0=\left(\nabla \cdot u,q\right)=\left({\nabla }_w\cdot Q_{h}u,q\right)=b\left(Q_{h}u,q\right)$ (25)

The difference of (25) and (9) yields the following equation

$b\left(e_h,q\right)=0$

for all $q\in W_h$ .

Which completes the proof of the lemma.

In the following, the proof process of Lemma 9 refers to reference [10] .

Lemma 9. If $\left(w\mathrm{,}\rho \mathrm{,}r\right)\in {\left[H^{r+1}\left(\Omega \right)\right]}^d\times H^r\left(\Omega \right)\times V_h$ and $1\le r\le k$ , with the

precondition of regular-shape ${\mathcal{T}}_h$ , we have the following estimation.

$\left|s\left(Q_{h}w\mathrm{,}v\right)\right|\le C{h}^r{\Vert w\Vert }_{r+1}\Vert \left|v\right|\Vert$

$\left|{\mathcal{l}}_w\left(v\right)\right|\le C{h}^r{\Vert w\Vert }_{r+1}\Vert \left|v\right|\Vert$

$\left|{\theta }_{\rho }\left(v\right)\right|\le C{h}^r{\Vert w\Vert }_{r+1}\Vert \left|v\right|\Vert$

6. Error Estimates

The following theorem is the main result of this paper.

Theorem 1. Let $\left(u\mathrm{,}p\right)\in {\left[H_0^1\left(\Omega \right)\cap H^{k+1}\left(\Omega \right)\right]}^d\times \left[L_0^2\left(\Omega \right)\cap H^k\left(\Omega \right)\right]$ and $\left(u_h\mathrm{,}{p}_h\right)\in V_h\times W_h$ be the solution of (1) and (9), respectively. the following error estimates is true.

$\begin{array}{l}\Vert Q_{h}u-u_h\Vert \le C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\displaystyle {\int }_0^t}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\right)\text{d}\tau \right)\\ \Vert Q_{h}u-u_h\Vert \le C{h}^k\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\right)\\ \Vert S_{h}p-p_h\Vert \le C{h}^{k+1}({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\\ \mathrm{<mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>}+{\displaystyle {\int }_0^t}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k+{\Vert u_{tt}\Vert }_{k+1}+{\Vert p_{tt}\Vert }_k\right)\text{d}\tau )\end{array}$

Proof. Let

$e_h=Q_{h}u-u_h=Q_{h}u-Q_h^{1}u+Q_h^{1}u-u_h=\theta +\eta$

$e_h\left(\cdot ,0\right)=\theta \left(\cdot ,0\right)=\eta \left(\cdot ,0\right)=0$

By the error of Equation (23), we have

$\begin{array}{l}\left({\theta }_t\mathrm{,}v\right)+a\left(\theta \mathrm{,}v\right)-b\left(v\mathrm{,}{S}_{h}p-S_h^{1}p\right)\\ ={\phi }_{u\mathrm{,}p}\left(v\right)-\left({\eta }_t\mathrm{,}v\right)-a\left(\eta \mathrm{,}v\right)+b\left(v\mathrm{,}{S}_h^{1}p-p_h\right)\end{array}$ (26)

Substituting (13) into (26), we obtain

$\left({\theta }_t,v\right)+a\left(\theta ,v\right)-b\left(v,S_{h}p-S_h^{1}p\right)=-\left({\eta }_t,v\right)$ (27)

Let $v=\theta =Q_{h}u-Q_h^{1}u$ , combing the Equation (25) and (14), we have

$b\left(\theta ,S_{h}p-S_h^{1}p\right)=0$

That is

$\left({\theta }_t,\theta \right)+a\left(\theta ,\theta \right)=-\left({\eta }_t,\theta \right)$

By Lemma 2 and Cauchy inequality, we have

$\frac{1}{2}\frac{\text{d}}{\text{d}t}{\Vert \theta \Vert }^2+{\Vert \left|\theta \right|\Vert }^2\le \Vert {\eta }_t\Vert \Vert \theta \Vert \le {\Vert {\eta }_t\Vert }^2+{\Vert \theta \Vert }^2$ (28)

By Gronwall Lemma, we have

$\frac{1}{2}\frac{\text{d}}{\text{d}t}{\Vert \theta \Vert }^2+{\Vert \left|\theta \right|\Vert }^2\le {\left(C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\right)\right)}^2$ (29)

By Cauchy inequality, we have

$\Vert \left|\theta \right|\Vert \le C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\right)$ (30)

Then take the integration about t of both side of Equation (28)

${\Vert \theta \left(\cdot \mathrm{,}t\right)\Vert }^2+2{\displaystyle {\int }_0^t}{\Vert \left|\theta \right|\Vert }^2\text{d}\tau \le {\Vert \theta \left(\cdot \mathrm{,0}\right)\Vert }^2+C{\left({\displaystyle {\int }_0^t}\Vert {\eta }_t\Vert \text{d}\tau \right)}^2+\frac{1}{4}\underset{\tau \le t}{{\displaystyle sup}}{\Vert \theta \left(\tau \right)\Vert }^2$

Since $\Vert \theta \left(\cdot ,0\right)\Vert =0$ , then

$\begin{array}{c}\Vert \theta \Vert \le \underset{\tau \le t}{{\displaystyle sup}}\Vert \theta \left(\tau \right)\Vert \le C{\displaystyle {\int }_0^t}\Vert {\eta }_t\Vert \text{d}\tau \\ \le C{h}^{k+1}{\displaystyle {\int }_0^t}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\right)\text{d}\tau \end{array}$ (31)

Combing the Equations (15), (29), (30) and triangle inequality, we have

$\Vert Q_{h}u-u_h\Vert \le C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\displaystyle {\int }_0^t}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\right)\text{d}\tau \right)$ (32)

$\Vert Q_{h}u-u_h\Vert \le C{h}^k\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k+{\Vert u_t\Vert }_{k+1}+{\Vert p_t\Vert }_k\right)$ (33)

Next, we proof the error estimate of pressure approximation $\Vert S_{h}p-p_h\Vert$ , by using error Equation (23), we have

$b\left(v\mathrm{,}{S}_{h}p-p_h\right)=\left(e_{h\mathrm{,}t}\mathrm{,}v\right)+a\left(e_h\mathrm{,}v\right)-{\phi }_{u\mathrm{,}p}(v)$

By using Lemma 2, Lemma 5 and Lemma 9, we obtain

$\begin{array}{c}b\left(v\mathrm{,}{S}_{h}p-p_h\right)\le \Vert e_{h\mathrm{,}t}\Vert \Vert v\Vert +\Vert \left|e_h\right|\Vert \Vert \left|v\right|\Vert +C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k\right)\Vert \left|v\right|\Vert \\ \le C\Vert e_{h\mathrm{,}t}\Vert \Vert v\Vert +\Vert \left|e_h\right|\Vert \Vert \left|v\right|\Vert +C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k\right)\Vert \left|v\right|\Vert \end{array}$

By Lemma 4, we have

$\Vert S_{h}p-p_h\Vert \le C\Vert e_{h\mathrm{,}t}\Vert +\Vert \left|e_h\right|\Vert +C{h}^{k+1}\left({\Vert u\Vert }_{k+1}+{\Vert p\Vert }_k\right)$ (34)