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:
In estimating convection terms, we will use Lemma 10 frequently without notice. We recall to obtain
The result in Lemma 13, , is essential to treat the next 2 convection terms. Invoking (33), we have
We note which comes from (30). Then we obtain
To attack , we first note Lemma 5 which is, for any , . If we apply
which is the result of Lemma 13, then we can conclude, in light of (34),
We now estimate using :
In light of Hölder inequality, (51) yields
Integral by parts leads to
In order to tackle , marking use of , we readily get
If we now apply inequality , then we can have
So we arrive at
Invoking (9), (37) and (60) lead to
Finally, the term becomes
For the new term , we argue as in Lemma 13 to arrive at
We now insert the above estimates into (93) to obtain
To evaluate errors of , we subtract two consecutive formulas (83) and choose by as follows:
where
We treat the term first by rewriting it as follows:
In estimating convection terms, we will use Lemma 10 frequently without notice. We recall and which is the result in Lemma 13 to obtain
In conjunction with which comes from (30), we rewrite as
Because we can readily get via Lemmas 5 and 13, we conclude, in conjunction with (60),
Before we attack , we note that the assumption is required to apply . We now imply and (51) to get
In light of Hölder inequality, becomes
Inserting the above estimates into (111) yields
Adding 2 equations (110) and (119) and then summing up from to lead up to (92) and complete the proof.
We now estimate optimal accuracy for velocity and temperature.
Lemma 17 (full rate of convergence for velocity and temperature). One denotes that and are the solutions of (16) and (24) with , respectively. And let and be solutions of (17) and (28) with . Let the exact solution of (1) be smooth enough and . If Assumptions 1 and 3–5 hold, then one has
Proof. We choose in (66) to obtain
where
We now estimate all the terms from to , respectively. The convection term can be rewritten as follows:
and we denote by , for , the four terms in the right-hand side. To estimate convection terms, we will use frequently Lemma 10 without notice. Using , we can readily get
Because and on boundary, we can use (32) and so we get
In light of , we can obtain
and so we can conclude
If we apply which can be derived from Lemmas 5 and 13, then we can get, by the help of Lemma 10 and Assumption 1,
In conjunction with (51), we can have
The definition and Assumption 1 give us
If we apply (127) again, then we arrive at
On the other hand, the truncation error term becomes
For the new term , we observe ; whence
Invoking and inserting the above estimates from and into (121) lead us to
On the other hand, we choose in (83) to obtain
where
In order to estimate , we note first and which come from Lemmas 5 and 13, respectively. Then Lemma 10 and the assumption yield
Since , we can readily get
Inserting the above estimates into (136) yields
In conjunction with the discrete Gronwall inequality, adding (135) and (140) and summing over from to lead to (120).
We now estimate the pressure error in . This hinges on the error estimates for the time derivative of velocity and temperature of Lemma 15.
Lemma 18 (pressure error estimate). Let the exact solution of (1) be smooth enough and . If Assumptions 1–5 hold, then one has
Proof. We first recall again inf-sup condition in Assumption 2. Consequently, it suffices to estimate in terms of . In conjunction with (10), we can rewrite (66) as We now proceed to estimate each term from to separately. We readily obtain Terms and can be dealt with thanks to the aid of Lemma 10 and as follows: In light of from Lemma 13, we can have Integrate by parts and Hölder inequality yield On the other hand, we have The new term can be bound by the Hölder inequality as Inserting the estimates from to back into (142) and employing the discrete inf-sup condition in Assumption 2, we obtain If we now square it, multiply it by , and sum over from to , then Lemmas 13, 15, and 17 derive (141).
5. Numerical Experiments
We finally document 3 computational performances of SGUM. The first is to check accuracy, and then the next 2 examples are physically relevant numerical simulations, the Bénard convection problem and the thermal driven cavity flow. We perform the last 2 examples under the same set within [2], but we conclude with different numerical simulation for the second test, the Bénard convection problem, from that of [2]. We impose Taylor-Hood () in all 3 experiments.
Example 19 (mesh analysis). In this first experiment, we choose square domain and impose forcing term the exact solution to become Table 1 is error decay with and . We conclude that the numerical accuracy of SGUM is optimal and consists with the result of Theorem 4.
Example 20 (Bénard convection). In order to explore the applicability of the SGUM, we consider the Bénard convection on the domain with forcing and . Figure 1 displays the initial and boundary conditions for velocity and temperature , as already studied in [2]. Figures 2–4 are simulations at with the nondimensional parameters , and . Figure 2 is the result for the case which is the same condition in [2], and so it displays similar behavior within [2] including 6 circulations in the velocity stream line. However, Figures 3 and 4, the higher resolution simulations with , and , display 8 circulations in the stream line. So we conclude that the high resolution result is correct simulation and thus Figure 2 and the result in [2] are not eventual simulation.
Example 21 (thermal driven cavity flow). We consider the thermal driven cavity flow in an enclosed square , as already studied in several papers [2, 6, 7]. The experiment is carried out with the same setting as in the work of Gresho et al. [6], which is shown in Figure 5. Figure 6 displays the evolution from rest to steady state .
(a)
(b)
(c)
(d)