- Email: [email protected]
A fully discrete stabilized mixed finite volume element formulation for the non-stationary conduction–convection problem
J. Math. Anal. Appl. 404 (2013) 71–85
Contents lists available at SciVerse ScienceDirect
Journal of Mathematical Analysis and Applications journal homepage: www.elsevier.com/locate/jmaa
A fully discrete stabilized mixed finite volume element formulation for the non-stationary conduction–convection problem✩ Zhendong Luo a,∗ , Hong Li b,∗ , Ping Sun c a
School of Mathematics and Physics, North China Electric Power University, Beijing 102206, China
b
School of Mathematical Sciences, Inner Mongolia University, Huhhot 010021, China
c
School of Mathematics and Computer Science, Guizhou Normal University, Guiyang 550001, China
article
info
Article history: Received 3 January 2013 Available online 13 March 2013 Submitted by G. Chen Keywords: Non-stationary conduction–convection problem Stabilized mixed finite volume element formulation Local Gauss integrals and parameter-free Error estimate
abstract At first, a semi-discrete formulation with respect to time for the non-stationary conduction–convection problem is recalled. Then, a fully discrete stabilized mixed finite volume element (SMFVE) formulation based on two local Gauss integrals and parameterfree is established directly from the semi-discrete formulation with respect to time for the non-stationary conduction–convection problem. Following this, the error estimates for the fully discrete SMFVE solutions are derived by means of the standard mixed finite element method. Finally, some numerical experiments are presented illustrating that the numerical errors are consistent with theoretical results, the degrees of freedom of the fully discrete SMFVE formulation are far fewer than those of the classical finite volume element (FVE) formulation without any stabilization, and its numerical solutions are more stable than those of the classical FVE formulation without any stabilization, thus validating that the fully discrete SMFVE formulation is feasible and efficient for finding the numerical solutions of the non-stationary conduction–convection problem. © 2013 Elsevier Inc. All rights reserved.
1. Introduction Let Ω ⊂ R2 be a bounded and connected polygonal domain. Consider the non-stationary conduction–convection problem including the velocity vector field and the pressure and temperature fields (see [11,17,18]). Problem I. Find u = (u1 , u2 ), p, and T such that, for tN > 0,
ut − ν ∆u + (u · ∇)u + ∇ p = T j , ∇ · u = 0 , Tt − γ0−1 ∆T + (u · ∇)T = 0, u (x, y, t ) = 0, T (x, y, t ) = ϕ(x, y, t ), u(x, y, 0) = 0, T (x, y, 0) = ψ(x, y)
(x, y, t ) ∈ Ω × (0, tN ), (x, y, t ) ∈ Ω × (0, tN ), (x, y, t ) ∈ Ω × (0, tN ), (x, y, t ) ∈ ∂ Ω × (0, tN ), (x, y) ∈ Ω ,
(1)
✩ This work is jointly supported by National Science Foundation of China (11271127, 11061009, and 11061021), Natural Science Foundation of Inner Mongolia (2012MS0106), Science Research Program of Guizhou (GJ[2011]2367), Science Research Program of Inner Mongolia Advanced Education (NJ10006), and Special Funds for Co-construction Project of Beijing and North China Electric Power University. ∗ Corresponding authors. E-mail addresses: [email protected] (Z. Luo), [email protected] (H. Li).
0022-247X/$ – see front matter © 2013 Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.jmaa.2013.03.001
72
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
where u = (u1 , u2 ) is the fluid velocity vector, p the pressure, T the temperature, tN the total time, j = (0, 1) the unit √ √ vector, ν = Pr /Re, Re the Reynolds number, Pr the Prandtl number, γ0 = RePr, and ϕ(x, y, t ) and ψ(x, y) are two given functions. For the sake of convenience and without loss of generality, we may as well suppose that ϕ(x, y, t ) = 0 in the following theoretical analysis. Since Problem I constitutes a system of nonlinear PDEs including the velocity vector field and the pressure and the temperature fields, there are no analytical solutions in general. One has to rely on numerical solutions (see, e.g., [9,17,18]). However, most of the existing papers use either the finite element (FE) method or finite difference (FD) schemes as discretization tools. The finite volume element (FVE) method (see [6,13,21]) is considered as one of the most effective numerical methods due to the following three advantages. First, it preserves the integral invariants of conservations of total energy and mass. Second, it is suitable for computations with complicated boundary and has higher accuracy than that of the FD method. Third, it has the same accuracy as FE method but is much more convenient and simpler than the FE method. It is also known as a box method (see [3]) or a generalized difference method (see [14,15]). It has been widely used to find numerical solutions of various types of PDEs, for example, second order elliptic equations, parabolic equations, Stokes equations, and the Navier–Stokes equations (see, e.g., [2–4,6–8,12–15,21–23]). Though a fully discrete FVE formulation without any stabilization for the non-stationary conduction–convection problem is proposed in [16], its trial function spaces of velocity and pressure have to satisfy an important convergence stability condition, i.e., the Babuška–Brezzi (BB) inequality (see [16,18]). Thus, there exist many difficulties for its theoretical analysis and actual applications. Therefore, it is very necessary to improve the method in [16]. To this end, we devote to establishing a fully discrete SMFVE formulation based on two local Gauss integrals and parameter-free for the nonstationary conduction–convection problem in this paper such that its trial function spaces of velocity and pressure satisfy the BB condition naturally. In addition, an optimizing reduced Petrov–Galerkin least squares mixed FE formulation based on residuals (see [19]) and a reduced mixed FE formulation without any stabilization (see [20]) have been established for the non-stationary conduction–convection problem, but they are completely different from the fully discrete SMFVE formulation in this paper since the FVE method is entirely different from and far more advantageous than the FE method as is mentioned above. The paper is organized as follows. In Section 2, we first recall the semi-discrete formulation with respect to time for the non-stationary conduction–convection problem, and then establish directly a fully discrete SMFVE formulation from the semi-discrete formulation with respect to time which could avoid the semi-discrete SMFVE formulation about spatial variables. It is a new type of study attempt for the non-stationary conduction–convection problem. In Section 3, the error estimates between the fully discrete SMFVE solutions and the exact solution are provided. In Section 4, some numerical experiments are presented illustrating that the errors of the fully discrete SMFVE solutions are consistent with the theoretical results obtained previously, the degrees of freedom of the fully discrete SMFVE formulation are far fewer than those of the classical fully discrete FVE formulation without any stabilization (see [16]), and the numerical solutions for the fully discrete SMFVE formulation are more stable than those for the classical fully discrete FVE formulation without any stabilization in [16], thus validating that the SMFVE formulation is feasible and efficient for finding the numerical solutions of the nonstationary conduction–convection problem. Section 5 provides main conclusions and some discussions. 2. Semi-discrete formulation about time and fully discrete SMFVE formulation for Problem I 2.1. Recall semi-discrete formulation about time for Problem I Sobolev spaces along with their properties used in this context are standard (see [1]). For example, H m (Ω ) = The 2 α α1 α2 v; 06|α|6m Ω |∂ v/(∂ x ∂ y )|0 dxdy < ∞ (where m > 0 is integer, α = (α1 , α2 ), α1 and α2 are two non-negative 1/2 α α1 α2 2 integers, and |α| = α1 + α2 ) with norm and semi-norm denoted by ∥v∥m = and 06|α|6m Ω |∂ v/(∂ x ∂ y )|0 dxdy 1/2 t α α1 α2 2 , respectively; L2 (L2 ) = v(x, y, t ); 0 N Ω v 2 (x, y, t )dxdydt < ∞ with |v|m = |α|=m Ω |∂ v/(∂ x ∂ y )|0 dxdy 1/2 t norm ∥v∥L2 (L2 ) = 0 N Ω v 2 dxdydt ; H01 (Ω ) = v ∈ H 1 (Ω ); v|∂ Ω = 0 . Let X = H01 (Ω )2 , M = L20 (Ω ) = q ∈ L2 (Ω ); Ω qdxdy = 0 , and W = H01 (Ω ). Put a(u, v ) = ν ∇ u · ∇ vdxdy, ∀u, v ∈ X , b(q, v ) = q divvdxdy, ∀v ∈ X , q ∈ M , Ω Ω 1 [(u∇ v ) · w − (u∇ w ) · v] dxdy, ∀u, v , w ∈ X , a1 (u, v , w ) = 2 Ω 1 [(u · ∇ T )φ − (u · ∇φ)T ] dxdy, ∀u ∈ X , ∀T , φ ∈ W , a2 (u, T , φ) = 2
Ω
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
d(T , φ) = γ0−1
Ω
∇ T · ∇φ dxdy,
73
∀T , φ ∈ W .
Thus, a mixed variational formulation for Problem I is written as follows. Problem II. Find (u, p, T ) ∈ H 1 (0, tN ; X )2 × L2 (0, tN ; M ) × H 1 (0, tN ; W ) such that, for almost all t ∈ (0, tN ),
(u , v ) + a(u, v ) + a1 (u, u, v ) − b(p, v ) = (T j , v ), t b(q, u) = 0, (Tt , φ) + d(T , φ) + a2 (u, T , φ) = 0, u(x, y, 0) = 0, T (x, y, 0) = ψ(x, y),
∀v ∈ X , ∀q ∈ M , ∀φ ∈ W , ( x, y ) ∈ Ω ,
(2)
where (·, ·) denotes inner product in L2 (Ω )2 or L2 (Ω ). There hold the following properties for above trilinear forms a1 (·, ·, ·) and a2 (·, ·, ·) (see [16,18]): a1 (u, v , w ) = −a1 (u, w , v ),
a1 (u, v , v ) = 0,
(3)
a2 (u, T , φ) = −a2 (u, φ, T ),
a2 (u, φ, φ) = 0,
∀u, v , w ∈ X , ∀u ∈ X , ∀T , φ ∈ W . There are the following properties for above bilinear forms a(·, ·), d(·, ·), and b(·, ·) (also see [16,18]):
(4)
a(v , v ) > ν|v |21 ,
sup
b(q, v )
v ∈X
2 1
∀u, v ∈ X ,
|d(T , φ)| 6 γ0 |T |1 |φ|1 , −1
∀φ ∈ W ;
> β∥q∥0 ,
|v |1
|a(u, v )| 6 ν|u|1 |v |1 ,
∀v ∈ X ;
d(φ, φ) > γ0 |φ| , −1
(5)
∀T , φ ∈ W ,
(6)
∀q ∈ M ,
(7)
where β is a constant. Define N0 =
sup
u,v ,w ∈X
a1 (u, v , w ) , |u|1 · |v |1 · |w |1
N˜ 0 =
sup
u∈X ,(T ,φ)∈W ×W
a2 (u, T , φ) . |u|1 · |T |1 · |φ|1
(8)
The following result is classical (see [16–18]). Theorem 1. If ψ ∈ L2 (Ω ), then the Problem II has at least a solution which, in addition, is unique provided that ∥ψ∥20 6
2ν 2 tN /(2N0 tN−1 exp(tN ) + νγ0 N˜ 02 ), and there are the following prior estimates
∥u∥20 + ν∥∇ u∥2L2 (L2 ) 6 tN2 ∥ψ∥20 exp(tN ),
∥T ∥20 + γ0−1 ∥∇ T ∥2L2 (L2 ) 6 ∥ψ∥20 .
Let N be the positive integer, k = tN /N denote the time step increment, tn = nk (0 6 n 6 N), and (un , pn , T n ) be the semi-discrete approximation of (u, p, T ) at tn = nk (n = 0, 1, . . . , N ) about time. If the differential quotients ut and Tt in Problem II at time t = tn are, respectively, approximated by means of the backward difference quotients ∂¯t un = (un − un−1 )/k and ∂¯t T n = (T n − T n−1 )/k, then the semi-discrete approximation scheme of Euler backward one step of Problem II with respect to time is written as follows. Problem III. Find (un , pn , T n ) ∈ X × M × W (1 6 n 6 N ) such that
(un , v ) + ka(un , v ) + ka1 (un−1 , un , v ) − kb(pn , v ) = k(T n j , v ) + (un−1 , v ), b(q, un ) = 0, (T n , φ) + kd(T n , φ) + ka2 (un−1 , T n , φ) = (T n−1 , φ), 0 u = 0, T 0 = ψ(x, y),
∀v ∈ X , ∀q ∈ M , ∀φ ∈ W , (x, y) ∈ Ω .
(9)
It has been proved that Problem III has a unique series of solutions {un , pn , T n }Nn=1 (see [16–18]) and the following results for Problem III are known and classical (see Theorems 5.6 and 5.7 in [18] or [16]). Theorem 2. Under the assumptions of Theorem 1, Problem III has a unique series of solutions (un , pn , T n ) ∈ X × M × W such that
∥un ∥20 + ν k
n
∥∇ ui ∥20 6 C0 ∥ψ∥20 , ∥T n ∥20 + γ0−1 k
i=1
∥u(tn ) −
∥ + kν
un 20
n
n
∥∇ T i ∥20 6 ∥ψ∥20 ,
∥pn ∥0 6 C ∥ψ∥0 ,
i=1
|u(ti ) − | + ∥p(tn ) − ui 21
∥ 6 Ck2 ∥ψ∥20 ,
pn 20
i=1
∥T (tn ) − T n ∥20 + kγ0−1
n
|T (ti ) − T i |21 6 Ck2 ∥ψ∥20 ,
1 6 n 6 N,
i =1
where C0 = 2nk exp(2nk) is a bounded real number, (u, p, T ) ∈ [H01 (Ω ) ∩ H 2 (Ω )]2 × [H 1 (Ω ) ∩ M ] × [H01 (Ω ) ∩ H 2 (Ω )] is the exact solution for the Problem I, and C is a constant independent of k.
74
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Fig. 1. Left chart is a triangle K partitioned into three sub-regions Kz . Right chart is a sample region with dotted lines indicating the corresponding control volume Vz .
2.2. Fully discrete SMFVE formulation for Problem I In order to establish the fully discrete SMFVE formulations for Problem II, it is necessary to introduce an FVE approximation for the spatial variables of Problem III (more details see [2,3,14,15]). Firstly, let ℑh = {K } be a quasi-uniform triangulation of Ω with h = max hK , where hK is the diameter of the triangle K ∈ ℑh (see [5,10,18]). In order to describe the SMFVE formulation, we introduce a dual partition ℑ∗h based on ℑh whose elements are called the control volumes. We construct the control volume in the same way as in [2,3,14,15]. Let zK = (xK , yk ) be the barycenter of K ∈ ℑh . We connect zK with line segments to the midpoints of the edges of K , thus partitioning K into three quadrilaterals Kz (z = (x, y) ∈ Zh (K ), where Zh (K ) are the vertices of K ). Then with each vertex z ∈ Zh = K ∈ℑh Zh (K ) we associate a control volume Vz , which consists of the union of the sub-regions Kz , sharing the vertex z. Finally, we obtain a group of control volumes covering the domain Ω , which is called a barycenter-type dual partition ℑ∗h of the triangulation ℑh (see Fig. 1). We denote the set of interior vertices of Zh by Zh◦ . Since the FE triangulation ℑh is quasi-uniform, the dual partition ℑ∗h is also quasi-uniform (see [2,3,5,10,14,15,18]), namely there exist two positive constants C1 and C2 , being independent of the spatial mesh size h and temporal mesh size k, such that C1 h2 6 mes(Vz ) 6 C2 h2 ,
∀Vz ∈ ℑ∗h .
(10)
Moreover, the barycenter-type dual partition can lead to relatively simple calculations. To this end, the trial function spaces Xh , Wh , and Mh of velocity, temperature, and pressure are, respectively, defined as follows: Xh = vh ∈ X ∩ C (Ω )2 ; vh |K ∈ [P1 (K )]2 , ∀K ∈ ℑh , Wh = wh ∈ W ∩ C (Ω ); wh |K ∈ P1 (K ), ∀K ∈ ℑh , Mh = {qh ∈ M ; qh |K ∈ P1 (K ), ∀K ∈ ℑh } ,
where P1 (K ) is the linear function space on K . Note that they are different from those in [16]. It is obvious that Xh ⊂ X = H01 (Ω )2 and Wh ⊂ W = H01 (Ω ). For (u, T ) ∈ X × W , let (Πh u, ρh T ) be the interpolation projection of (u, T ) onto the trial function spaces Xh × Wh . Then, due to the interpolation theory of Sobolev spaces (see [5,10,14,16,18]), there hold the following error estimates
|u − Πh u|m 6 Ch2−m |u|2 ,
∀u ∈ H 2 (Ω )2 , m = 0, 1,
(11)
|T − ρh T |m 6 Ch
∀T ∈ H (Ω ), m = 0, 1,
(12)
2−m
|T |2 ,
2
where C in this context indicates a positive constant which is possibly different at different occurrences, being independent of the spatial mesh size h and temporal mesh size k. ˜ h of velocity and temperature are, respectively, chosen as follows: The test spaces X˜ h and W X˜ h = vh ∈ L2 (Ω )2 ; vh |Vz ∈ [P0 (Vz )]2 (Vz ∩ ∂ Ω = ∅), vh |Vz = 0 (Vz ∩ ∂ Ω ̸= ∅), ∀Vz ∈ ℑ∗h , ˜ h = wh ∈ L2 (Ω ); wh |Vz ∈ P0 (Vz ) (Vz ∩ ∂ Ω = ∅), wh |Vz = 0 (Vz ∩ ∂ Ω ̸= ∅), ∀Vz ∈ ℑ∗h , W
(13)
where P0 (Vz ) is the constant function space on Vz . In fact, they can be spanned by the following basis functions
φz (x, y) =
1, (x, y) ∈ Vz , 0, elsewhere,
∀z ∈ Zh◦ .
(14)
˜ h , i.e., For (u, w) ∈ X × W , let (Πh∗ u, ρh∗ w) be the interpolation projection of (u, w) onto the test spaces X˜ h × W (Πh∗ u, ρh∗ w) =
(u(z ), w(z ))φz .
(15)
z ∈Zh◦
By the interpolation theory of Sobolev spaces (see [5,10,14,16,18]), we have that
∥u − Πh∗ u∥0 6 Ch|u|1 ;
∥w − ρh∗ w∥0 6 Ch|w|1 .
(16)
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
75
Due to the loss of continuity of the piecewise constant vector functions in X˜ h and the piecewise constant functions in ˜ h on the boundary of two neighboring elements, namely X˜ h ̸⊂ Xh and W ˜ h ̸⊂ Wh , the bilinear forms a(u, v ), b(v , p), and W d(T , φ) and trilinear forms a1 (u, v , w ) and a2 (u, T , φ) in (2) must be revised. Noting that divu = 0, by Green’s formula, there hold the following formulas
Ω
∆u · vdxdy =
∇ p · vdxdy = Ω
Ω
Ω
Ω
∆u · vdxdy = Vz
Vz ∈ℑ∗ h
pv · nds −
∆T φ dxdy =
∆T φ dxdy =
Vz ∈ℑ∗ h
(u · ∇)T φ dxdy =
Vz ∈ℑ∗ h
(v ∇ u) · nds −
∇ u · ∇ vdxdy,
(17)
Vz
Vz ∈ℑ∗ h
pdivvdxdy,
(18)
Vz
Vz ∈ℑ∗ h
(u · ∇)v · wdxdy =
Vz ∈ℑ∗ h
Vz
Vz ∈ℑ∗ h
∂ Vz
Vz ∈ℑ∗ h
∂ Vz
Vz ∈ℑ∗ h
∂ Vz
(∇ T φ) · nds −
Vz ∈ℑ∗ h
(u · ∇)v · wdxdy = Vz
Vz ∈ℑ∗ h
(u · ∇)T φ dxdy = Vz
∂ Vz
Vz ∈ℑ∗ h
∂ Vz
∇ T ∇φ dxdy,
(19)
Vz
(u · n)(v · w )ds −
Vz ∈ℑ∗ h
T φ(u · n)ds −
Vz ∈ℑ∗ h
(u∇ w ) · vdxdy, Vz
(u · ∇φ)T dxdy, Vz
(20)
(21)
where ∂ Vz denotes the line integrals, with the counter clockwise direction, on the boundary ∂ Vz of the dual element; n = (n1 , n2 ) is the unit outer normal vector to ∂ Vz . So bilinear forms a(u, v ), b(p, v ), and d(T , ϕ) and trilinear forms a1 (u, v , w ) and a2 (u, T , φ) in (2) are, respectively, rewritten as
a(u, v ) = ν
∇ u · ∇ vdxdy − ∂ Vz
Vz
Vz ∈ℑ∗ h
b(p, v ) = −
d(T , φ) = γ0
−1
∇ T ∇φ dxdy − γ0
−1
Vz ∈ℑ∗ h
(u∇ w ) · vdxdy + Vz
Vz ∈ℑ∗ h
a2 (u, v , φ) = −
pdivvdxdy ,
Vz
Vz ∈ℑ∗ h
Vz ∈ℑ∗ h
(22)
(23)
Vz
a1 (u, v , w ) = −
(v ∇ u) · nds ,
pv · nds −
∂ Vz
Vz ∈ℑ∗ h
(u · ∇φ)T dxdy +
(∇ T φ) · nds,
(24)
(u · n)(v · w )ds,
(25)
T φ(u · n)ds.
(26)
Vz ∈ℑ∗ h
Vz
∂ Vz
∂ Vz
Vz ∈ℑ∗ h
∂ Vz
˜ h are, respectively, the piecewise constant vector function space and the piecewise constant function space Since X˜ h and W with the characteristic functions of the dual elements Vz as the basis functions, there hold the following expressions a(u, v ) = −ν
Vz ∈ℑ∗ h
a1 (u, v , w ) =
(u · n)(v · w )ds,
∂ Vz
Vz ∈ℑ∗ h
d(T , φ) = −γ0−1
T φ(u · n)ds,
∂ Vz
Vz ∈ℑ∗ h
∀v ∈ X˜ h ;
b(p, v ) = −
Vz ∈ℑ∗ h
Vz ∈ℑ∗ h
a2 (u, T , φ) =
∂ Vz
(v ∇ u) · nds,
∂ Vz
φ∇ T · nds,
∂ Vz
pv · nds,
∀v ∈ X˜ h ;
(27)
∀w ∈ X˜ h ; (28)
˜ h; ∀φ ∈ W ˜ h. ∀φ ∈ W
Then a fully discrete SMFVE formulation for Problem II is written as follows.
(29)
76
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Problem IV. Find (unh , pnh , Thn ) ∈ Xh × Mh × Wh (1 6 n 6 N ) such that
(∂¯ un , Π ∗ v ) + ah (unh , Πh∗ vh ) + a1h (unh−1 , unh , Πh∗ vh ) + bh (pnh , Πh∗ vh ) = (Thn j , Πh∗ vh ), t h nh h b(qh , uh ) + Dh (pnh , qh ) = 0, (∂¯t T n , ρ ∗ wh ) + dh (Thn , ρh∗ wh ) + a2h (unh−1 , Thn , ρh∗ wh ) = 0, 0 h h0 uh = 0, Th = ρh ψ(x, y),
∀vh ∈ Xh , ∀qh ∈ Mh , ∀wh ∈ Wh , (x, y) ∈ Ω ,
(30)
where ah (unh , Πh∗ vh ) = −ν
∂ Vz
Vz ∈ℑ∗ h
(vh (z )∇ unh ) · nds;
dh (Thn , ρh∗ wh ) = −γ0−1
∂ Vz
wh (z )
wh ( z )
K ,2
K ∈ℑh
∂ Vz
qh nds;
(31)
(32)
∇ T · nds; ∂ Vz
Vz ∈ℑ∗ h
Dh (pnh , qh ) = ε
vh ( z )
(unh−1 · n)(unh · vh (z ))ds;
Vz ∈ℑ∗ h
Vz ∈ℑ∗ h
a2h (unh−1 , Thn , ρh∗ wh ) =
Vz ∈ℑ∗ h
a1h (unh−1 , unh , Πh∗ vh ) =
bh (qh , Πh∗ vh ) =
∂ Vz
pnh qh dxdy −
Thn (unh−1 · n)ds;
(33)
K ,1
pnh qh dxdy ,
(34)
here ε is a positive real number and K ,i g (x, y)dxdy indicate an appropriate Gauss integral on K that is exact for polynomials of degree i (i = 1, 2) and g (x, y) = ph qh is a polynomial of degree not more than i (i = 1, 2). Thus, for all test functions qh ∈ Mh , the trial function ph ∈ Mh must be piecewise constant when i = 1. Consequently, we define the L2 -projection ˆ h such that ∀p ∈ L2 (Ω ) satisfies operator ϱh : L2 (Ω ) → W
(p, qh ) = (ϱh p, qh ),
ˆ h, ∀qh ∈ W
(35)
ˆ h ⊂ L2 (Ω ) denotes the piecewise constant space associated with ℑh . The projection operator ϱh has the following where W properties (see [18,5]). ∥ϱh p∥0 6 C ∥p∥0 ,
∀p ∈ L2 (Ω ),
∥p − ϱh p∥0 6 Ch∥p∥1 ,
(36)
∀p ∈ H (Ω ). 1
(37)
Now, using the definition of ϱh , we can rewrite the bilinear form Dh (·, ·) as follows: Dh (ph , qh ) = ε(ph − ϱh ph , qh ) = ε(ph − ϱh ph , qh − ϱh qh ).
(38)
Remark 1. Problem IV is also referred to as Euler backward one step fully discrete SMFVE formulation. 3. Error estimates of solutions for fully discrete SMFVE formulation of Problem II In order to derive the existence, the uniqueness, the stability, and the error estimates of the solutions for fully discrete SMFVE formulation, i.e., Problem IV, it is necessary to introduce some preliminary lemmas. ¯ )2 and Wh ⊂ C (Ω ¯ ) as well (17)–(29), from [2,12,14,16,22,24], we have the following two lemmas. Due to Xh ⊂ C (Ω Lemma 3. There hold the following results: ah (uh , Πh∗ vh ) = a(uh , vh ), a1h (vh , uh , Πh∗ wh ) = a1 (vh , uh , wh ), ∗ a1h (vh , uh , Πh uh ) = 0, ∀uh , vh , wh ∈ Xh ,
(39)
dh (Th , ρh∗ φh ) = d(Th , φh ), a2h (uh , Th , ρh∗ φh ) = a2 (uh , Th , φh ), ∗ a2h (uh , Th , ρh Th ) = 0, ∀Th , φh ∈ Wh , ∀uh ∈ Xh ,
(40)
bh (ph , Πh∗ vh ) = −b(ph , vh ),
(41)
∀vh ∈ Xh , ∀ph ∈ Mh .
Further, ah (uh , Πh vh ) and dh (Th , ρh wh ) are all symmetric, bounded, and positive definite, i.e., ∗
∗
ah (uh , Πh∗ vh ) = ah (vh , Πh∗ uh ),
∀uh , vh ∈ Xh ,
(42)
dh (Th , ρh wh ) = dh (wh , ρh Th ),
∀Th , wh ∈ Wh ,
(43)
∗
∗
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
77
and there exist positive constants h0 , C0 , and C˜ 0 such that, when 0 < h 6 h0 , ah (uh , Πh∗ uh ) > ν|uh |21 ,
|ah (uh , Πh∗ vh )| 6 C0 ∥uh ∥1 ∥vh ∥1 ,
dh (Th , ρh Th ) > γ0 | | , ∗
−1
|dh (Th , ρh wh )| 6 C˜ 0 ∥Th ∥1 ∥wh ∥1 , ∗
Th 21
∀uh , vh ∈ Xh , ∀Th , wh ∈ Wh .
(44) (45)
Lemma 4. There holds the following statement:
(uh , Πh∗ vh ) = (vh , Πh∗ uh ),
∀ u h , vh ∈ X h .
(46)
For any u ∈ H m (Ω )2 (m = 0, 1) and vh ∈ Xh ,
|(u, vh ) − (u, Πh∗ vh )| 6 Chm+n ∥u∥m ∥vh ∥n ,
n = 0, 1.
(47)
Set |||uh |||0 = (uh , Πh∗ uh )1/2 , then ||| · |||0 is equivalent to ∥ · ∥0 on Xh , namely there exist two positive constants C3 and C4 such that C3 ∥uh ∥0 6 |||uh |||0 6 C4 ∥uh ∥0 ,
∀uh ∈ Xh .
(48)
Remark 2. For scalar functions, i.e., if uh and vh in Xh are, respectively, substituted with wh and Th in Wh , the results of Lemma 4 hold (see Theorems 3.2.1 and 5.1.5 in [14]). The following discrete Gronwall Lemma (see [16,18]) is useful for the proofs of the existence, the uniqueness, the stability, and the error estimates of the solutions for Problem IV. Lemma 5 (Discrete Gronwall Lemma). If {an }, {bn }, and {cn } are three positive sequences, and {cn } is monotone, they satisfy
¯ an + bn 6 cn + λ
n −1
ai ,
λ¯ > 0, a0 + b0 6 c0 ,
(49)
i =0
then
¯ an + bn 6 cn exp(nλ),
n > 0.
(50)
By Lemma 3 and using the same approach as the proof of Theorem 4.1 in [2], there holds the following inequality sup
(vh ,qh )∈Uh ×Mh
(unh , Πh∗ vh ) + k[a(unh , vh ) − b(pnh , vh ) + b(qh , unh ) + Dh (pnh , qh )] ∥vh ∥1 + ∥qh ∥0
˜ unh ∥0 + k∥∇ unh ∥0 + k∥pnh ∥0 ), > β(∥
∀(unh , pnh ) ∈ Uh × Mh ,
(51)
where β˜ is independent of h and k. For Problem IV, there hold the following results of the existence, the uniqueness, and the stability of its solutions. Theorem 6. Under the hypotheses of Theorem 2, there exists a unique series of solutions (unh , pnh , Thn )(n = 1, 2, . . . , N ) of the fully discrete SMFVE formulation, i.e., Problem IV satisfying
∥unh ∥20 + ∥Thn ∥20 + k2 ∥pnh ∥20 + k
n (∥uih ∥21 + ∥Thi ∥21 ) 6 C˜ ∥ψ∥20 + ∥ψ∥40 ,
(52)
i=1
where C˜ is a constant only depending on C3 , C4 , and ν . Proof. By Lemma 3, the third equation of (30) and its initial condition are rewritten as the following system of equations:
D (Thn , wh ) = G(wh ), Th0 = ρh ψ(x, y),
∀wh ∈ Wh , n = 1, 2, . . . , N , (x, y) ∈ Ω ,
(53)
where D (Thn , wh ) = (Thn , ρh∗ wh ) + kd(Thn , wh ) + ka2 (unh−1 , Thn , wh ) and G(wh ) = (Thn−1 , ρh∗ wh ). For given uhn−1 and Thn−1 , D (·, ·) and G(·) are obviously bounded bilinear and linear functions, respectively. Moreover, by using (4) and Lemmas 3 and 4, we have that D (Thn , Thn ) = (Thn , ρh∗ Thn ) + kd(Thn , Thn ) + ka2 (unh−1 , Thn , Thn ) = (Thn , ρh∗ Thn ) + kd(Thn , Thn ) > α1 ∥Thn ∥21 (where α1 = min{C3 , kγ0−1 }), namely D (·, ·) is the positive definite bilinear function. Therefore, the system of Eq. (53) has a unique series of solutions {Thn }Nn=1 . Taking wh = Thn in (53), by using Lemma 3, Hölder inequality, and Cauchy inequality, we obtain that
|||Thn |||20 + kγ0−1 ∥∇ Thn ∥20 6 D (Thn , Thn ) = G(Thn ) 6
1 2
2
(|||Thn |||20 + |||Thn−1 |||0 ),
n = 1, 2, . . . , N .
(54)
78
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Further, we obtain that 2
|||Thn |||20 + 2kγ0−1 ∥∇ Thn ∥20 6 |||Thn−1 |||0 ,
n = 1, 2, . . . , N .
(55)
Summing (55) from 1 to n and using ∥Th0 ∥0 = ∥ρh ψ(x, y)∥0 6 ∥ψ(x, y)∥0 yields that
∥Thn ∥20 + k
n
∥∇ Thi ∥20 6 C0 ∥ψ∥20 ,
n = 1, 2, . . . , N ,
(56)
i=1
where C0 = C4 min 1/(2γ0 ), C3−1 . For known Thn , the first and second equations as well as their initial condition in Problem IV have a unique series of solutions (unh , pnh )(n = 1, 2, . . . , N ) by means of the theory of mixed FE methods (see [18,5]) and inequality (51). Taking vh = unh in the first equation of Problem IV and qh = pnh in the second equation of Problem IV, by using Lemmas 3 and 4, Hölder inequality, and Cauchy inequality, we obtain that
1 2
2
(|||unh |||20 − |||unh−1 |||0 ) + kν|unh |21 + kε∥pnh − ρh pnh ∥20 = k(Thn j , Πh∗ unh ) 6
k 2ν
∥Thn ∥2−1 +
kν 2
|unh |21 .
(57)
Further, we obtain that 2
|||unh |||20 − |||unh−1 |||0 + kν|unh |21 6 kν −1 ∥Thn ∥2−1 .
(58)
Summing (58) from 1 to n yields that
|||unh |||20 + kν
n
|uih |21 6 kν −1
i =1
n
∥Thi ∥2−1 6 kν −1
i =1
n
∥Thi ∥20 6 ν −1 C0 ∥ψ∥20 .
(59)
i =1
By the first equation of Problem IV and using (51), Hölder inequality, and (59), we obtain that
β˜ k∥pnh ∥0 6
sup
(vh ,qh )∈Uh ×Mh
6 sup
(unh , Πh∗ vh ) + k[a(unh , vh ) − b(pnh , vh ) + b(qh , unh ) + Dh (pnh , qh )] ∥vh ∥1 + ∥qh ∥0
k(Thn j , Πh∗ vh ) − ka1 (unh−1 , unh , Πh∗ vh ) + (unh−1 , Πh∗ vh )
∥v h ∥1
vh ∈Uh
6 k∥
Thn −1
∥
+ kN0 ∥∇
∥ ∥∇ unh ∥0 + ∥uhn−1 ∥−1 .
unh−1 0
(60)
Combining (60) and (59) with (56) yields (52) and completes the proof of Theorem 6. For (u , p ) ∈ X × M and n
n
unh
∈ Xh , put
Ah ((Sh un , Qh pn ); (vh , qh )) = a(Sh un , vh ) − b(Qh pn , vh ) + b(qh , Sh un ) + a1 (uhn−1 , Sh un , vh ),
∀(vh , qh ) ∈ Xh × Mh ,
(61)
A((un , pn ); (vh , qh )) = a(un , vh ) − b(pn , vh ) + b(qh , un ) + a1 (un−1 , un , vh ),
∀(vh , qh ) ∈ Xh × Mh .
By the FE methods (see, e.g., [16,18]) for the non-stationary Navier–Stokes equations, there holds the following Lemma 7. Lemma 7. Let (Sh , Qh ) : X × M → Xh × Mh be the stabilized Navier–Stokes projection, i.e., for (un , pn ) ∈ X × M, there exist (Sh un , Qh pn ) (n = 0, 1, . . . , N) such that kAh ((Sh un , Qh pn ); (vh , qh )) + (Sh un − Sh un−1 , vh ) + kDh (Qh pn , qh ) = kA((un , pn ); (vh , qh )) + (un − un−1 , vh ),
∀(vh , qh ) ∈ Xh × Mh ,
(62)
∥Sh un ∥1 + ∥Qh pn ∥0 6 C (∥un ∥1 + ∥pn ∥0 ).
(63)
If (u , p ) ∈ H (Ω ) × H (Ω ) and h = O(k), there hold the following equalities n
n
2
2
1
∥un − Sh un ∥20 + k2 ∥pn − Qh pn ∥20 + kν
n i =1
∥ui − Sh ui ∥21 6 Ch5
n (∥ui ∥22 + ∥pi ∥21 ),
0 6 n 6 N.
(64)
i =1
And if (un , pn , T n ) ∈ H 2 (Ω )2 × H 1 (Ω ) × H 2 (Ω ) (n = 1, 2, . . . , N) are the solutions to Problem III, then there exists a positive constant C independent of the discretization parameters such that
∥un − Sh un ∥20 + k2 ∥pn − Qh pn ∥20 + kν
n i =1
∥ui − Sh ui ∥21 6 Ch4 ∥ψ∥22 ,
0 6 n 6 N.
(65)
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
79
By the classical FE methods (see, e.g., [10,18]) for elliptic equations, we have the following Lemma 8. Lemma 8. Let Rh : W → Wh be a generalized Ritz projection, i.e., for given uhn−1 ∈ Xh , T n−1 ∈ W , and Thn−1 ∈ Wh , and for any T n ∈ W (n = 1, 2, . . . , N), there exist Rh T n ∈ Wh (n = 1, 2, . . . , N) such that
(Rh T n , wh ) + kd(Rh T n , wh ) + ka2 (unh−1 , Rh T n , wh ) − (Rh T n−1 , wh ) = (T n , wh ) + kd(T n , wh ) + ka2 (un−1 , T n , wh ) − (T n−1 , wh ),
∀wh ∈ Wh , n = 1, 2, . . . , N .
(66)
If (un , pn , T n ) (n = 1, 2, . . . , N) are the solutions of Problem III and T n ∈ H 2 (Ω ) ∩ W , there hold the following inequalities
∥Rh T n ∥0 + k1/2 ∥Rh T n ∥1 6 C ∥T n ∥1 , ∥Rh T n − T n ∥20 + kγ0−1
n
n = 0, 1, 2, . . . , N ,
|Rh T i − T i |21 6 Ch4 ∥ψ∥22 ,
(67)
n = 0, 1, 2, . . . , N .
(68)
i =1
Theorem 9. Let (u, p, T ) be the solution to Problem II and (unh , pnh , Thn ) the solution to Problem IV. Under the hypotheses of Theorems 1 and 2, if p0h = p0 = 0 (or p0h = Qh p0 ), h = O(k), u0 ∈ H 2 (Ω )2 , and ψ ∈ H 1 (Ω ), then there hold that
∥u(tn ) − unh ∥0 + ∥T (tn ) − Thn ∥0 + k[∥p(tn ) − pnh ∥0 + ∥u(tn ) − unh ∥1 + ∥T (tn ) − Thn ∥1 ] 6 Ch2 (∥ψ∥1 ) , n = 1, 2, . . . , N .
(69)
Proof. Subtracting Problem IV from Problem III taking v = vh and q = qh , using Lemmas 3 and 4, we obtain the system of error equations:
n (u − unh , vh ) + (unh − Πh∗ unh , vh − Πh∗ vh ) + ka(un − unh , vh ) + ka1 (un−1 , un , vh ) − ka1 (uhn−1 , unh , vh ) −kb(pn − pnh , vh ) = k(j (T n − Thn ), vh ) − k(jThn , Πh∗ vh − vh ) + (un−1 − uhn−1 , vh ) − (unh−1 , Πh∗ vh − vh ), ∀vh ∈ Xh , n = 1, 2, . . . , N , b(qh , un − unh ) − ε(pnh − ϱh pnh , qh − ϱh qh ) = 0, ∀qh ∈ Mh , n = 1, 2, . . . , N , (T n − Thn , wh ) − (Thn , ρh∗ wh − wh ) + kd(T n − Thn , wh ) + ka2 (un−1 , T n , wh ) − ka2 (uhn−1 , Thn , wh ) = (T n−1 − Thn−1 , wh ) − (Thn−1 , ρh∗ wh − wh ), ∀wh ∈ Wh , n = 1, 2, . . . , N , 0 u − u0h = 0, T 0 − Th0 = ψ(x, y) − ρh ψ(x, y), (x, y) ∈ Ω .
(70)
Let En = Sh un − unh . By using (62) and the system of error equations (70), we obtain that
∥En ∥20 + kν|En |21 = (En , En ) + ka(En , En ) = [(Sh un − un , En ) + ka(Sh un − un , En )] + [(un − unh , En ) + ka(un − unh , En )] = [(Sh un−1 − un−1 , En ) + kb(Qh pn − pn , Eh ) + ka1 (un−1 , un , En ) − ka1 (uhn−1 , Sh un , En )] + [kb(pn − pnh , En ) − ka1 (un−1 , un , En ) + ka1 (unh−1 , unh , En ) + (unh , Πh∗ En − En ) + k(j (T n − Thn ), En ) − k(jThn , Πh∗ En − En ) + (un−1 − uhn−1 , En ) − (uhn−1 , Πh∗ En − En )] = (En−1 , En ) + k(j (T n − Thn ), En ) − k(jThn , Πh∗ En − En ) + (unh − unh−1 , Πh∗ En − En ) + kb(Qh pn − pnh , Eh ).
(71)
By using Hölder inequality, Cauchy inequality, (47), (65), Lemmas 3 and 4, and Taylor’s formula, if h = O(k) we get that
|(En−1 , En )| 6 ∥En−1 ∥0 ∥En ∥0 6
1 2
1
∥En−1 ∥20 + ∥En ∥20 ,
(72)
2
|k(j (T n − Thn ), En ) − k(jThn , Πh∗ En − En )| 6 Ck∥T n − Thn ∥0 ∥En ∥0 + Ckh2 ∥Thn ∥0 ∥∇ En ∥0 νk 6 Ck∥T n − Thn ∥20 + Ck∥En ∥20 + Ckh4 + |En |21 , 4
|(unh − uhn−1 , Πh∗ En − En )| = (Sh un − un + un − un−1 + un−1 − Sh un−1 + En−1 − En , Πh∗ En − En ) 6 Ch(∥Sh un − un ∥20 + h∥un − un−1 ∥21 + ∥un−1 νk − Sh un−1 ∥20 + ∥En−1 ∥20 + ∥En ∥20 ) + |En |21 4 νk 6 Ch(∥En ∥20 + ∥En−1 ∥20 ) + Ch5 ∥ψ∥22 + Ck2 h3 ∥ut ∥2L∞ (H 1 ) + |En |21 . 4
(73)
(74)
80
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Note that b(qh , Sh un − un ) = −kε(Qh pn − ϱh (Qh pn ), qh − ϱh (qh )). By the property of operator ϱh and the second equation of (70), we have that b(Qh pn − pnh , Sh un − unh ) = b(Qh pn − pnh , Sh un − un ) + b(Qh pn − pnh , un − unh )
ε = − (Qh pn − pnh − ϱh (Qh pn − pnh ), Qh pn − pnh − ϱh (Qh pn − pnh )) 2
ε − (Qh pn−1 − pnh−1 − ϱh (Qh pn−1 − pnh−1 ), Qh pn − pnh − ϱh (Qh pn − pnh )) 2
ε
ε
6 − ∥Qh pn − pnh − ϱh (Qh pn − pnh )∥20 + ∥Qh pn−1 − ph 4 4
n −1
− ϱh (Qh pn−1 − phn−1 )∥20 .
(75)
Combining (72)–(75) with (71) yields that
ε ε ∥En ∥20 + kν|En |21 + ∥Qh pn − pnh − ϱh (Qh pn − pnh )∥20 − ∥Qh pn−1 − pnh−1 − ϱh (Qh pn−1 − pnh−1 )∥20 2
6∥
En−1 20
2
4
n
∥ + Ckh + Ck∥T −
∥ + Ck(∥ ∥ + ∥En−1 ∥20 ).
Thn 20
En 20
(76)
Summing (76) from 1 to n, if h is sufficiently small such that Ch 6 1/2 in (76) and p0h = p0 = 0 (or p0h = Qh p0 ), we obtain that
∥En ∥20 + kν
n
|Ei |21 6 Cnkh4 + Ck
i =1
n
∥T i − Thi ∥20 + Ck
i=1
n
∥Ei ∥20 .
(77)
i =0
If k is sufficiently small such that Ck 6 1/2 in (77), we obtain from (77) that
∥En ∥20 + 2kν
n
|Ei |21 6 Ch4 + Ck
i=1
n
∥T i − Thi ∥20 + Ck
n −1
i=1
∥Ei ∥20 .
(78)
i =0
Applying Gronwall Lemma to (78) yields that
∥ ∥ + kν En 20
n
En 21
4
| | 6C h +k
i =1
n
i
∥T −
Thi 20
∥
exp(Ckn).
(79)
i=1
By using triangle inequality, (64), and Lemma 7, we get that n
∥u −
∥ + kν
unh 20
n
i
|u −
uih 21
4
| 6C h +k
i=1
n
i
∥T −
Thi 20
∥
.
(80)
i =1
By using (51) and the system of error equations (70), we have that
β˜ k∥Qh pn − pnh ∥0 6
sup
(vh ,qh )∈Uh ×Mh
(Sh un − unh , Πh∗ vh )
+ k a(Sh un − unh , vh ) − b(Qh pn − pnh , vh ) + b(qh , Sh un − unh ) + Dh (pnh , qh ) /(∥vh ∥1 + ∥qh ∥0 ) |b(Qh pn − pn , vh )| |b(pn − pnh , vh )| + k sup ∥v h ∥1 ∥v h ∥1 vh ∈Uh vh ∈Uh −1 n n n n n n ∗ 6 C (k∥Qh p − p ∥0 + ∥u − uh ∥0 ) + k k (u − uh , Πh vh ) − (un−1 − unh−1 , Πh∗ vh ) 6 C ∥Sh un − unh ∥0 + k sup
+ a(un − unh , vh ) + a1 (un−1 , un − unh , vh ) + a1 (un−1 − unh−1 , unh , vh ) − (j (T n − Thn ), Πh∗ vh ) /∥vh ∥1 n −1 6 Ck ∥Qh pn − pn ∥0 + k−1 ∥un − unh ∥0 + ∥un−1 − uh ∥0 + |un−1 − unh−1 |1 + |un − unh |1 + ∥T n − Thn ∥0 .
(81)
Combining (80) with (81) yields that
∥Qh pn − pnh ∥0 6 Ch ∥ψ∥2 + ∥T n − Thn ∥0 + k
n i =1
1/2 ∥T i − Thi ∥20
.
(82)
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
81
By using triangle inequality, (65) and (82), we get that
∥pn − pnh ∥0 6 Ch ∥ψ∥2 + ∥T n − Thn ∥0 + k
n
1/2 . ∥T i − Thi ∥20
(83)
i =1
Let en = Rh T n − Thn . By using the system of error equations (70), (4), (66)–(68) in Lemma 8, Lemma 3, and (16), we obtain that
∥en ∥20 + kγ0−1 |en |21 = (en , en ) + kd(en , en ) = [(Rh T n − T n , en ) + kd(Rh T n − T n , en )] + [(T n − Thn , ρh∗ en ) + kd(T n − Thn , en )] = [(Rh T n−1 − T n−1 , en ) + ka2 (un−1 , T n , en ) − ka2 (unh−1 , Rh T n , en )] + [(Thn , ρh∗ en − en ) + ka2 (unh−1 , Thn , en ) − ka2 (un−1 , T n , en ) + (T n−1 − Thn−1 , en ) − (Thn−1 , ρh∗ en − en )] = (en−1 , en ) + (Thn − Thn−1 , ρh∗ en − en ) 6 ∥en−1 ∥0 ∥en ∥0 + Ch(∥en ∥0 + ∥Rh T n − T n ∥0 + h∥T n − T n−1 ∥1
+ ∥T n−1 − Rh T n−1 ∥0 + ∥en−1 ∥0 )|en |1 6
1 2
k ∥en−1 ∥20 + ∥en ∥20 + Ch h4 + k2 h2 + ∥en ∥20 + ∥en−1 ∥20 + |en |21 . 2γ0
(84)
Further, we get from (84) that
∥en ∥20 + kγ0−1 |en |21 6 Ch h4 + k2 h2 + ∥en ∥20 + ∥en−1 ∥20 + ∥en−1 ∥20 .
(85)
Summing (85) from 1 to n and using (68) and (12) yield that
∥en ∥20 + kγ0−1
n
|ei |21 6 Cnh(h4 + k2 ) + ∥e0 ∥20 + Ck
i=1
n
∥ei ∥20
i =1
6 C (h4 + k2 h2 ) + Ck
n
∥ei ∥20 + C ∥Rh ψ − ψ∥20 + C ∥ψ − ρh ψ∥20
i =1
6 C (h4 + k2 h2 ) + Ck
n
∥ei ∥20 .
(86)
i =1
If k is sufficiently small such that Ck 6 1/2 in (86), we obtain from (86) that
∥en ∥20 + kγ0−1
n
|ei |21 6 Ch4 + Ck
i=1
n −1
∥ei ∥20 .
(87)
i =0
Applying the Gronwall Lemma to (87) yields that
∥en ∥20 + kγ0−1
n
|ei |21 6 Ch4 exp(Ck) 6 Ch4 .
(88)
i=1
By using triangle inequality, (88), and (68), we obtain that
∥T n − Thn ∥20 + kγ0−1
n
|T i − Thi |21 6 Ch4 .
(89)
i =1
Combining (79) with (89) yields that
∥un − unh ∥20 + kν
n
|ui − uih |21 6 Ch4 .
(90)
i =1
Combining (89) with (83) yields that
∥pn − pnh ∥0 6 Ch. Combining (90) and (91) with (89) yields (69) which completes the proof of Theorem 9.
(91)
82
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Combining Theorem 2 with Theorem 9 yields the following conclusion. Theorem 10. Under the hypotheses of Theorems 2 and 9, there hold the error estimates between the solution (u, p, T ) to Problem II and the solutions (unh , pnh , Thn )(n = 1, 2, . . . , N ) to Problem IV as follows
∥u(tn ) − unh ∥20 + ∥T (tn ) − Thn ∥20 + k
n (|u(ti ) − uih |21 + |T (ti ) − Thi |21 ) 6 C (k2 + h4 ),
(92)
i =1
∥pn − pnh ∥0 6 C (k + h),
(93)
where C is a constant which is only dependent on γ0 , ν , Ω , ψ , and maximum total time upper bound tN , but independent of k and h. 4. Some numerical experiments In this section, we present some numerical experiments with a physical model of cavity flow and Reynolds number Re = 7.1 × 103 and Pr = 0.71 by the fully discrete SMFVE Formulation, i.e., Problem IV to validate that the numerical computational results are consistent with the theoretical conclusions. Moreover, it is shown that the fully discrete SMFVE formulation based on two local Gauss integrals and parameter-free is feasible and efficient for finding the numerical solutions of the non-stationary conduction–convection problem and its numerical solutions are more stable than those of the classical fully discrete FVE formulation without any stabilization in [16]. Especially, its degrees of freedom are far fewer than those of the classical fully discrete FVE formulation without any stabilization in [16]. 4.1. Numerical experiment for the fully discrete SMFVE formulation
¯ = [0, 1] × [0, 1]. We first divide the cavity into 100 × 100 = 10 000 small Let the side length of the cavity be 1 and Ω squares with side length △x = △y = 0.01, and then link diagonal of the√square to divide each square into two triangles in the same direction, which constitutes triangularizations ℑh with h = 2 × 10−2 . The dual decomposition ℑ∗h is taken as barycenter dual decomposition, namely the barycenter of the right triangle K ∈ ℑh is taken as the node of the dual decomposition. We take the time step increment as k = 0.01. Let the initial value and the boundary value of u = (u1 , u2 ) be equal to 0. And let T = 0 on left and lower boundary of the cavity, ∂ T /∂ y = 0 on upper boundary of the cavity, and T = 4y(1 − y) on the right boundary of the cavity. The initial value of T is 0. We find a SMFVE numerical solution (unh , pnh , Thn ) by means of Problem IV when n = 200 (i.e., t = 2), whose unh , pnh , and Thn are depicted graphically on the left charts in Figs. 2–4, respectively. The degrees of freedom for the fully discrete SMFVE Formulation, i.e., Problem IV on each time level are about 4Np = 4 × 104 (where Np is the number of nodal points in ℑh , see [5,10,18]). 4.2. Numerical experiment for the classical FVE formulation without any stabilization
¯ = [0, 1] × [0, 1]. We also divide the cavity into 100 × 100 = 10000 small Let still the side length of the cavity be 1 and Ω squares with side length △x = △y = 0.01 and link diagonal √ of the square to divide each square into two triangles in the same direction and to constitute triangularizations ℑh (h = 2 × 10−2 ). In order to make BB condition satisfied,√ it is necessary ˜ h with h = 2 × 10−2 /2. to link three midpoints on sides of each triangle in ℑh , which constitutes triangularizations ℑ ˜ ∗h is also taken as barycenter dual decomposition, namely the barycenter of the triangle K ∈ ℑ˜ h The dual decomposition ℑ is taken as the node of the dual decomposition. The time step increment, the initial values and the boundary values are the same as those in Section 4.1. The classical fully discrete FVE formulation without any stabilization is written as follows (see [16]). Problem V. Find (unh , pnh , Thn ) ∈ X h × M h × W h (1 6 n 6 N ) such that
(∂¯ un , Π ∗ v ) + ah (unh , Πh∗ vh ) + a1h (unh−1 , unh , Πh∗ vh ) + bh (pnh , Πh∗ vh ) = (jThn , Πh∗ vh ), t h nh h b(qh , uh ) = 0, (∂¯t Thn , ρh∗ wh ) + dh (Thn , ρh∗ wh ) + a2h (unh−1 , Thn , ρh∗ wh ) = 0, 0 uh = 0, Th0 = ρh ψ(x, y),
∀vh ∈ X h , ∀qh ∈ M h , ∀wh ∈ W h , (x, y) ∈ Ω ,
where all other signs are the same as those in Problem IV except for trial function spaces X h , W h , and M h of velocity, temperature, and pressure being, respectively, defined as follows in order to make BB condition satisfied (see [16]):
˜h , X h = vh ∈ X ∩ C (Ω )2 ; vh |K ∈ P12 (K ), ∀K ∈ ℑ
˜h , W h = wh ∈ W ∩ C (Ω ); wh |K ∈ P1 (K ), ∀K ∈ ℑ M h = qh ∈ L20 (Ω ); qh |K ∈ P0 (K ), ∀K ∈ ℑh .
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
83
Fig. 2. When Re = 7.1 × 103 and Pr = 0.71, the left chart and the right chart are the SMFVE solution and the classical FVE solution of the velocity u at the time level t = 2, respectively.
Fig. 3. When Re = 7.1 × 103 and Pr = 0.71, the left chart and the right chart are the SMFVE solution and the classical FVE solution of the temperature T at the time level t = 2, respectively.
Thus, the degrees of freedom for the classical fully discrete FVE Formulation, i.e., Problem V on each time level are about 3Np + NK + 3Ns ≈ 14Np = 14 × 104 (where Np , NK , and Ns are the numbers of nodal points, elements, and sides in ℑh , respectively) since Np : NK : Ns ≈ 1 : 2 : 3 (see [5,10,18]). They are 3.5 times those for the fully discrete SMFVE Formulation. We also find a FVE numerical solution (unh , pnh , Thn ) to Problem V when n = 200 (i.e., t = 2), whose unh , pnh , and Thn are depicted graphically on the right charts in Figs. 2–4, respectively. The curves of the upper and lower charts in Fig. 5 are the relative errors (log 10) of the fully discrete SMFVE solutions and the classical fully discrete FVE solutions at time t ∈ (0, 2] with respect to L2 -norm, respectively. Due to the truncation error accumulation in computational process, the errors of numerical solutions appear increase (see Fig. 5), but the truncation error accumulation for the fully discrete SMFVE formulation Problem IV is far smaller than that for the classical fully discrete FVE formulation with any stabilization and the relative errors (which does not exceed 2×10−2 ) of SMFVE numerical solutions are far smaller than those of classical FVE solutions (also see Fig. 5). Moreover, it is shown that the numerical computational results are consistent with the theoretical conclusions for the fully discrete SMFVE formulation, since the theoretical and numerical errors achieve O(10−2 ). Since the fully discrete SMFVE formulation adopts the stabilization based on two local Gauss integrals and parameter-free, its numerical solutions are more stable than those of the fully discrete FVE formulation without any stabilization (also see Figs. 2–4). For example, the chart of velocity u of the classical FVE solution at t = 2 has obvious chaos (see the right chart in Fig. 2) so that its numerical computational effectiveness has no far better than that of SMFVE formulation (see the left chart in Fig. 2). Especially, the fully discrete SMFVE formulation can reduce great many
84
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Fig. 4. When Re = 7.1 × 103 and Pr = 0.71, the left chart and the right chart are the SMFVE solution and the classical FVE solution of the pressure p at the time level t = 2, respectively.
Fig. 5. When Re = 7.1 × 103 and Pr = 0.71, the upper chart and the lower chart are the relative errors (log 10) of the SMFVE solutions and the classical FVE solution of the velocity u, the temperature T , and the pressure p at the time t ∈ (0, 2], respectively.
of degrees of freedom (its degrees of freedom are by a factor of 3.5 times fewer than those for the classical fully discrete FVE formulation in [16]). Therefore, the fully discrete SMFVE formulation is far better than the classical fully discrete FVE formulation without any stabilization in [16].
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
85
5. Conclusions and discussions We have established directly a fully discrete SMFVE formulation from the semi-discrete formulation with respect to time for the non-stationary conduction–convection problem. Thus, we could avoid the semi-discrete SMFVE formulation about spatial variables. We have also provided the existence and stabilization of the fully discrete SMFVE solutions and the error estimates between the fully discrete SMFVE solutions and the accurate solution. We have also given some numerical experiments illustrating that the numerical computational results are consistent with the theoretical results obtained previously, validating that the fully discrete SMFVE formulation is feasible and efficient for finding the numerical solutions of the non-stationary conduction–convection problem, and showing that the errors of the fully discrete SMFVE solutions are far smaller than those of the classical fully discrete FVE solutions, its degrees of freedom are far fewer than those of the classical fully discrete FVE formulation without any stabilization in [16], and the numerical solutions for the fully discrete SMFVE formulation are more stable than those for the classical fully discrete FVE formulation without any stabilization in [16]. Especially, the trial function spaces of velocity and pressure for the fully discrete SMFVE formulation satisfy the BB condition naturally, which is a main improvement for method in [16]. Consequently, the fully discrete SMFVE formulation are far better than the classical fully discrete FVE formulation without any stabilization in [16]. Though an FVE method for time-dependent Stokes equations (see [22]) and a penalty FVE method for transient Navier–Stokes equations (see [12]) have been provided, the non-stationary conduction–convection problem as indicated in [16] is different from either of the above. The Stokes equations are linear, while our conduction–convection problem is nonlinear. Moreover, the non-stationary conduction–convection problem includes a nonlinear energy equation of the temperature field with the velocity vector field coupled, which is more complex than the transient Navier–Stokes equations without energy equation in [12]. Therefore, the study of the SMFVE method for the non-stationary conduction–convection problem has far more difficulties, more important, more serviceable, and more challenging than those for the timedependent Stokes equations and the transient Navier–Stokes equations. It is the improvement and innovation for the existing methods (see, e.g., [12,16,22] or others). Acknowledgments The authors cordially thank all the Reviewers, Editors, and Professor & Editor in Chief Goong Chen for their very valuable suggestions that led to improvements of our manuscript. References [1] R.A. Adams, Sobolev Spaces, Academic Press, New York, 1975. [2] J. An, P. Sun, Z.D. Luo, X.M Huang, A stabilized fully discrete finite volume element formulation for non-stationary Stokes equations, Math. Numer. Sin. 33 (2) (2011) 213–224. [3] R.E. Bank, D.J. Rose, Some error estimates for the box methods, SIAM J. Numer. Anal. 24 (4) (1987) 777–787. [4] P. Blanc, R. Eymerd, R. Herbin, A error estimate for finite volume methods for the Stokes equations on equilateral triangular meshes, Numer. Methods Partial Differential Equations 20 (2004) 907–918. [5] F. Brezzi, M. Fortin, Mixed and Hybrid Finite Element Methods, Springer–Verlag, New York, 1991. [6] Z. Cai, S. McCormick, On the accuracy of the finite volume element method for diffusion equations on composite grid, SIAM J. Numer. Anal. 27 (3) (1990) 636–655. [7] P. Chatzipantelidis, R.D. Lazarrov, V. Thome, Error estimates for a finite volume element method for parabolic equations in convex in polygonal domains, Numer. Methods Partial Differential Equations 20 (2004) 650–674. [8] S.H. Chou, D.Y. Kwak, A covolume method based on rotated bilinears for the generalized Stokes problem, SIAM J. Numer. Anal. 35 (1998) 494–507. [9] M.A. Christon, P.M. Gresho, S.B. Sutton, Computational predictability of time-dependent natural convection flows in enclosures (including a benchmark solution), in: MIT Special Issue on Thermal Convection, Internat. J. Numer. Methods Fluids 40 (8) (2002) 953–980 (special issue). [10] P.G. Ciarlet, The Finite Element Method for Elliptic Problems, North-Holland, Amsterdam, 1978. [11] E. DiBenedetto, A. Friedman, Conduction–convection problems with change of phase, J. Differential Equations 62 (1986) 129–185. [12] G.L. He, Y.N. He, Z.X. Chen, A penalty finite volume method for the transient Navier–Stokes equations, Appl. Numer. Math. 58 (2008) 1583–1613. [13] W.P. Jones, K.P. Menziest, Analysis of the cell-centred finite volume method for the diffusion equation, J. Comput. Phys. 165 (2000) 45–68. [14] R.H. Li, Z.Y. Chen, W. Wu, Generalized Difference Methods for Differential Equations-Numerical Analysis of Finite Volume Methods, in: Monographs and Textbooks in Pure and Applied Mathematics, vol. 226, Marcel Dekker Inc., New York, 2000. [15] Y. Li, R.H. Li, Generalized difference methods on arbitrary quadrilateral networks, J. Comput. Math. 17 (1999) 653–672. [16] H. Li, Z.D. Luo, P. Sun, J. An, A finite volume element formulation and error analysis for the non-stationary conduction–convection problem, J. Math. Anal. Appl. 396 (2012) 864–879. [17] Z.D. Luo, The mixed finite element method for the non-stationary conduction–convection problems, Math. Numer. Sin. 20 (1) (1998) 69–88. [18] Z.D. Luo, The bases and Applications of Mixed Finite Element Methods, Chinese Science Press, Beijing, 2006. [19] Z.D. Luo, J. Chen, I.M. Navon, J. Zhu, An optimizing reduced PLSMFE formulation for non-stationary conduction–convection problems, Internat. J. Numer. Methods. Fluids 60 (2009) 409–436. [20] Z.D. Luo, Z.H. X, J. Chen, A reduced MFE formulation based on POD for the non-stationary conduction–convection problems, Acta Math. sci. 31 (5) (2011) 765–1785. [21] E. Süli, Convergence of finite volume schemes for Poissons equation on nonuniform meshes, SIAM J. Numer. Anal. 28 (5) (1991) 1419–1430. [22] M. Yang, H.L. Song, A postprocessing finite volume method for time-dependent Stokes equations, Appl. Numer. Math. 59 (2009) 1922–1932. [23] X. Ye, On the relation between finite volume and finite element methods applied to the Stokes equations, Numer. Methods Partial Differential Equations 17 (2001) 440–453. [24] M. Zlámal, Finite element methods for nonlinear parabolic equations, RAIRO Anal. Numer. 11 (1) (1977) 93–107.
Contents lists available at SciVerse ScienceDirect
Journal of Mathematical Analysis and Applications journal homepage: www.elsevier.com/locate/jmaa
A fully discrete stabilized mixed finite volume element formulation for the non-stationary conduction–convection problem✩ Zhendong Luo a,∗ , Hong Li b,∗ , Ping Sun c a
School of Mathematics and Physics, North China Electric Power University, Beijing 102206, China
b
School of Mathematical Sciences, Inner Mongolia University, Huhhot 010021, China
c
School of Mathematics and Computer Science, Guizhou Normal University, Guiyang 550001, China
article
info
Article history: Received 3 January 2013 Available online 13 March 2013 Submitted by G. Chen Keywords: Non-stationary conduction–convection problem Stabilized mixed finite volume element formulation Local Gauss integrals and parameter-free Error estimate
abstract At first, a semi-discrete formulation with respect to time for the non-stationary conduction–convection problem is recalled. Then, a fully discrete stabilized mixed finite volume element (SMFVE) formulation based on two local Gauss integrals and parameterfree is established directly from the semi-discrete formulation with respect to time for the non-stationary conduction–convection problem. Following this, the error estimates for the fully discrete SMFVE solutions are derived by means of the standard mixed finite element method. Finally, some numerical experiments are presented illustrating that the numerical errors are consistent with theoretical results, the degrees of freedom of the fully discrete SMFVE formulation are far fewer than those of the classical finite volume element (FVE) formulation without any stabilization, and its numerical solutions are more stable than those of the classical FVE formulation without any stabilization, thus validating that the fully discrete SMFVE formulation is feasible and efficient for finding the numerical solutions of the non-stationary conduction–convection problem. © 2013 Elsevier Inc. All rights reserved.
1. Introduction Let Ω ⊂ R2 be a bounded and connected polygonal domain. Consider the non-stationary conduction–convection problem including the velocity vector field and the pressure and temperature fields (see [11,17,18]). Problem I. Find u = (u1 , u2 ), p, and T such that, for tN > 0,
ut − ν ∆u + (u · ∇)u + ∇ p = T j , ∇ · u = 0 , Tt − γ0−1 ∆T + (u · ∇)T = 0, u (x, y, t ) = 0, T (x, y, t ) = ϕ(x, y, t ), u(x, y, 0) = 0, T (x, y, 0) = ψ(x, y)
(x, y, t ) ∈ Ω × (0, tN ), (x, y, t ) ∈ Ω × (0, tN ), (x, y, t ) ∈ Ω × (0, tN ), (x, y, t ) ∈ ∂ Ω × (0, tN ), (x, y) ∈ Ω ,
(1)
✩ This work is jointly supported by National Science Foundation of China (11271127, 11061009, and 11061021), Natural Science Foundation of Inner Mongolia (2012MS0106), Science Research Program of Guizhou (GJ[2011]2367), Science Research Program of Inner Mongolia Advanced Education (NJ10006), and Special Funds for Co-construction Project of Beijing and North China Electric Power University. ∗ Corresponding authors. E-mail addresses: [email protected] (Z. Luo), [email protected] (H. Li).
0022-247X/$ – see front matter © 2013 Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.jmaa.2013.03.001
72
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
where u = (u1 , u2 ) is the fluid velocity vector, p the pressure, T the temperature, tN the total time, j = (0, 1) the unit √ √ vector, ν = Pr /Re, Re the Reynolds number, Pr the Prandtl number, γ0 = RePr, and ϕ(x, y, t ) and ψ(x, y) are two given functions. For the sake of convenience and without loss of generality, we may as well suppose that ϕ(x, y, t ) = 0 in the following theoretical analysis. Since Problem I constitutes a system of nonlinear PDEs including the velocity vector field and the pressure and the temperature fields, there are no analytical solutions in general. One has to rely on numerical solutions (see, e.g., [9,17,18]). However, most of the existing papers use either the finite element (FE) method or finite difference (FD) schemes as discretization tools. The finite volume element (FVE) method (see [6,13,21]) is considered as one of the most effective numerical methods due to the following three advantages. First, it preserves the integral invariants of conservations of total energy and mass. Second, it is suitable for computations with complicated boundary and has higher accuracy than that of the FD method. Third, it has the same accuracy as FE method but is much more convenient and simpler than the FE method. It is also known as a box method (see [3]) or a generalized difference method (see [14,15]). It has been widely used to find numerical solutions of various types of PDEs, for example, second order elliptic equations, parabolic equations, Stokes equations, and the Navier–Stokes equations (see, e.g., [2–4,6–8,12–15,21–23]). Though a fully discrete FVE formulation without any stabilization for the non-stationary conduction–convection problem is proposed in [16], its trial function spaces of velocity and pressure have to satisfy an important convergence stability condition, i.e., the Babuška–Brezzi (BB) inequality (see [16,18]). Thus, there exist many difficulties for its theoretical analysis and actual applications. Therefore, it is very necessary to improve the method in [16]. To this end, we devote to establishing a fully discrete SMFVE formulation based on two local Gauss integrals and parameter-free for the nonstationary conduction–convection problem in this paper such that its trial function spaces of velocity and pressure satisfy the BB condition naturally. In addition, an optimizing reduced Petrov–Galerkin least squares mixed FE formulation based on residuals (see [19]) and a reduced mixed FE formulation without any stabilization (see [20]) have been established for the non-stationary conduction–convection problem, but they are completely different from the fully discrete SMFVE formulation in this paper since the FVE method is entirely different from and far more advantageous than the FE method as is mentioned above. The paper is organized as follows. In Section 2, we first recall the semi-discrete formulation with respect to time for the non-stationary conduction–convection problem, and then establish directly a fully discrete SMFVE formulation from the semi-discrete formulation with respect to time which could avoid the semi-discrete SMFVE formulation about spatial variables. It is a new type of study attempt for the non-stationary conduction–convection problem. In Section 3, the error estimates between the fully discrete SMFVE solutions and the exact solution are provided. In Section 4, some numerical experiments are presented illustrating that the errors of the fully discrete SMFVE solutions are consistent with the theoretical results obtained previously, the degrees of freedom of the fully discrete SMFVE formulation are far fewer than those of the classical fully discrete FVE formulation without any stabilization (see [16]), and the numerical solutions for the fully discrete SMFVE formulation are more stable than those for the classical fully discrete FVE formulation without any stabilization in [16], thus validating that the SMFVE formulation is feasible and efficient for finding the numerical solutions of the nonstationary conduction–convection problem. Section 5 provides main conclusions and some discussions. 2. Semi-discrete formulation about time and fully discrete SMFVE formulation for Problem I 2.1. Recall semi-discrete formulation about time for Problem I Sobolev spaces along with their properties used in this context are standard (see [1]). For example, H m (Ω ) = The 2 α α1 α2 v; 06|α|6m Ω |∂ v/(∂ x ∂ y )|0 dxdy < ∞ (where m > 0 is integer, α = (α1 , α2 ), α1 and α2 are two non-negative 1/2 α α1 α2 2 integers, and |α| = α1 + α2 ) with norm and semi-norm denoted by ∥v∥m = and 06|α|6m Ω |∂ v/(∂ x ∂ y )|0 dxdy 1/2 t α α1 α2 2 , respectively; L2 (L2 ) = v(x, y, t ); 0 N Ω v 2 (x, y, t )dxdydt < ∞ with |v|m = |α|=m Ω |∂ v/(∂ x ∂ y )|0 dxdy 1/2 t norm ∥v∥L2 (L2 ) = 0 N Ω v 2 dxdydt ; H01 (Ω ) = v ∈ H 1 (Ω ); v|∂ Ω = 0 . Let X = H01 (Ω )2 , M = L20 (Ω ) = q ∈ L2 (Ω ); Ω qdxdy = 0 , and W = H01 (Ω ). Put a(u, v ) = ν ∇ u · ∇ vdxdy, ∀u, v ∈ X , b(q, v ) = q divvdxdy, ∀v ∈ X , q ∈ M , Ω Ω 1 [(u∇ v ) · w − (u∇ w ) · v] dxdy, ∀u, v , w ∈ X , a1 (u, v , w ) = 2 Ω 1 [(u · ∇ T )φ − (u · ∇φ)T ] dxdy, ∀u ∈ X , ∀T , φ ∈ W , a2 (u, T , φ) = 2
Ω
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
d(T , φ) = γ0−1
Ω
∇ T · ∇φ dxdy,
73
∀T , φ ∈ W .
Thus, a mixed variational formulation for Problem I is written as follows. Problem II. Find (u, p, T ) ∈ H 1 (0, tN ; X )2 × L2 (0, tN ; M ) × H 1 (0, tN ; W ) such that, for almost all t ∈ (0, tN ),
(u , v ) + a(u, v ) + a1 (u, u, v ) − b(p, v ) = (T j , v ), t b(q, u) = 0, (Tt , φ) + d(T , φ) + a2 (u, T , φ) = 0, u(x, y, 0) = 0, T (x, y, 0) = ψ(x, y),
∀v ∈ X , ∀q ∈ M , ∀φ ∈ W , ( x, y ) ∈ Ω ,
(2)
where (·, ·) denotes inner product in L2 (Ω )2 or L2 (Ω ). There hold the following properties for above trilinear forms a1 (·, ·, ·) and a2 (·, ·, ·) (see [16,18]): a1 (u, v , w ) = −a1 (u, w , v ),
a1 (u, v , v ) = 0,
(3)
a2 (u, T , φ) = −a2 (u, φ, T ),
a2 (u, φ, φ) = 0,
∀u, v , w ∈ X , ∀u ∈ X , ∀T , φ ∈ W . There are the following properties for above bilinear forms a(·, ·), d(·, ·), and b(·, ·) (also see [16,18]):
(4)
a(v , v ) > ν|v |21 ,
sup
b(q, v )
v ∈X
2 1
∀u, v ∈ X ,
|d(T , φ)| 6 γ0 |T |1 |φ|1 , −1
∀φ ∈ W ;
> β∥q∥0 ,
|v |1
|a(u, v )| 6 ν|u|1 |v |1 ,
∀v ∈ X ;
d(φ, φ) > γ0 |φ| , −1
(5)
∀T , φ ∈ W ,
(6)
∀q ∈ M ,
(7)
where β is a constant. Define N0 =
sup
u,v ,w ∈X
a1 (u, v , w ) , |u|1 · |v |1 · |w |1
N˜ 0 =
sup
u∈X ,(T ,φ)∈W ×W
a2 (u, T , φ) . |u|1 · |T |1 · |φ|1
(8)
The following result is classical (see [16–18]). Theorem 1. If ψ ∈ L2 (Ω ), then the Problem II has at least a solution which, in addition, is unique provided that ∥ψ∥20 6
2ν 2 tN /(2N0 tN−1 exp(tN ) + νγ0 N˜ 02 ), and there are the following prior estimates
∥u∥20 + ν∥∇ u∥2L2 (L2 ) 6 tN2 ∥ψ∥20 exp(tN ),
∥T ∥20 + γ0−1 ∥∇ T ∥2L2 (L2 ) 6 ∥ψ∥20 .
Let N be the positive integer, k = tN /N denote the time step increment, tn = nk (0 6 n 6 N), and (un , pn , T n ) be the semi-discrete approximation of (u, p, T ) at tn = nk (n = 0, 1, . . . , N ) about time. If the differential quotients ut and Tt in Problem II at time t = tn are, respectively, approximated by means of the backward difference quotients ∂¯t un = (un − un−1 )/k and ∂¯t T n = (T n − T n−1 )/k, then the semi-discrete approximation scheme of Euler backward one step of Problem II with respect to time is written as follows. Problem III. Find (un , pn , T n ) ∈ X × M × W (1 6 n 6 N ) such that
(un , v ) + ka(un , v ) + ka1 (un−1 , un , v ) − kb(pn , v ) = k(T n j , v ) + (un−1 , v ), b(q, un ) = 0, (T n , φ) + kd(T n , φ) + ka2 (un−1 , T n , φ) = (T n−1 , φ), 0 u = 0, T 0 = ψ(x, y),
∀v ∈ X , ∀q ∈ M , ∀φ ∈ W , (x, y) ∈ Ω .
(9)
It has been proved that Problem III has a unique series of solutions {un , pn , T n }Nn=1 (see [16–18]) and the following results for Problem III are known and classical (see Theorems 5.6 and 5.7 in [18] or [16]). Theorem 2. Under the assumptions of Theorem 1, Problem III has a unique series of solutions (un , pn , T n ) ∈ X × M × W such that
∥un ∥20 + ν k
n
∥∇ ui ∥20 6 C0 ∥ψ∥20 , ∥T n ∥20 + γ0−1 k
i=1
∥u(tn ) −
∥ + kν
un 20
n
n
∥∇ T i ∥20 6 ∥ψ∥20 ,
∥pn ∥0 6 C ∥ψ∥0 ,
i=1
|u(ti ) − | + ∥p(tn ) − ui 21
∥ 6 Ck2 ∥ψ∥20 ,
pn 20
i=1
∥T (tn ) − T n ∥20 + kγ0−1
n
|T (ti ) − T i |21 6 Ck2 ∥ψ∥20 ,
1 6 n 6 N,
i =1
where C0 = 2nk exp(2nk) is a bounded real number, (u, p, T ) ∈ [H01 (Ω ) ∩ H 2 (Ω )]2 × [H 1 (Ω ) ∩ M ] × [H01 (Ω ) ∩ H 2 (Ω )] is the exact solution for the Problem I, and C is a constant independent of k.
74
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Fig. 1. Left chart is a triangle K partitioned into three sub-regions Kz . Right chart is a sample region with dotted lines indicating the corresponding control volume Vz .
2.2. Fully discrete SMFVE formulation for Problem I In order to establish the fully discrete SMFVE formulations for Problem II, it is necessary to introduce an FVE approximation for the spatial variables of Problem III (more details see [2,3,14,15]). Firstly, let ℑh = {K } be a quasi-uniform triangulation of Ω with h = max hK , where hK is the diameter of the triangle K ∈ ℑh (see [5,10,18]). In order to describe the SMFVE formulation, we introduce a dual partition ℑ∗h based on ℑh whose elements are called the control volumes. We construct the control volume in the same way as in [2,3,14,15]. Let zK = (xK , yk ) be the barycenter of K ∈ ℑh . We connect zK with line segments to the midpoints of the edges of K , thus partitioning K into three quadrilaterals Kz (z = (x, y) ∈ Zh (K ), where Zh (K ) are the vertices of K ). Then with each vertex z ∈ Zh = K ∈ℑh Zh (K ) we associate a control volume Vz , which consists of the union of the sub-regions Kz , sharing the vertex z. Finally, we obtain a group of control volumes covering the domain Ω , which is called a barycenter-type dual partition ℑ∗h of the triangulation ℑh (see Fig. 1). We denote the set of interior vertices of Zh by Zh◦ . Since the FE triangulation ℑh is quasi-uniform, the dual partition ℑ∗h is also quasi-uniform (see [2,3,5,10,14,15,18]), namely there exist two positive constants C1 and C2 , being independent of the spatial mesh size h and temporal mesh size k, such that C1 h2 6 mes(Vz ) 6 C2 h2 ,
∀Vz ∈ ℑ∗h .
(10)
Moreover, the barycenter-type dual partition can lead to relatively simple calculations. To this end, the trial function spaces Xh , Wh , and Mh of velocity, temperature, and pressure are, respectively, defined as follows: Xh = vh ∈ X ∩ C (Ω )2 ; vh |K ∈ [P1 (K )]2 , ∀K ∈ ℑh , Wh = wh ∈ W ∩ C (Ω ); wh |K ∈ P1 (K ), ∀K ∈ ℑh , Mh = {qh ∈ M ; qh |K ∈ P1 (K ), ∀K ∈ ℑh } ,
where P1 (K ) is the linear function space on K . Note that they are different from those in [16]. It is obvious that Xh ⊂ X = H01 (Ω )2 and Wh ⊂ W = H01 (Ω ). For (u, T ) ∈ X × W , let (Πh u, ρh T ) be the interpolation projection of (u, T ) onto the trial function spaces Xh × Wh . Then, due to the interpolation theory of Sobolev spaces (see [5,10,14,16,18]), there hold the following error estimates
|u − Πh u|m 6 Ch2−m |u|2 ,
∀u ∈ H 2 (Ω )2 , m = 0, 1,
(11)
|T − ρh T |m 6 Ch
∀T ∈ H (Ω ), m = 0, 1,
(12)
2−m
|T |2 ,
2
where C in this context indicates a positive constant which is possibly different at different occurrences, being independent of the spatial mesh size h and temporal mesh size k. ˜ h of velocity and temperature are, respectively, chosen as follows: The test spaces X˜ h and W X˜ h = vh ∈ L2 (Ω )2 ; vh |Vz ∈ [P0 (Vz )]2 (Vz ∩ ∂ Ω = ∅), vh |Vz = 0 (Vz ∩ ∂ Ω ̸= ∅), ∀Vz ∈ ℑ∗h , ˜ h = wh ∈ L2 (Ω ); wh |Vz ∈ P0 (Vz ) (Vz ∩ ∂ Ω = ∅), wh |Vz = 0 (Vz ∩ ∂ Ω ̸= ∅), ∀Vz ∈ ℑ∗h , W
(13)
where P0 (Vz ) is the constant function space on Vz . In fact, they can be spanned by the following basis functions
φz (x, y) =
1, (x, y) ∈ Vz , 0, elsewhere,
∀z ∈ Zh◦ .
(14)
˜ h , i.e., For (u, w) ∈ X × W , let (Πh∗ u, ρh∗ w) be the interpolation projection of (u, w) onto the test spaces X˜ h × W (Πh∗ u, ρh∗ w) =
(u(z ), w(z ))φz .
(15)
z ∈Zh◦
By the interpolation theory of Sobolev spaces (see [5,10,14,16,18]), we have that
∥u − Πh∗ u∥0 6 Ch|u|1 ;
∥w − ρh∗ w∥0 6 Ch|w|1 .
(16)
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
75
Due to the loss of continuity of the piecewise constant vector functions in X˜ h and the piecewise constant functions in ˜ h on the boundary of two neighboring elements, namely X˜ h ̸⊂ Xh and W ˜ h ̸⊂ Wh , the bilinear forms a(u, v ), b(v , p), and W d(T , φ) and trilinear forms a1 (u, v , w ) and a2 (u, T , φ) in (2) must be revised. Noting that divu = 0, by Green’s formula, there hold the following formulas
Ω
∆u · vdxdy =
∇ p · vdxdy = Ω
Ω
Ω
Ω
∆u · vdxdy = Vz
Vz ∈ℑ∗ h
pv · nds −
∆T φ dxdy =
∆T φ dxdy =
Vz ∈ℑ∗ h
(u · ∇)T φ dxdy =
Vz ∈ℑ∗ h
(v ∇ u) · nds −
∇ u · ∇ vdxdy,
(17)
Vz
Vz ∈ℑ∗ h
pdivvdxdy,
(18)
Vz
Vz ∈ℑ∗ h
(u · ∇)v · wdxdy =
Vz ∈ℑ∗ h
Vz
Vz ∈ℑ∗ h
∂ Vz
Vz ∈ℑ∗ h
∂ Vz
Vz ∈ℑ∗ h
∂ Vz
(∇ T φ) · nds −
Vz ∈ℑ∗ h
(u · ∇)v · wdxdy = Vz
Vz ∈ℑ∗ h
(u · ∇)T φ dxdy = Vz
∂ Vz
Vz ∈ℑ∗ h
∂ Vz
∇ T ∇φ dxdy,
(19)
Vz
(u · n)(v · w )ds −
Vz ∈ℑ∗ h
T φ(u · n)ds −
Vz ∈ℑ∗ h
(u∇ w ) · vdxdy, Vz
(u · ∇φ)T dxdy, Vz
(20)
(21)
where ∂ Vz denotes the line integrals, with the counter clockwise direction, on the boundary ∂ Vz of the dual element; n = (n1 , n2 ) is the unit outer normal vector to ∂ Vz . So bilinear forms a(u, v ), b(p, v ), and d(T , ϕ) and trilinear forms a1 (u, v , w ) and a2 (u, T , φ) in (2) are, respectively, rewritten as
a(u, v ) = ν
∇ u · ∇ vdxdy − ∂ Vz
Vz
Vz ∈ℑ∗ h
b(p, v ) = −
d(T , φ) = γ0
−1
∇ T ∇φ dxdy − γ0
−1
Vz ∈ℑ∗ h
(u∇ w ) · vdxdy + Vz
Vz ∈ℑ∗ h
a2 (u, v , φ) = −
pdivvdxdy ,
Vz
Vz ∈ℑ∗ h
Vz ∈ℑ∗ h
(22)
(23)
Vz
a1 (u, v , w ) = −
(v ∇ u) · nds ,
pv · nds −
∂ Vz
Vz ∈ℑ∗ h
(u · ∇φ)T dxdy +
(∇ T φ) · nds,
(24)
(u · n)(v · w )ds,
(25)
T φ(u · n)ds.
(26)
Vz ∈ℑ∗ h
Vz
∂ Vz
∂ Vz
Vz ∈ℑ∗ h
∂ Vz
˜ h are, respectively, the piecewise constant vector function space and the piecewise constant function space Since X˜ h and W with the characteristic functions of the dual elements Vz as the basis functions, there hold the following expressions a(u, v ) = −ν
Vz ∈ℑ∗ h
a1 (u, v , w ) =
(u · n)(v · w )ds,
∂ Vz
Vz ∈ℑ∗ h
d(T , φ) = −γ0−1
T φ(u · n)ds,
∂ Vz
Vz ∈ℑ∗ h
∀v ∈ X˜ h ;
b(p, v ) = −
Vz ∈ℑ∗ h
Vz ∈ℑ∗ h
a2 (u, T , φ) =
∂ Vz
(v ∇ u) · nds,
∂ Vz
φ∇ T · nds,
∂ Vz
pv · nds,
∀v ∈ X˜ h ;
(27)
∀w ∈ X˜ h ; (28)
˜ h; ∀φ ∈ W ˜ h. ∀φ ∈ W
Then a fully discrete SMFVE formulation for Problem II is written as follows.
(29)
76
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Problem IV. Find (unh , pnh , Thn ) ∈ Xh × Mh × Wh (1 6 n 6 N ) such that
(∂¯ un , Π ∗ v ) + ah (unh , Πh∗ vh ) + a1h (unh−1 , unh , Πh∗ vh ) + bh (pnh , Πh∗ vh ) = (Thn j , Πh∗ vh ), t h nh h b(qh , uh ) + Dh (pnh , qh ) = 0, (∂¯t T n , ρ ∗ wh ) + dh (Thn , ρh∗ wh ) + a2h (unh−1 , Thn , ρh∗ wh ) = 0, 0 h h0 uh = 0, Th = ρh ψ(x, y),
∀vh ∈ Xh , ∀qh ∈ Mh , ∀wh ∈ Wh , (x, y) ∈ Ω ,
(30)
where ah (unh , Πh∗ vh ) = −ν
∂ Vz
Vz ∈ℑ∗ h
(vh (z )∇ unh ) · nds;
dh (Thn , ρh∗ wh ) = −γ0−1
∂ Vz
wh (z )
wh ( z )
K ,2
K ∈ℑh
∂ Vz
qh nds;
(31)
(32)
∇ T · nds; ∂ Vz
Vz ∈ℑ∗ h
Dh (pnh , qh ) = ε
vh ( z )
(unh−1 · n)(unh · vh (z ))ds;
Vz ∈ℑ∗ h
Vz ∈ℑ∗ h
a2h (unh−1 , Thn , ρh∗ wh ) =
Vz ∈ℑ∗ h
a1h (unh−1 , unh , Πh∗ vh ) =
bh (qh , Πh∗ vh ) =
∂ Vz
pnh qh dxdy −
Thn (unh−1 · n)ds;
(33)
K ,1
pnh qh dxdy ,
(34)
here ε is a positive real number and K ,i g (x, y)dxdy indicate an appropriate Gauss integral on K that is exact for polynomials of degree i (i = 1, 2) and g (x, y) = ph qh is a polynomial of degree not more than i (i = 1, 2). Thus, for all test functions qh ∈ Mh , the trial function ph ∈ Mh must be piecewise constant when i = 1. Consequently, we define the L2 -projection ˆ h such that ∀p ∈ L2 (Ω ) satisfies operator ϱh : L2 (Ω ) → W
(p, qh ) = (ϱh p, qh ),
ˆ h, ∀qh ∈ W
(35)
ˆ h ⊂ L2 (Ω ) denotes the piecewise constant space associated with ℑh . The projection operator ϱh has the following where W properties (see [18,5]). ∥ϱh p∥0 6 C ∥p∥0 ,
∀p ∈ L2 (Ω ),
∥p − ϱh p∥0 6 Ch∥p∥1 ,
(36)
∀p ∈ H (Ω ). 1
(37)
Now, using the definition of ϱh , we can rewrite the bilinear form Dh (·, ·) as follows: Dh (ph , qh ) = ε(ph − ϱh ph , qh ) = ε(ph − ϱh ph , qh − ϱh qh ).
(38)
Remark 1. Problem IV is also referred to as Euler backward one step fully discrete SMFVE formulation. 3. Error estimates of solutions for fully discrete SMFVE formulation of Problem II In order to derive the existence, the uniqueness, the stability, and the error estimates of the solutions for fully discrete SMFVE formulation, i.e., Problem IV, it is necessary to introduce some preliminary lemmas. ¯ )2 and Wh ⊂ C (Ω ¯ ) as well (17)–(29), from [2,12,14,16,22,24], we have the following two lemmas. Due to Xh ⊂ C (Ω Lemma 3. There hold the following results: ah (uh , Πh∗ vh ) = a(uh , vh ), a1h (vh , uh , Πh∗ wh ) = a1 (vh , uh , wh ), ∗ a1h (vh , uh , Πh uh ) = 0, ∀uh , vh , wh ∈ Xh ,
(39)
dh (Th , ρh∗ φh ) = d(Th , φh ), a2h (uh , Th , ρh∗ φh ) = a2 (uh , Th , φh ), ∗ a2h (uh , Th , ρh Th ) = 0, ∀Th , φh ∈ Wh , ∀uh ∈ Xh ,
(40)
bh (ph , Πh∗ vh ) = −b(ph , vh ),
(41)
∀vh ∈ Xh , ∀ph ∈ Mh .
Further, ah (uh , Πh vh ) and dh (Th , ρh wh ) are all symmetric, bounded, and positive definite, i.e., ∗
∗
ah (uh , Πh∗ vh ) = ah (vh , Πh∗ uh ),
∀uh , vh ∈ Xh ,
(42)
dh (Th , ρh wh ) = dh (wh , ρh Th ),
∀Th , wh ∈ Wh ,
(43)
∗
∗
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
77
and there exist positive constants h0 , C0 , and C˜ 0 such that, when 0 < h 6 h0 , ah (uh , Πh∗ uh ) > ν|uh |21 ,
|ah (uh , Πh∗ vh )| 6 C0 ∥uh ∥1 ∥vh ∥1 ,
dh (Th , ρh Th ) > γ0 | | , ∗
−1
|dh (Th , ρh wh )| 6 C˜ 0 ∥Th ∥1 ∥wh ∥1 , ∗
Th 21
∀uh , vh ∈ Xh , ∀Th , wh ∈ Wh .
(44) (45)
Lemma 4. There holds the following statement:
(uh , Πh∗ vh ) = (vh , Πh∗ uh ),
∀ u h , vh ∈ X h .
(46)
For any u ∈ H m (Ω )2 (m = 0, 1) and vh ∈ Xh ,
|(u, vh ) − (u, Πh∗ vh )| 6 Chm+n ∥u∥m ∥vh ∥n ,
n = 0, 1.
(47)
Set |||uh |||0 = (uh , Πh∗ uh )1/2 , then ||| · |||0 is equivalent to ∥ · ∥0 on Xh , namely there exist two positive constants C3 and C4 such that C3 ∥uh ∥0 6 |||uh |||0 6 C4 ∥uh ∥0 ,
∀uh ∈ Xh .
(48)
Remark 2. For scalar functions, i.e., if uh and vh in Xh are, respectively, substituted with wh and Th in Wh , the results of Lemma 4 hold (see Theorems 3.2.1 and 5.1.5 in [14]). The following discrete Gronwall Lemma (see [16,18]) is useful for the proofs of the existence, the uniqueness, the stability, and the error estimates of the solutions for Problem IV. Lemma 5 (Discrete Gronwall Lemma). If {an }, {bn }, and {cn } are three positive sequences, and {cn } is monotone, they satisfy
¯ an + bn 6 cn + λ
n −1
ai ,
λ¯ > 0, a0 + b0 6 c0 ,
(49)
i =0
then
¯ an + bn 6 cn exp(nλ),
n > 0.
(50)
By Lemma 3 and using the same approach as the proof of Theorem 4.1 in [2], there holds the following inequality sup
(vh ,qh )∈Uh ×Mh
(unh , Πh∗ vh ) + k[a(unh , vh ) − b(pnh , vh ) + b(qh , unh ) + Dh (pnh , qh )] ∥vh ∥1 + ∥qh ∥0
˜ unh ∥0 + k∥∇ unh ∥0 + k∥pnh ∥0 ), > β(∥
∀(unh , pnh ) ∈ Uh × Mh ,
(51)
where β˜ is independent of h and k. For Problem IV, there hold the following results of the existence, the uniqueness, and the stability of its solutions. Theorem 6. Under the hypotheses of Theorem 2, there exists a unique series of solutions (unh , pnh , Thn )(n = 1, 2, . . . , N ) of the fully discrete SMFVE formulation, i.e., Problem IV satisfying
∥unh ∥20 + ∥Thn ∥20 + k2 ∥pnh ∥20 + k
n (∥uih ∥21 + ∥Thi ∥21 ) 6 C˜ ∥ψ∥20 + ∥ψ∥40 ,
(52)
i=1
where C˜ is a constant only depending on C3 , C4 , and ν . Proof. By Lemma 3, the third equation of (30) and its initial condition are rewritten as the following system of equations:
D (Thn , wh ) = G(wh ), Th0 = ρh ψ(x, y),
∀wh ∈ Wh , n = 1, 2, . . . , N , (x, y) ∈ Ω ,
(53)
where D (Thn , wh ) = (Thn , ρh∗ wh ) + kd(Thn , wh ) + ka2 (unh−1 , Thn , wh ) and G(wh ) = (Thn−1 , ρh∗ wh ). For given uhn−1 and Thn−1 , D (·, ·) and G(·) are obviously bounded bilinear and linear functions, respectively. Moreover, by using (4) and Lemmas 3 and 4, we have that D (Thn , Thn ) = (Thn , ρh∗ Thn ) + kd(Thn , Thn ) + ka2 (unh−1 , Thn , Thn ) = (Thn , ρh∗ Thn ) + kd(Thn , Thn ) > α1 ∥Thn ∥21 (where α1 = min{C3 , kγ0−1 }), namely D (·, ·) is the positive definite bilinear function. Therefore, the system of Eq. (53) has a unique series of solutions {Thn }Nn=1 . Taking wh = Thn in (53), by using Lemma 3, Hölder inequality, and Cauchy inequality, we obtain that
|||Thn |||20 + kγ0−1 ∥∇ Thn ∥20 6 D (Thn , Thn ) = G(Thn ) 6
1 2
2
(|||Thn |||20 + |||Thn−1 |||0 ),
n = 1, 2, . . . , N .
(54)
78
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Further, we obtain that 2
|||Thn |||20 + 2kγ0−1 ∥∇ Thn ∥20 6 |||Thn−1 |||0 ,
n = 1, 2, . . . , N .
(55)
Summing (55) from 1 to n and using ∥Th0 ∥0 = ∥ρh ψ(x, y)∥0 6 ∥ψ(x, y)∥0 yields that
∥Thn ∥20 + k
n
∥∇ Thi ∥20 6 C0 ∥ψ∥20 ,
n = 1, 2, . . . , N ,
(56)
i=1
where C0 = C4 min 1/(2γ0 ), C3−1 . For known Thn , the first and second equations as well as their initial condition in Problem IV have a unique series of solutions (unh , pnh )(n = 1, 2, . . . , N ) by means of the theory of mixed FE methods (see [18,5]) and inequality (51). Taking vh = unh in the first equation of Problem IV and qh = pnh in the second equation of Problem IV, by using Lemmas 3 and 4, Hölder inequality, and Cauchy inequality, we obtain that
1 2
2
(|||unh |||20 − |||unh−1 |||0 ) + kν|unh |21 + kε∥pnh − ρh pnh ∥20 = k(Thn j , Πh∗ unh ) 6
k 2ν
∥Thn ∥2−1 +
kν 2
|unh |21 .
(57)
Further, we obtain that 2
|||unh |||20 − |||unh−1 |||0 + kν|unh |21 6 kν −1 ∥Thn ∥2−1 .
(58)
Summing (58) from 1 to n yields that
|||unh |||20 + kν
n
|uih |21 6 kν −1
i =1
n
∥Thi ∥2−1 6 kν −1
i =1
n
∥Thi ∥20 6 ν −1 C0 ∥ψ∥20 .
(59)
i =1
By the first equation of Problem IV and using (51), Hölder inequality, and (59), we obtain that
β˜ k∥pnh ∥0 6
sup
(vh ,qh )∈Uh ×Mh
6 sup
(unh , Πh∗ vh ) + k[a(unh , vh ) − b(pnh , vh ) + b(qh , unh ) + Dh (pnh , qh )] ∥vh ∥1 + ∥qh ∥0
k(Thn j , Πh∗ vh ) − ka1 (unh−1 , unh , Πh∗ vh ) + (unh−1 , Πh∗ vh )
∥v h ∥1
vh ∈Uh
6 k∥
Thn −1
∥
+ kN0 ∥∇
∥ ∥∇ unh ∥0 + ∥uhn−1 ∥−1 .
unh−1 0
(60)
Combining (60) and (59) with (56) yields (52) and completes the proof of Theorem 6. For (u , p ) ∈ X × M and n
n
unh
∈ Xh , put
Ah ((Sh un , Qh pn ); (vh , qh )) = a(Sh un , vh ) − b(Qh pn , vh ) + b(qh , Sh un ) + a1 (uhn−1 , Sh un , vh ),
∀(vh , qh ) ∈ Xh × Mh ,
(61)
A((un , pn ); (vh , qh )) = a(un , vh ) − b(pn , vh ) + b(qh , un ) + a1 (un−1 , un , vh ),
∀(vh , qh ) ∈ Xh × Mh .
By the FE methods (see, e.g., [16,18]) for the non-stationary Navier–Stokes equations, there holds the following Lemma 7. Lemma 7. Let (Sh , Qh ) : X × M → Xh × Mh be the stabilized Navier–Stokes projection, i.e., for (un , pn ) ∈ X × M, there exist (Sh un , Qh pn ) (n = 0, 1, . . . , N) such that kAh ((Sh un , Qh pn ); (vh , qh )) + (Sh un − Sh un−1 , vh ) + kDh (Qh pn , qh ) = kA((un , pn ); (vh , qh )) + (un − un−1 , vh ),
∀(vh , qh ) ∈ Xh × Mh ,
(62)
∥Sh un ∥1 + ∥Qh pn ∥0 6 C (∥un ∥1 + ∥pn ∥0 ).
(63)
If (u , p ) ∈ H (Ω ) × H (Ω ) and h = O(k), there hold the following equalities n
n
2
2
1
∥un − Sh un ∥20 + k2 ∥pn − Qh pn ∥20 + kν
n i =1
∥ui − Sh ui ∥21 6 Ch5
n (∥ui ∥22 + ∥pi ∥21 ),
0 6 n 6 N.
(64)
i =1
And if (un , pn , T n ) ∈ H 2 (Ω )2 × H 1 (Ω ) × H 2 (Ω ) (n = 1, 2, . . . , N) are the solutions to Problem III, then there exists a positive constant C independent of the discretization parameters such that
∥un − Sh un ∥20 + k2 ∥pn − Qh pn ∥20 + kν
n i =1
∥ui − Sh ui ∥21 6 Ch4 ∥ψ∥22 ,
0 6 n 6 N.
(65)
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
79
By the classical FE methods (see, e.g., [10,18]) for elliptic equations, we have the following Lemma 8. Lemma 8. Let Rh : W → Wh be a generalized Ritz projection, i.e., for given uhn−1 ∈ Xh , T n−1 ∈ W , and Thn−1 ∈ Wh , and for any T n ∈ W (n = 1, 2, . . . , N), there exist Rh T n ∈ Wh (n = 1, 2, . . . , N) such that
(Rh T n , wh ) + kd(Rh T n , wh ) + ka2 (unh−1 , Rh T n , wh ) − (Rh T n−1 , wh ) = (T n , wh ) + kd(T n , wh ) + ka2 (un−1 , T n , wh ) − (T n−1 , wh ),
∀wh ∈ Wh , n = 1, 2, . . . , N .
(66)
If (un , pn , T n ) (n = 1, 2, . . . , N) are the solutions of Problem III and T n ∈ H 2 (Ω ) ∩ W , there hold the following inequalities
∥Rh T n ∥0 + k1/2 ∥Rh T n ∥1 6 C ∥T n ∥1 , ∥Rh T n − T n ∥20 + kγ0−1
n
n = 0, 1, 2, . . . , N ,
|Rh T i − T i |21 6 Ch4 ∥ψ∥22 ,
(67)
n = 0, 1, 2, . . . , N .
(68)
i =1
Theorem 9. Let (u, p, T ) be the solution to Problem II and (unh , pnh , Thn ) the solution to Problem IV. Under the hypotheses of Theorems 1 and 2, if p0h = p0 = 0 (or p0h = Qh p0 ), h = O(k), u0 ∈ H 2 (Ω )2 , and ψ ∈ H 1 (Ω ), then there hold that
∥u(tn ) − unh ∥0 + ∥T (tn ) − Thn ∥0 + k[∥p(tn ) − pnh ∥0 + ∥u(tn ) − unh ∥1 + ∥T (tn ) − Thn ∥1 ] 6 Ch2 (∥ψ∥1 ) , n = 1, 2, . . . , N .
(69)
Proof. Subtracting Problem IV from Problem III taking v = vh and q = qh , using Lemmas 3 and 4, we obtain the system of error equations:
n (u − unh , vh ) + (unh − Πh∗ unh , vh − Πh∗ vh ) + ka(un − unh , vh ) + ka1 (un−1 , un , vh ) − ka1 (uhn−1 , unh , vh ) −kb(pn − pnh , vh ) = k(j (T n − Thn ), vh ) − k(jThn , Πh∗ vh − vh ) + (un−1 − uhn−1 , vh ) − (unh−1 , Πh∗ vh − vh ), ∀vh ∈ Xh , n = 1, 2, . . . , N , b(qh , un − unh ) − ε(pnh − ϱh pnh , qh − ϱh qh ) = 0, ∀qh ∈ Mh , n = 1, 2, . . . , N , (T n − Thn , wh ) − (Thn , ρh∗ wh − wh ) + kd(T n − Thn , wh ) + ka2 (un−1 , T n , wh ) − ka2 (uhn−1 , Thn , wh ) = (T n−1 − Thn−1 , wh ) − (Thn−1 , ρh∗ wh − wh ), ∀wh ∈ Wh , n = 1, 2, . . . , N , 0 u − u0h = 0, T 0 − Th0 = ψ(x, y) − ρh ψ(x, y), (x, y) ∈ Ω .
(70)
Let En = Sh un − unh . By using (62) and the system of error equations (70), we obtain that
∥En ∥20 + kν|En |21 = (En , En ) + ka(En , En ) = [(Sh un − un , En ) + ka(Sh un − un , En )] + [(un − unh , En ) + ka(un − unh , En )] = [(Sh un−1 − un−1 , En ) + kb(Qh pn − pn , Eh ) + ka1 (un−1 , un , En ) − ka1 (uhn−1 , Sh un , En )] + [kb(pn − pnh , En ) − ka1 (un−1 , un , En ) + ka1 (unh−1 , unh , En ) + (unh , Πh∗ En − En ) + k(j (T n − Thn ), En ) − k(jThn , Πh∗ En − En ) + (un−1 − uhn−1 , En ) − (uhn−1 , Πh∗ En − En )] = (En−1 , En ) + k(j (T n − Thn ), En ) − k(jThn , Πh∗ En − En ) + (unh − unh−1 , Πh∗ En − En ) + kb(Qh pn − pnh , Eh ).
(71)
By using Hölder inequality, Cauchy inequality, (47), (65), Lemmas 3 and 4, and Taylor’s formula, if h = O(k) we get that
|(En−1 , En )| 6 ∥En−1 ∥0 ∥En ∥0 6
1 2
1
∥En−1 ∥20 + ∥En ∥20 ,
(72)
2
|k(j (T n − Thn ), En ) − k(jThn , Πh∗ En − En )| 6 Ck∥T n − Thn ∥0 ∥En ∥0 + Ckh2 ∥Thn ∥0 ∥∇ En ∥0 νk 6 Ck∥T n − Thn ∥20 + Ck∥En ∥20 + Ckh4 + |En |21 , 4
|(unh − uhn−1 , Πh∗ En − En )| = (Sh un − un + un − un−1 + un−1 − Sh un−1 + En−1 − En , Πh∗ En − En ) 6 Ch(∥Sh un − un ∥20 + h∥un − un−1 ∥21 + ∥un−1 νk − Sh un−1 ∥20 + ∥En−1 ∥20 + ∥En ∥20 ) + |En |21 4 νk 6 Ch(∥En ∥20 + ∥En−1 ∥20 ) + Ch5 ∥ψ∥22 + Ck2 h3 ∥ut ∥2L∞ (H 1 ) + |En |21 . 4
(73)
(74)
80
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Note that b(qh , Sh un − un ) = −kε(Qh pn − ϱh (Qh pn ), qh − ϱh (qh )). By the property of operator ϱh and the second equation of (70), we have that b(Qh pn − pnh , Sh un − unh ) = b(Qh pn − pnh , Sh un − un ) + b(Qh pn − pnh , un − unh )
ε = − (Qh pn − pnh − ϱh (Qh pn − pnh ), Qh pn − pnh − ϱh (Qh pn − pnh )) 2
ε − (Qh pn−1 − pnh−1 − ϱh (Qh pn−1 − pnh−1 ), Qh pn − pnh − ϱh (Qh pn − pnh )) 2
ε
ε
6 − ∥Qh pn − pnh − ϱh (Qh pn − pnh )∥20 + ∥Qh pn−1 − ph 4 4
n −1
− ϱh (Qh pn−1 − phn−1 )∥20 .
(75)
Combining (72)–(75) with (71) yields that
ε ε ∥En ∥20 + kν|En |21 + ∥Qh pn − pnh − ϱh (Qh pn − pnh )∥20 − ∥Qh pn−1 − pnh−1 − ϱh (Qh pn−1 − pnh−1 )∥20 2
6∥
En−1 20
2
4
n
∥ + Ckh + Ck∥T −
∥ + Ck(∥ ∥ + ∥En−1 ∥20 ).
Thn 20
En 20
(76)
Summing (76) from 1 to n, if h is sufficiently small such that Ch 6 1/2 in (76) and p0h = p0 = 0 (or p0h = Qh p0 ), we obtain that
∥En ∥20 + kν
n
|Ei |21 6 Cnkh4 + Ck
i =1
n
∥T i − Thi ∥20 + Ck
i=1
n
∥Ei ∥20 .
(77)
i =0
If k is sufficiently small such that Ck 6 1/2 in (77), we obtain from (77) that
∥En ∥20 + 2kν
n
|Ei |21 6 Ch4 + Ck
i=1
n
∥T i − Thi ∥20 + Ck
n −1
i=1
∥Ei ∥20 .
(78)
i =0
Applying Gronwall Lemma to (78) yields that
∥ ∥ + kν En 20
n
En 21
4
| | 6C h +k
i =1
n
i
∥T −
Thi 20
∥
exp(Ckn).
(79)
i=1
By using triangle inequality, (64), and Lemma 7, we get that n
∥u −
∥ + kν
unh 20
n
i
|u −
uih 21
4
| 6C h +k
i=1
n
i
∥T −
Thi 20
∥
.
(80)
i =1
By using (51) and the system of error equations (70), we have that
β˜ k∥Qh pn − pnh ∥0 6
sup
(vh ,qh )∈Uh ×Mh
(Sh un − unh , Πh∗ vh )
+ k a(Sh un − unh , vh ) − b(Qh pn − pnh , vh ) + b(qh , Sh un − unh ) + Dh (pnh , qh ) /(∥vh ∥1 + ∥qh ∥0 ) |b(Qh pn − pn , vh )| |b(pn − pnh , vh )| + k sup ∥v h ∥1 ∥v h ∥1 vh ∈Uh vh ∈Uh −1 n n n n n n ∗ 6 C (k∥Qh p − p ∥0 + ∥u − uh ∥0 ) + k k (u − uh , Πh vh ) − (un−1 − unh−1 , Πh∗ vh ) 6 C ∥Sh un − unh ∥0 + k sup
+ a(un − unh , vh ) + a1 (un−1 , un − unh , vh ) + a1 (un−1 − unh−1 , unh , vh ) − (j (T n − Thn ), Πh∗ vh ) /∥vh ∥1 n −1 6 Ck ∥Qh pn − pn ∥0 + k−1 ∥un − unh ∥0 + ∥un−1 − uh ∥0 + |un−1 − unh−1 |1 + |un − unh |1 + ∥T n − Thn ∥0 .
(81)
Combining (80) with (81) yields that
∥Qh pn − pnh ∥0 6 Ch ∥ψ∥2 + ∥T n − Thn ∥0 + k
n i =1
1/2 ∥T i − Thi ∥20
.
(82)
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
81
By using triangle inequality, (65) and (82), we get that
∥pn − pnh ∥0 6 Ch ∥ψ∥2 + ∥T n − Thn ∥0 + k
n
1/2 . ∥T i − Thi ∥20
(83)
i =1
Let en = Rh T n − Thn . By using the system of error equations (70), (4), (66)–(68) in Lemma 8, Lemma 3, and (16), we obtain that
∥en ∥20 + kγ0−1 |en |21 = (en , en ) + kd(en , en ) = [(Rh T n − T n , en ) + kd(Rh T n − T n , en )] + [(T n − Thn , ρh∗ en ) + kd(T n − Thn , en )] = [(Rh T n−1 − T n−1 , en ) + ka2 (un−1 , T n , en ) − ka2 (unh−1 , Rh T n , en )] + [(Thn , ρh∗ en − en ) + ka2 (unh−1 , Thn , en ) − ka2 (un−1 , T n , en ) + (T n−1 − Thn−1 , en ) − (Thn−1 , ρh∗ en − en )] = (en−1 , en ) + (Thn − Thn−1 , ρh∗ en − en ) 6 ∥en−1 ∥0 ∥en ∥0 + Ch(∥en ∥0 + ∥Rh T n − T n ∥0 + h∥T n − T n−1 ∥1
+ ∥T n−1 − Rh T n−1 ∥0 + ∥en−1 ∥0 )|en |1 6
1 2
k ∥en−1 ∥20 + ∥en ∥20 + Ch h4 + k2 h2 + ∥en ∥20 + ∥en−1 ∥20 + |en |21 . 2γ0
(84)
Further, we get from (84) that
∥en ∥20 + kγ0−1 |en |21 6 Ch h4 + k2 h2 + ∥en ∥20 + ∥en−1 ∥20 + ∥en−1 ∥20 .
(85)
Summing (85) from 1 to n and using (68) and (12) yield that
∥en ∥20 + kγ0−1
n
|ei |21 6 Cnh(h4 + k2 ) + ∥e0 ∥20 + Ck
i=1
n
∥ei ∥20
i =1
6 C (h4 + k2 h2 ) + Ck
n
∥ei ∥20 + C ∥Rh ψ − ψ∥20 + C ∥ψ − ρh ψ∥20
i =1
6 C (h4 + k2 h2 ) + Ck
n
∥ei ∥20 .
(86)
i =1
If k is sufficiently small such that Ck 6 1/2 in (86), we obtain from (86) that
∥en ∥20 + kγ0−1
n
|ei |21 6 Ch4 + Ck
i=1
n −1
∥ei ∥20 .
(87)
i =0
Applying the Gronwall Lemma to (87) yields that
∥en ∥20 + kγ0−1
n
|ei |21 6 Ch4 exp(Ck) 6 Ch4 .
(88)
i=1
By using triangle inequality, (88), and (68), we obtain that
∥T n − Thn ∥20 + kγ0−1
n
|T i − Thi |21 6 Ch4 .
(89)
i =1
Combining (79) with (89) yields that
∥un − unh ∥20 + kν
n
|ui − uih |21 6 Ch4 .
(90)
i =1
Combining (89) with (83) yields that
∥pn − pnh ∥0 6 Ch. Combining (90) and (91) with (89) yields (69) which completes the proof of Theorem 9.
(91)
82
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Combining Theorem 2 with Theorem 9 yields the following conclusion. Theorem 10. Under the hypotheses of Theorems 2 and 9, there hold the error estimates between the solution (u, p, T ) to Problem II and the solutions (unh , pnh , Thn )(n = 1, 2, . . . , N ) to Problem IV as follows
∥u(tn ) − unh ∥20 + ∥T (tn ) − Thn ∥20 + k
n (|u(ti ) − uih |21 + |T (ti ) − Thi |21 ) 6 C (k2 + h4 ),
(92)
i =1
∥pn − pnh ∥0 6 C (k + h),
(93)
where C is a constant which is only dependent on γ0 , ν , Ω , ψ , and maximum total time upper bound tN , but independent of k and h. 4. Some numerical experiments In this section, we present some numerical experiments with a physical model of cavity flow and Reynolds number Re = 7.1 × 103 and Pr = 0.71 by the fully discrete SMFVE Formulation, i.e., Problem IV to validate that the numerical computational results are consistent with the theoretical conclusions. Moreover, it is shown that the fully discrete SMFVE formulation based on two local Gauss integrals and parameter-free is feasible and efficient for finding the numerical solutions of the non-stationary conduction–convection problem and its numerical solutions are more stable than those of the classical fully discrete FVE formulation without any stabilization in [16]. Especially, its degrees of freedom are far fewer than those of the classical fully discrete FVE formulation without any stabilization in [16]. 4.1. Numerical experiment for the fully discrete SMFVE formulation
¯ = [0, 1] × [0, 1]. We first divide the cavity into 100 × 100 = 10 000 small Let the side length of the cavity be 1 and Ω squares with side length △x = △y = 0.01, and then link diagonal of the√square to divide each square into two triangles in the same direction, which constitutes triangularizations ℑh with h = 2 × 10−2 . The dual decomposition ℑ∗h is taken as barycenter dual decomposition, namely the barycenter of the right triangle K ∈ ℑh is taken as the node of the dual decomposition. We take the time step increment as k = 0.01. Let the initial value and the boundary value of u = (u1 , u2 ) be equal to 0. And let T = 0 on left and lower boundary of the cavity, ∂ T /∂ y = 0 on upper boundary of the cavity, and T = 4y(1 − y) on the right boundary of the cavity. The initial value of T is 0. We find a SMFVE numerical solution (unh , pnh , Thn ) by means of Problem IV when n = 200 (i.e., t = 2), whose unh , pnh , and Thn are depicted graphically on the left charts in Figs. 2–4, respectively. The degrees of freedom for the fully discrete SMFVE Formulation, i.e., Problem IV on each time level are about 4Np = 4 × 104 (where Np is the number of nodal points in ℑh , see [5,10,18]). 4.2. Numerical experiment for the classical FVE formulation without any stabilization
¯ = [0, 1] × [0, 1]. We also divide the cavity into 100 × 100 = 10000 small Let still the side length of the cavity be 1 and Ω squares with side length △x = △y = 0.01 and link diagonal √ of the square to divide each square into two triangles in the same direction and to constitute triangularizations ℑh (h = 2 × 10−2 ). In order to make BB condition satisfied,√ it is necessary ˜ h with h = 2 × 10−2 /2. to link three midpoints on sides of each triangle in ℑh , which constitutes triangularizations ℑ ˜ ∗h is also taken as barycenter dual decomposition, namely the barycenter of the triangle K ∈ ℑ˜ h The dual decomposition ℑ is taken as the node of the dual decomposition. The time step increment, the initial values and the boundary values are the same as those in Section 4.1. The classical fully discrete FVE formulation without any stabilization is written as follows (see [16]). Problem V. Find (unh , pnh , Thn ) ∈ X h × M h × W h (1 6 n 6 N ) such that
(∂¯ un , Π ∗ v ) + ah (unh , Πh∗ vh ) + a1h (unh−1 , unh , Πh∗ vh ) + bh (pnh , Πh∗ vh ) = (jThn , Πh∗ vh ), t h nh h b(qh , uh ) = 0, (∂¯t Thn , ρh∗ wh ) + dh (Thn , ρh∗ wh ) + a2h (unh−1 , Thn , ρh∗ wh ) = 0, 0 uh = 0, Th0 = ρh ψ(x, y),
∀vh ∈ X h , ∀qh ∈ M h , ∀wh ∈ W h , (x, y) ∈ Ω ,
where all other signs are the same as those in Problem IV except for trial function spaces X h , W h , and M h of velocity, temperature, and pressure being, respectively, defined as follows in order to make BB condition satisfied (see [16]):
˜h , X h = vh ∈ X ∩ C (Ω )2 ; vh |K ∈ P12 (K ), ∀K ∈ ℑ
˜h , W h = wh ∈ W ∩ C (Ω ); wh |K ∈ P1 (K ), ∀K ∈ ℑ M h = qh ∈ L20 (Ω ); qh |K ∈ P0 (K ), ∀K ∈ ℑh .
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
83
Fig. 2. When Re = 7.1 × 103 and Pr = 0.71, the left chart and the right chart are the SMFVE solution and the classical FVE solution of the velocity u at the time level t = 2, respectively.
Fig. 3. When Re = 7.1 × 103 and Pr = 0.71, the left chart and the right chart are the SMFVE solution and the classical FVE solution of the temperature T at the time level t = 2, respectively.
Thus, the degrees of freedom for the classical fully discrete FVE Formulation, i.e., Problem V on each time level are about 3Np + NK + 3Ns ≈ 14Np = 14 × 104 (where Np , NK , and Ns are the numbers of nodal points, elements, and sides in ℑh , respectively) since Np : NK : Ns ≈ 1 : 2 : 3 (see [5,10,18]). They are 3.5 times those for the fully discrete SMFVE Formulation. We also find a FVE numerical solution (unh , pnh , Thn ) to Problem V when n = 200 (i.e., t = 2), whose unh , pnh , and Thn are depicted graphically on the right charts in Figs. 2–4, respectively. The curves of the upper and lower charts in Fig. 5 are the relative errors (log 10) of the fully discrete SMFVE solutions and the classical fully discrete FVE solutions at time t ∈ (0, 2] with respect to L2 -norm, respectively. Due to the truncation error accumulation in computational process, the errors of numerical solutions appear increase (see Fig. 5), but the truncation error accumulation for the fully discrete SMFVE formulation Problem IV is far smaller than that for the classical fully discrete FVE formulation with any stabilization and the relative errors (which does not exceed 2×10−2 ) of SMFVE numerical solutions are far smaller than those of classical FVE solutions (also see Fig. 5). Moreover, it is shown that the numerical computational results are consistent with the theoretical conclusions for the fully discrete SMFVE formulation, since the theoretical and numerical errors achieve O(10−2 ). Since the fully discrete SMFVE formulation adopts the stabilization based on two local Gauss integrals and parameter-free, its numerical solutions are more stable than those of the fully discrete FVE formulation without any stabilization (also see Figs. 2–4). For example, the chart of velocity u of the classical FVE solution at t = 2 has obvious chaos (see the right chart in Fig. 2) so that its numerical computational effectiveness has no far better than that of SMFVE formulation (see the left chart in Fig. 2). Especially, the fully discrete SMFVE formulation can reduce great many
84
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
Fig. 4. When Re = 7.1 × 103 and Pr = 0.71, the left chart and the right chart are the SMFVE solution and the classical FVE solution of the pressure p at the time level t = 2, respectively.
Fig. 5. When Re = 7.1 × 103 and Pr = 0.71, the upper chart and the lower chart are the relative errors (log 10) of the SMFVE solutions and the classical FVE solution of the velocity u, the temperature T , and the pressure p at the time t ∈ (0, 2], respectively.
of degrees of freedom (its degrees of freedom are by a factor of 3.5 times fewer than those for the classical fully discrete FVE formulation in [16]). Therefore, the fully discrete SMFVE formulation is far better than the classical fully discrete FVE formulation without any stabilization in [16].
Z. Luo et al. / J. Math. Anal. Appl. 404 (2013) 71–85
85
5. Conclusions and discussions We have established directly a fully discrete SMFVE formulation from the semi-discrete formulation with respect to time for the non-stationary conduction–convection problem. Thus, we could avoid the semi-discrete SMFVE formulation about spatial variables. We have also provided the existence and stabilization of the fully discrete SMFVE solutions and the error estimates between the fully discrete SMFVE solutions and the accurate solution. We have also given some numerical experiments illustrating that the numerical computational results are consistent with the theoretical results obtained previously, validating that the fully discrete SMFVE formulation is feasible and efficient for finding the numerical solutions of the non-stationary conduction–convection problem, and showing that the errors of the fully discrete SMFVE solutions are far smaller than those of the classical fully discrete FVE solutions, its degrees of freedom are far fewer than those of the classical fully discrete FVE formulation without any stabilization in [16], and the numerical solutions for the fully discrete SMFVE formulation are more stable than those for the classical fully discrete FVE formulation without any stabilization in [16]. Especially, the trial function spaces of velocity and pressure for the fully discrete SMFVE formulation satisfy the BB condition naturally, which is a main improvement for method in [16]. Consequently, the fully discrete SMFVE formulation are far better than the classical fully discrete FVE formulation without any stabilization in [16]. Though an FVE method for time-dependent Stokes equations (see [22]) and a penalty FVE method for transient Navier–Stokes equations (see [12]) have been provided, the non-stationary conduction–convection problem as indicated in [16] is different from either of the above. The Stokes equations are linear, while our conduction–convection problem is nonlinear. Moreover, the non-stationary conduction–convection problem includes a nonlinear energy equation of the temperature field with the velocity vector field coupled, which is more complex than the transient Navier–Stokes equations without energy equation in [12]. Therefore, the study of the SMFVE method for the non-stationary conduction–convection problem has far more difficulties, more important, more serviceable, and more challenging than those for the timedependent Stokes equations and the transient Navier–Stokes equations. It is the improvement and innovation for the existing methods (see, e.g., [12,16,22] or others). Acknowledgments The authors cordially thank all the Reviewers, Editors, and Professor & Editor in Chief Goong Chen for their very valuable suggestions that led to improvements of our manuscript. References [1] R.A. Adams, Sobolev Spaces, Academic Press, New York, 1975. [2] J. An, P. Sun, Z.D. Luo, X.M Huang, A stabilized fully discrete finite volume element formulation for non-stationary Stokes equations, Math. Numer. Sin. 33 (2) (2011) 213–224. [3] R.E. Bank, D.J. Rose, Some error estimates for the box methods, SIAM J. Numer. Anal. 24 (4) (1987) 777–787. [4] P. Blanc, R. Eymerd, R. Herbin, A error estimate for finite volume methods for the Stokes equations on equilateral triangular meshes, Numer. Methods Partial Differential Equations 20 (2004) 907–918. [5] F. Brezzi, M. Fortin, Mixed and Hybrid Finite Element Methods, Springer–Verlag, New York, 1991. [6] Z. Cai, S. McCormick, On the accuracy of the finite volume element method for diffusion equations on composite grid, SIAM J. Numer. Anal. 27 (3) (1990) 636–655. [7] P. Chatzipantelidis, R.D. Lazarrov, V. Thome, Error estimates for a finite volume element method for parabolic equations in convex in polygonal domains, Numer. Methods Partial Differential Equations 20 (2004) 650–674. [8] S.H. Chou, D.Y. Kwak, A covolume method based on rotated bilinears for the generalized Stokes problem, SIAM J. Numer. Anal. 35 (1998) 494–507. [9] M.A. Christon, P.M. Gresho, S.B. Sutton, Computational predictability of time-dependent natural convection flows in enclosures (including a benchmark solution), in: MIT Special Issue on Thermal Convection, Internat. J. Numer. Methods Fluids 40 (8) (2002) 953–980 (special issue). [10] P.G. Ciarlet, The Finite Element Method for Elliptic Problems, North-Holland, Amsterdam, 1978. [11] E. DiBenedetto, A. Friedman, Conduction–convection problems with change of phase, J. Differential Equations 62 (1986) 129–185. [12] G.L. He, Y.N. He, Z.X. Chen, A penalty finite volume method for the transient Navier–Stokes equations, Appl. Numer. Math. 58 (2008) 1583–1613. [13] W.P. Jones, K.P. Menziest, Analysis of the cell-centred finite volume method for the diffusion equation, J. Comput. Phys. 165 (2000) 45–68. [14] R.H. Li, Z.Y. Chen, W. Wu, Generalized Difference Methods for Differential Equations-Numerical Analysis of Finite Volume Methods, in: Monographs and Textbooks in Pure and Applied Mathematics, vol. 226, Marcel Dekker Inc., New York, 2000. [15] Y. Li, R.H. Li, Generalized difference methods on arbitrary quadrilateral networks, J. Comput. Math. 17 (1999) 653–672. [16] H. Li, Z.D. Luo, P. Sun, J. An, A finite volume element formulation and error analysis for the non-stationary conduction–convection problem, J. Math. Anal. Appl. 396 (2012) 864–879. [17] Z.D. Luo, The mixed finite element method for the non-stationary conduction–convection problems, Math. Numer. Sin. 20 (1) (1998) 69–88. [18] Z.D. Luo, The bases and Applications of Mixed Finite Element Methods, Chinese Science Press, Beijing, 2006. [19] Z.D. Luo, J. Chen, I.M. Navon, J. Zhu, An optimizing reduced PLSMFE formulation for non-stationary conduction–convection problems, Internat. J. Numer. Methods. Fluids 60 (2009) 409–436. [20] Z.D. Luo, Z.H. X, J. Chen, A reduced MFE formulation based on POD for the non-stationary conduction–convection problems, Acta Math. sci. 31 (5) (2011) 765–1785. [21] E. Süli, Convergence of finite volume schemes for Poissons equation on nonuniform meshes, SIAM J. Numer. Anal. 28 (5) (1991) 1419–1430. [22] M. Yang, H.L. Song, A postprocessing finite volume method for time-dependent Stokes equations, Appl. Numer. Math. 59 (2009) 1922–1932. [23] X. Ye, On the relation between finite volume and finite element methods applied to the Stokes equations, Numer. Methods Partial Differential Equations 17 (2001) 440–453. [24] M. Zlámal, Finite element methods for nonlinear parabolic equations, RAIRO Anal. Numer. 11 (1) (1977) 93–107.