- About this Journal ·
- Abstracting and Indexing ·
- Aims and Scope ·
- Annual Issues ·
- Article Processing Charges ·
- Articles in Press ·
- Author Guidelines ·
- Bibliographic Information ·
- Citations to this Journal ·
- Contact Information ·
- Editorial Board ·
- Editorial Workflow ·
- Free eTOC Alerts ·
- Publication Ethics ·
- Reviewers Acknowledgment ·
- Submit a Manuscript ·
- Subscription Information ·
- Table of Contents

Journal of Applied Mathematics

Volume 2013 (2013), Article ID 372906, 21 pages

http://dx.doi.org/10.1155/2013/372906

## Fully Discrete Finite Element Approximation for the Stabilized Gauge-Uzawa Method to Solve the Boussinesq Equations

Department of Mathematics, Kangwon National University, Chuncheon 200-701, Republic of Korea

Received 9 February 2013; Accepted 11 March 2013

Academic Editor: Jinyun Yuan

Copyright © 2013 Jae-Hong Pyo. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

The stabilized Gauge-Uzawa method (SGUM), which is a 2nd-order projection type algorithm used to solve Navier-Stokes equations, has been newly constructed in the work of Pyo, 2013. In this paper, we apply the SGUM to the evolution Boussinesq equations, which model the thermal driven motion of incompressible fluids. We prove that SGUM is unconditionally stable, and we perform error estimations on the fully discrete finite element space via variational approach for the velocity, pressure, and temperature, the three physical unknowns. We conclude with numerical tests to check accuracy and physically relevant numerical simulations, the Bénard convection problem and the thermal driven cavity flow.

#### 1. Introduction

The stabilized Gauge-Uzawa method (SGUM) is a 2nd-order projection type method to solve the evolution Navier-Stokes equations. In this paper, we extend SGUM to the evolution Boussinesq equations: given a bounded polygon in with or , with initial conditions , in , vanishing Dirichlet boundary conditions , on , and pressure mean-value . The forcing functions and are given, and is the vector of gravitational acceleration. The nondimensional numbers and are reciprocal of the Reynolds and Prandtl numbers, respectively, whereas is the Grashof number. The Boussinesq system (1) describes fluid motion due to density differences which are in turn induced by temperature gradients: hot and thus less dense fluid tends to rise against gravity and cooler fluid falls in its place. The simplest governing equations are thus the Navier-Stokes equations for motion of an incompressible fluid, with forcing due to buoyancy and the heat equation for diffusion and transport of heat. Density differences are thus ignored altogether except for buoyancy.

The projection type methods are representative solvers for the incompressible flows, and the Gauge-Uzawa method is a typical projection method. The Gauge-Uzawa method was constructed in [1] to solve Navier-Stokes equations and extended to more complicated problems, the Boussinesq equations in [2] and the nonconstant density Navier-Stokes equations in [3]. However, most of studies for the Gauge-Uzawa method have been limited only for the first-order accuracy backward Euler time marching algorithm. The second-order Gauge-Uzawa method using BDF2 scheme was introduced in [4] and proved superiority for accuracy on the normal mode space, but we could not get any theoretical proof via energy estimate even stability and we suffer from weak stability performance on the numerical test. Recently, we construct SGUM in [5] which is unconditionally stable for semidiscrete level to solve the Navier-Stokes equations. The goal of this paper is to extend SGUM to the Boussinesq equations (1), which model the motion of an incompressible viscous fluid due to thermal effects [6, 7]. We will estimate errors and stability on the fully discrete finite element space. The main difficulties in the fully discrete estimation arise from losing the cancellation law due to the failing of the divergence free condition of the discrete velocity function. The strategy of projection type methods computes first an artificial velocity and then decomposes it to divergence free velocity and curl free functions. However, the divergence free condition cannot be preserved in discrete finite element space, and so the cancellation law (12) can not be satisfied any more. In order to solve this difficulty, we impose the discontinuous velocity on across interelement boundaries to make fulfill discrete divergence free velocity (12) automatically. We will discuss this issue at Remark 2 below. This discontinuity makes it difficult to treat nonlinear term and to apply the integration by parts, because the discontinuous solution is not included in . So we need to hire technical skills in proof of this paper.

One more remarkable discovery is in the second numerical test at the last section which is the Bénard convection problem with the same setting in [2]. In this performance, we newly find out that the number of circulations depends on the time step size . We obtain similar simulation within [2] for a relatively larger , but the behavior is changed for the small . So we conclude that the numerical result of the Bénard convection problem in [2] is not an eventual solution and a more smaller time marching step is required to get the desired simulation.

We will denote as the time marching size. Also we will use as the difference of 2 consecutive functions for example, for any sequence function ,

In order to introduce the finite element discretization, we need further notations. Let be the Sobolev space with derivatives in , set and , where or , and denote by the subspace of of functions with vanishing mean value. We indicate with the norm in and with the inner product in . Let be a shape-regular quasiuniform partition of of mesh size into closed elements [8–10]. The vector and scalar finite element spaces are where , , and are spaces of polynomials with degree bounded uniformly with respect to [9, 10]. We stress that the space is composed of continuous functions to ensure the crucial equality Using the following discrete counterpart of the form we now introduce the fully discrete SGUM to solve the evolution Boussinesq equations (1).

*Algorithm 1 (the fully discrete stabilized Gauge-Uzawa method). *Compute , , and via any method satisfying Assumption 4 below, and set and . Repeat for .

*Step 1. *Set and , and then find as the solution of

*Step 2. *Find as the solution of

*Step 3. *Update and according to

*Step 4. *Update pressure by

*Step 5. *Find as the solution of

*Remark 2 (discontinuity of ). * We note that is a discontinuous function across inter-element boundaries and that, in light of (7) and (8), is discrete divergence free in the sense that

We now summarize the results of this paper along with organization. We introduce appropriate Assumptions 1–5 in Section 2 and introduce well-known lemmas. In Section 3, we prove the stability result.

Theorem 3 (stability). *The SGUM is unconditionally stable in the sense that, for all , the following a priori bound holds:
*

We then will carry out the following optimal error estimates through several lemmas in Section 4.

Theorem 4 (error estimates). * Suppose the exact solution of (1) is smooth enough and . If Assumptions 1, 3, 4, and 5 mentioned later hold, then the errors of Algorithm 1 will be bound by
**
Moreover, if Assumption 2 also hold, then one has
*

We note that the condition in Theorem 4 can be omitted for the linearized Boussinesq equations (see Remark 16). Finally, we perform numerical tests in Section 5 to check accuracy and physically relevant numerical simulations, the Bénard convection problem and the thermal driven cavity flow.

#### 2. Preliminaries

In this section, we introduce 5 assumptions and well known lemmas to use in proof of main theorems. We resort to a duality argument for which is the stationary Stokes system with vanishing boundary condition , as well as Poisson’s equation with boundary condition .

We now state a basic assumption about .

*Assumption 1 (regularity of and ). *The unique solutions of (16) and of (17) satisfy
We remark that the validity of Assumption 1 is known if is of class [11, 12], or if is a two-dimensional convex polygon [13] and is generally believed for convex polyhedral [12].

We impose the following properties for relations between the spaces and .

*Assumption 2 (discrete inf-sup). *There exists a constant such that

*Assumption 3 (shape regularity and quasiuniformity [8–10]). *There exists a constant such that the ratio between the diameter of an element and the diameter of the largest ball contained in is bounded uniformly by , and is comparable with the meshsize for all .

In order to launch Algorithm 1, we need to set via any first-order methods which hold the following conditions.

*Assumption 4 (the setting of the first step values). *Let be the exact solution of (1) at . The first step value satisfies

*Assumption 5 (approximability [8–10]). *For each , there exist approximations such that

The following elementary but crucial relations are derived in [14].

Lemma 5 (inverse inequality). *If denotes the Clément interpolant, then
*

Lemma 6 (div-grad relation). *If , then
*

Let now indicate the finite element solution of (16); namely,

Then we can find the well-known lemmas in [8–10].

Lemma 7 (error estimates for mixed FEM). *Let be the solutions of (16), and, be the Stokes projection defined by (24), respectively. If Assumptions 1, 2, 3, and 5 hold, then
*

*Proof. *Inequality (25) is standard [8–10]. To establish (26), we just deal with the -norm, since the other can be treated similarly. If denotes the Clément interpolant, then and
as a consequence of an inverse estimate. This completes the proof.

*Remark 8 ( stability of ). *The bound is a simple by-product of (16). To see this, we add and subtract , use the stability of in , and observe that (25) implies .

Lemma 9 (error estimates for Poisson's equation [8, 10]). *Let be the solution of Poisson's equation (17), and, be its finite element approximation
**
If Assumptions 1, 3, and 5 hold, then there exists a positive constant satisfying
*

We finally state without proof several properties of the nonlinear form . In view of (5), we have the following properties of for all : Applying Sobolev imbedding lemma yields the following useful results.

Lemma 10 (bounds on nonlinear convection [1, 12]). *Let with , and let . Then
**
In addition,
*

*Remark 11 (Treatment of convection term). *As we mention in Remark 2, is discontinuous across inter-element boundaries and so . Thus, we can't directly apply (32) anymore to treat the convection term. To solve this difficulty, we apply (34) together with (26) and inverse inequality Lemma 5.

We will use the following algebraic identities frequently to treat time derivative terms.

Lemma 12 (inner product of time derivative terms). * For any sequence , one has
*

#### 3. Proof of Stability

We show that the SGUM is unconditionally stable via a standard energy method in this section. We start to prove stability with rewriting the momentum equation (6) by using (8) and (10) as follows: We now choose and apply (35) to obtain where We give thanks to (30) for eliminating the convection term. In light of , (12) and (36) yield Before we estimate , we evaluate an inequality via choosing in (9) to get Lemma 6 derives , and so (9) and (37) lead us to Clearly, we have We now attack with which comes from (8) and (12). Then we arrive at In conjunction with which comes from (30), if we choose in (11) and use (35), then we obtain Inserting back into (39), adding with (46), and then summing over from to lead to (13) by help of discrete Gronwall inequality and the equality . So we finish the proof of Theorem 3.

#### 4. Error Estimates

We prove Theorem 4 which is error estimates for SGUM of Algorithm 1. This proof is carried out through several lemmas. We start to prove this theorem with defining to be the Stokes projection of the true solution at time . It means that is the solution of We also define as the solution of And we denote notations , From Lemma 7, we can deduce In conjunction with the definition of in (26), we can derive We now carry out error evaluation by comparing (65) with (6)–(10) and then by comparing (82) and (11). We derive strong estimates of order , and this result is instrumental in proving weak estimates of order for the errors Then, in conjunction with (50), we can readily get the same accuracy for the errors Additionally, we denote We readily obtain the following properties: Moreover, from (12), Whence we deduce crucial orthogonality properties: We also point out that, owing to Lemma 6, defined in (9) satisfies In conjunction with , Assumption 5 leads to

We now estimate the first-order accuracy for velocity and temperature in Lemma 13 and then the 2nd-order accuracy for time derivative of velocity and temperature in Lemma 15. The result of Lemma 13 is instrumental to treat the convection term in proof of Lemma 15. We will use Lemmas 13 and 15 to prove optimal error accuracy in Lemma 17. Finally, we will prove pressure error estimate in Lemma 18.

Lemma 13 (reduced rate of convergence for velocity and temperature). *Suppose the exact solution of (1) is smooth enough. If Assumptions 3–5 hold, then the velocity and temperature error functions satisfy
*

*Proof. *We resort to the Taylor theorem to write (1) as follows:
where
are the truncation errors. In conjunction with the definition of the Stokes projection , we readily get a weak formulation for (62), :
We replace the pressure term in (6) with (10) and then subtract it from (65) to obtain
Choosing in (66) and using (58) and (35) yield
where
We now estimate all the terms from to , respectively. To tackle , we first add and subtract to obtain
Because of , which comes from (30), the second term of can be replaced by
If we apply Lemma 10, then we can readily obtain
Since we have according to (52), we arrive at
In light of (51) and (58), becomes
In order to estimate and , we note that the cancellation law (12) gives . Then yields
In conjunction with the definition , can be evaluated by
If we now apply inequality , then we can get
So we arrive at
In light of and (37), (9) and (59) yield
Also we readily get
The Hölder inequality and (60) yield
Inserting the above estimates into (67) gives
On the other hand, the definition (48) of leads to a weak formula of (63) as, for all ,
We now subtract (11) from (82) to derive
Choosing yields
where
To estimate , we note which comes from (30), then and yield
We use which is (60) to attack :
Inserting the above estimates into (84) gives
Adding 2 inequalities (81) and (88) and summing over from to imply
Because of , we have , and thus, Assumption 4 yields and the first 4 terms in the right-hand side can be bound by Assumption 4 and properties and which are directly deduced from the conditions in Algorithm 1. Moreover, the remaining terms can be treated by the discrete Gronwall lemma. Finally, in conjunction with (58), we conclude the desired result and complete this proof.

*Remark 14 (optimal estimation). *In order to get optimal accuracy, we must get rid of the terms of and by applying duality argument in Lemma 17. To do this, we first evaluate the errors for time derivative of velocity and temperature in Lemma 15. Thus, we need to evaluate optimal initial errors for the case , and so we have to recompute again (77). We start to rewrite as
In light of Assumption 4, we arrive at

We now start to estimate errors for time derivative of velocity.

Lemma 15 (error estimate for time derivative of velocity and temperature). *Suppose the exact solution of (1) is smooth enough and . If Assumptions 3–5 hold, then the time derivative velocity and temperature error functions satisfy
*

*Remark 16 (the condition ). *The assumption requires to control convection terms which are used at only (100) and (117), so we can omit this condition for the linearized Boussinesq equations. However, this condition cannot be removed for nonlinear equation case, because (99) must be bounded by .

*Proof. *We subtract two consecutive formulas of (66) and impose to obtain
where
We now estimate each term from to separately. The convection term can be rewritten as follows: