Abstract
In this work, a stationary Stokes flow with thermal effects is studied both mathematically and numerically. First, existence, uniqueness and regularity of the weak solution of the problem are established. Next, finite element approximation to the problem, based on a fixed point algorithm, is proposed. Then, an error estimate between continuous solution and discrete one is obtained. Finally, some numerical tests are presented to confirm the theoretical results.
finite element analysis; thermally coupled Stokes problem; incompressible quasi-Newtonian flow
Finite element analysis of a coupled thermally dependent viscosity flow problem
Jiang Zhu; Abimael F.D. Loula; José Karam F.; João N.C. Guerreiro
Laboratório Nacional de Computação Científica, MCT, Avenida Getúlio Vargas 333, 25651-075 Petrópolis, RJ, Brazil, E-mails: jiang@lncc.br / aloc@lncc.br / jkfi@lncc.br / joao@lncc.br
ABSTRACT
In this work, a stationary Stokes flow with thermal effects is studied both mathematically and numerically. First, existence, uniqueness and regularity of the weak solution of the problem are established. Next, finite element approximation to the problem, based on a fixed point algorithm, is proposed. Then, an error estimate between continuous solution and discrete one is obtained. Finally, some numerical tests are presented to confirm the theoretical results.
Mathematical subject classification: 35J, 35Q, 65N30, 76.
Key words: finite element analysis, thermally coupled Stokes problem, incompressible quasi-Newtonian flow.
1 Introduction
In modelling engineering problems describing incompressible quasi Newtonian flows with viscous heating we need to consider the following thermally coupled Stokes problem (see for instance [1-4], and the references therein):
where u: W ® d is the velocity, p: W ® is the pressure, q : W ® is the temperature, W is a bounded open subset of d, d = 2 or 3, G its boundary. The viscosity µ is a function of q , µ = µ ( q ). D(u) = ( Ñu + ÑuT) is the strain rate tensor, and |D(u)|2 is the second invariant of D(u).
Problems of this type have received especial attention recently. Mathematical analysis of this class of problems can be found, for example, in [2, 5]. In [5], a convergence result for an iterative method was obtained under very strong regularity hypothesis. To our knowledge, there is no general result on the numerical analysis of problem (1.1). For mathematical and numerical analyses of simpler problems consisting of nonlinear coupled systems of two scalar elliptic equations, we refer to [6-12].
Complete mathematical and numerical analyses to a coupled nonlinear system of scalar elliptic equations are presented in [11]. In the present work, we will extend the analyses presented in [11] to problem (1.1). We admit for simplicity homogeneous boundary conditions and assume that the coupling function m Î C () is bounded, i.e., there exist constants K2> K1 > 0 such that, for all x Î ,
We first establish existence, uniqueness and regularity of the weak solution of problem (1.1). Then, we apply a fixed point algorithm and propose a finite element approximation. We prove the convergence of the fixed point algorithm and derive error estimates for the discrete iterative solutions. Finally, we present some numerical results to confirm the predicted rates of convergence of the finite element approximations and to illustrate the influences of nonhomogeneous boundary conditions and of the source term f on existence and stability of solution for a two-dimensional model.
2 Variational formulation
Let Wm,r(W) denote the standard Sobolev space with its norm , for m > 0 and 1 < r < ¥ . We write Hm( W ) = Wm,2( W ) when r = 2, with the norm , and Lr( W ) = W0,r( W ) when m = 0, with the norm . ( W ) is the closure of the space ( W ) for the norm . Vector variables are, in general, denoted with bold face. We denote also Wm,r( W ) = [Wm,r( W )]d, ( W ) = [( W )]d, Hm( W ) = [Hm( W )]d, ( W ) = [( W )]d, and Lr( W ) = [Lr( W )]d.
Throughout this work, we assume that f Î L2( W ), then the variational formulation of problem (1.1) can be defined as:
where
(· ,·) denotes the duality between Lr( W )d and Lr' ( W )d, d = 1, 2, 3, r' is the dual number of r. ( W ) = {q Î Lr( W ) | òW q = 0}. Introducing the space:
we associate with (2.1) the following problem:
The classical Korn's inequality implies that the norm is equivalent to the norm in space ( W ). Condition (1.2) implies that
and
b(v, q) satisfies the inf-sup condition, i.e. there exists a constant b > 0 such that (cf. [13])
3 Existence, uniqueness and regularity
Definition 1. We denote by Sr for r Î (1, ¥ ) the class of regular subsets G in dfor which the Stokes operator maps (G) = { (G) | Ñ · v = 0 } onto W-1,r(G).
Remark 1. For r Î (1, ¥ ), a bounded C1 domain or a bounded Lipschitz domain with sufficiently small Lipschitz constant depending on d and r, is of class Sr [14].
From now on, we assume W is of class Sr for some r > 2. For 1 < s < r, we define Ms > 1 by
Remark 2. Similarly to [15] and [16], we can see that Ms < ¥ and especially M2 = 1.
Similarly to [12] and [15], we can prove
Lemma 1.For any given q , if u Î ( W ) satisfies (2.1.i, ii) (or (2.6.i)), then there exist s Î (2, min{6, r}] satisfying
and a constant Cs > 0 defined by
such that u Î ( W ) and the following estimate holds
If s defined in Lemma 1 is such that
we have
Thus, by the density of ( W ) Ç L¥ ( W ) in ( W ), problem (2.1) can be written equivalently as:
And the associated problem (2.6) can also be written equivalently as:
We now prove existence of a solution to problem (3.8). For any given x Î L2( W ), we denote by ux Î V the solution of
and define by qx Î ( W ) the solution of
By (3.10), (1.2) and (3.6), we have
where the constant C0 depends only on W , K1 and K2.
Let BR be the ball in L2( W ) defined by
Then, the map T defined by
is compact since ( W ) is compactly imbedded in L2( W ), and satisfies that T(L2( W )) Ì BR. We only need show that the map T is continuous, then the solvability of problem (3.8) comes from the Schauder Fixed Point Theorem. Furthermore, by Lemma I.2.1 in [13], there exists p Î ( W ) such that (u, p, q ) solves problem (3.7).
To show the continuity of the map T, let xj® x in L2( W ), by Lemma 1, the corresponding solutions {} of
satisfy
So, there is a subsequence denoted by {} such that
The uniqueness of the solution of (3.9) implies that (the whole sequence)
Then, noticing (3.5), we can see that
Since (3.10) is a standard elliptic problem, the continuity of the map T is well known.
Theorem 1 (Existence).If s defined in Lemma 1 satisfies (3.5), then problems (2.1) and (2.6) are equivalent to problems (3.7) and (3.8), respectively. Problem (3.7) has a solution (u, p, q ) whereas (u, q ) is a solution of problem (3.8), and the following estimates hold:
and
where
Proof. It is only needed to prove (3.20). In fact, by the Sobolev inequality, we have
where C > 0 is a constant depending only on W , K1, K2 and s.
To study the uniqueness of the problem, we need to assume that the function µ is Lipschitz continuous, i.e., there is a Lipschitz constant L, for any x 1, x 2Î , such that
Suppose (3.8) has two solutions (u1, q1) and (u2, q2), and let = u1 - u2 and = q 1 - q2. Then, by (3.8), we have " v Î V
and
By (3.24), (3.23), (3.19) and the Sobolev inequality, for
we have
where C > 0 is a constant dependent on W , K1, K2, L and s.
Let h = in (3.25), and considering the Hölder inequality, the Sobolev inequality, (3.30), (3.27) and (3.19), we have
where is a constant directly proportional to K2 and L, inversely to K1, and dependent on W and s as well.
Therefore, if
then it holds that = 0, which implies that = 0 by (3.27). With the above result, we can state:
Theorem 2 (Uniqueness). If conditions (3.23), (3.26) and (3.29) hold. Then, problem (3.8) (or (3.7)) has a unique solution.
Theorem 3 (Regularity). If
and if s defined in Lemma 1 satisfies that
then we have
and
where C > 0 is a constant dependent on W , K1, K2, L and s.
Proof. From (3.30), (3.19) and (3.20), we get
Thus, by Theorem I.5.4 in [13], we have
Then,
4 A fixed point algorithm
From the numerical point of view, it is interesting to introduce an iterative scheme to solve problems (3.7) and (3.8). The scheme proposed in this section is based on a fixed point algorithm.
For an arbitrary q0, and n = 1, 2, ... , we can get an iterative solution of problem (3.7) {(un,pn, qn)} by:
And an iterative solution of problem (3.8) {(un, qn)} can be obtained by:
Theorem 4.The solutions {(un,pn, qn)} of (4.1) and {(un, qn)} of (4.2) satisfy
where s andare same as in Theorem 1. Moreover, with the same assumptions of Theorem 3, we have
and
Proof. The proof is similar to that of Theorem 1 and 3.
Theorem 5. If problem (3.7) has a unique solution (u,p, q ), then the sequence {(un,pn, qn)} defined by (4.1) converges in ( W )×( W ) ×( W ) to (u,p, q); and the sequence {(un, qn)} defined by (4.2) converges in V×( W ) to (u, q).
Proof. We only give a proof for problem (4.1), since the proof for problem (4.2) is similar. If this Theorem is not true, then there exist some small constant e0 > 0 and an infinite subsequence of {( un, qn )}, denoted by {()}, such that
On the other hand, Theorem 4 implies that {()} Î ( W )×( W ) ×( W ). Since the space ( W )×( W ) ×( W ) is compact in ( W )×( W ) ×( W ), and noticing the fact that any limit of {(un,pn, qn)} should satisfy (3.7) and that problem (3.7) has a unique solution. Then we can get a subsequence of {()} which converges to ( u,p, q ) in ( W )×( W )×( W ), which leads a contradiction to (4.7). Hence, we complete the proof. To further study u- un, p - pn and q - qn, by (3.7) and (4.1), and using similar arguments to (3.27) and (3.28), we can deduce that
where is the same as in (3.28). The estimate for p - pn comes from (3.7), (4.1) and the inf-sup condition (2.9),
Thus, we have
Theorem 6. If condition (3.29) holds, then, the fixed point algorithms (4.1) and (4.2) work with the linear convergence rate, and the following estimates hold:
where.
5 Finite element approximation
For simplicity we assume that W is a polygonal (or polyhedral) domain discretized by a quasi uniform mesh of Ne triangles (or tetrahedrons) or convex quadrilaterals (or hexahedrons), with mesh parameter h. Let Sh be the Lagrangian finite element space of C0( W ) piecewise linear polynomials, and = Sh Ç ( W ). Let also Xh Ì H1( W ) and Qh Ì L2( W ) be two finite element spaces, and = X Ç ( W ) and = Qh Ç ( W ) such that the following hypotheses hold:
Hypothesis H1 (Approximation property of ). There exists an operator PXÎ (H2( W ); Xh) Ç (H2( W ) Ç ( W ); ) such that:
Hypothesis H2 (Approximation property of Qh). There exists an operator P QÎ (L2( W );Qh) such that:
Hypothesis H3 (Uniform inf-sup condition). For each qh Î there exists a vh Î such that:
with a constant C > 0 independent of h, qh and vh.
Remark 3. Hypothesis H3 is equivalent to the discrete inf-sup condition, i.e. (cf. [13])
The Galerkin approximation to problem (4.1) reads:
Given as an approximation of q 0, for n = 1, 2, ... , {()} can be calculated by:
Since Ñ · Î ( W ), then (5.6.ii) is equivalent to
Thus, we can define the space:
and the problem associated with (5.6) is:
To analyze problem (5.6), we introduce the standard Stokes and elliptic projections () Î × × defined by
It is well known that
Lemma 2. There exists a constant C > 0 independent of h and n such that the following estimates hold:
The following inverse properties of the finite element spaces and are useful.
Lemma 3. For any h hÎ , we have
where
and M is a constant independent of h.
Proof. See [11]
For the errors of and , we have:
Lemma 4. There exists a constant C > 0 only dependent on W , K1, K2and L such that
Proof. By (5.6.i) and (4.1.i), we have, " vh Î
Let vh = . By noticing that b(, qh) = 0, " qh Î and the fact that (cf. [13]) , we can get
Thus, Lemma 2 leads to (5.16).
Lemma 5. There exists a constant C > 0 only dependent on W , K1, K2, L and b* such that
Proof. By (5.6.ii) and (4.1.ii), we have, " vh Î
By Remark 3, Hypothesis H3 yields:
In the last inequality, we applied Lemmas 2 and 4.
Lemma 6. The following estimate
holds with C > 0 only dependent on W , K1, K2and L, and M(h) defined by (5.15).
Proof. By (5.6.iii), (5.11) and (4.1), we have
By the Hölder inequality and the Sobolev inequality,
Since a2-b2 = (a-b)2 + 2b(a-b), then R2 can be split into:
By the Hölder inequality and the inverse inequality (5.14),
To estimate R22, we have
Combining (5.19)-(5.23), and noticing that
we obtain (5.18).
Let us now make an inductive hypothesis: for sufficiently small h,
In fact, when n = 1, we can choose as the standard elliptic projection of q 0, thus, M(h) < 1 for sufficiently small h. If (5.25) holds for n-1, then, we have
where is a constant only dependent on W , K1, K2 and L. If
and for sufficiently small h such that
or
then, we can obtain
Thus, choosing h sufficiently small such that
where C depends on M, then, the inductive hypothesis (5.25) holds. Meanwhile, we have
Theorem 7. Let () be the solution of problem (5.6) and () the solution of problem (5.9). Then, for sufficiently small h satisfying (5.29), () converges to the solution (un,pn, q n) of problem (4.1), and the following error estimates hold:
where C is a constant independent of n, h and f, and( f ) < 1 is defined by (5.27).
Now, let C * = max{, } where and are defined by (3.28) and (5.26) respectively. Thus,
implies (3.29) and (5.27). Hence, we get the main result
Theorem 8. If condition (5.32) holds, then problem (3.7) has a unique solution ( u,p, q ) and ( u, q ) is also the solution of problem (3.8), the finite element solution sequence {()} of (5.6) (where {()} solves problem (5.9)) converges to ( u,p, q ) and the following estimates hold, for sufficiently small h satisfying (5.29),
where C is a constant independent of n, h and f, and M * (f) < 1 is defined by (5.32).
6 Numerical results
The iterative solution method described in the previous section has been tested in two dimensions considering two temperature dependence functions. In the first case we considered a bounded µ (s) given by
This case has been performed only to confirm the convergence estimates obtained here by prescribing homogeneous boundary conditions as adopted in the analysis. Exact solution has been assumed to be the results obtained with a refined mesh of 4096 elements, quadratic for the velocity field and linear for the pressures combined with linear temperature approximation elements. The plots of Fig. 1 confirm the first order convergence rate for the gradient of the velocities and of the temperatures as obtained in the analysis. Fig. 1 also shows the second order convergence of the velocity and the temperature fields as expected but not demonstrated here. Results for the convergence in terms of the number of iterations versus values of f are shown in Fig. 2. We note that for a critical value of f the number of iterations increases largely, but it still remains practically independent of f before the critical value.
The second example deals with an axisymmetric flow in a cylindrical tube. In this study we consider the Arrhenius function
with H the activation energy. The tube is 10 cm long with a diameter of 4 cm. It is a developing hydrodynamics and thermal problem. A Dirichlet boundary condition for the velocity has been prescribed at the entrance section, u = {1.0; 0.0} cm/s, homogeneous Neumann condition has been set at the end of the tube and non slip at the wall.
For the temperature, two cases have been considered. In the first case, a 500K has been prescribed at the entrance and along the cylinder wall. Changes in the temperature field are due to the coupling viscous dissipation, with higher temperatures in the core of the tube. Of course, the solution of the uncoupled analogue of this problem is a uniform temperature of 500K all over the domain. The corresponding temperature profiles are shown in Fig. 3.
The second case corresponds to what is known as extended or modified Graetz problem [17], after L. Graetz (1885) studies [18] on developed hydrodynamics-developing thermal fields for fluids passing a flat plate. The model problem consists of a fluid at 700K, entering the same tube with the same velocity boundary conditions of the first case studied, cooled by a prescribed temperature of 500K at the tube wall. Changes in the temperature field are determined for the uncoupled and for the coupled situations. In this problem the coupled effects do not alter the velocity field but they are strong enough to affect the temperature field. The activation energy has been set to 5q min. The classical temperature profiles (isothermal lines) obtained for the uncoupled case are shown in Fig. 4. For the coupled situation, the isothermal lines obtained are depicted in Fig. 5, exhibiting a hotter region along the core of the duct, as expected.
Concerning convergence of the iterative method, 120 iterations have been required for the convergence of the coupled modified Graetz problem while in the coupled problem considered in the first case convergence has been achieved with only 45 iterations.
7 Conclusions
In this paper, we gave the complete mathematical analyses, such as existence, uniqueness and regularity of the weak solution of problem (1.1). The uniqueness is conditional and dependent on the source term f which should be small. We applied a fixed point algorithm to solve the nonlinear problem and proposed a finite element approximation. We proved the convergence of the fixed point algorithm and derive error estimates for the discrete iterative solutions. We got condition (5.32) for f which guarantees not only the uniqueness of the weak solution but also the convergence. Finally, we presented numerical implementations for a two-dimensional model, using Taylor-Hood elements for the velocity-pressure and bilinear elements for the temperature, to confirm the predicted rates of convergence of the finite element approximations and to illustrate the influences of nonhomogeneous boundary conditions and of the source term f on existence and stability of solution.
REFERENCES
[1] P. Avenas, J.F. Agassant and J.Ph. Sergent. La mise en forme des matières plastiques. Lavoisier, 1982.
[2] J. Baranger and A. Mikelic, Stationary solutions to a quasi-newtonian flow with viscous heating. M3AS: Math. Model. Meth. Appl. Sci., 5 (1995), 725-738.
[3] J. Bass, Thermoelasticity. McGraw-Hill Encyclopedia of Physics, Parker SP (ed). McGraw-Hill, New York, 1982.
[4] J. Karam F., Thermally Coupled Flow of Non-Newtonian Fluids in Conjugated Problems. Doctoral Thesis, COPPE, Federal University of Rio de Janeiro, 1996 (in Portuguese).
[5] S. Wardi, A convergence result for an iterative method for the equations of a stationary quasi-newtonian flow with temperature dependent viscosity. M2AN: Modél. Math. Anal. Numér., 32 (1998), 391-404.
[6] G. Cimatti, A bound for the temperature in the thermistor problem. IMA J. Appl. Math., 40 (1988), 15-22.
[7] G. Cimatti, Remarks on existence and uniqueness for the thermistor problem under mixed boundary conditions. Quart. Appl. Math., 47 (1989), 117-121.
[8] G. Cimatti and G. Prodi, Existence results for a nonlinear elliptic system modelling a temperature dependent electrical resistor. Ann. Mat. Pura Appl., 152 (1989), 227-236.
[9] C.M. Elliott and S. Larsson, A finite element model for the time-dependent Joule heating problem. Math. Comp., 64 (1995), 1433-1453.
[10] S.D. Howison, J.F. Rodrigues and M. Shillor, Stationary solutions to the thermistor problem. J. Math. Anal. Appl., 174 (1993), 573-588.
[11] A.F.D. Loula and J. Zhu, Finite element analysis of a coupled nonlinear system. Comput. Appl. Math., 20 (2001), 321-339.
[12] J. Zhu and A.F.D. Loula, Mixed finite element analysis of a thermally nonlinear coupled problem. Numer. Methods Partial Differential Eq., 22 (2006), 180-196.
[13] V. Girault and P.A. Raviart, Finite Element Methods for Navier-Stokes Equations. Springer-Verlag: Berlin, 1986.
[14] G.P. Galdi, C.G. Simader and H. Sohr, On the Stokes problem in Lipschitz domains. Ann. Mat. pura appl., CLXVII (1994), 147-163.
[15] K. Groger, A W1,p - estimate for solutions to mixed boundary value problems for second order elliptic differential equations. Math. Ann., 283 (1989), 679-687.
[16] N.G. Meyers, An Lp-estimate for the gradient of solutions of second order elliptic divergence equations. Ann. Sc. Norm. Sup. Pisa, 17 (1963), 189-206.
[17] R.K. Shah and A.L. London, Laminar Flow Forced Convection in Ducts. Academic Press, NY, 1978.
[18] L. Graetz, Uber die Wärmeleitungsfähigkeit von Flüssigkeiten (On Thermal Conductivity of Liquids), Part 1: Ann. Phys. Chem., 18 (1883), 79-94; Part 2: Ann. Phys. Chem., 25 (1885), 337-357.
Received: 01/XI/05. Accepted: 13/XII/05.
This work was supported partially by Brazilian Research Council (CNPq) and Research Foundation of Rio de Janeiro State (FAPERJ).
#632/05.
- [1] P. Avenas, J.F. Agassant and J.Ph. Sergent. La mise en forme des matičres plastiques Lavoisier, 1982.
- [2] J. Baranger and A. Mikelic, Stationary solutions to a quasi-newtonian flow with viscous heating. M3AS: Math. Model. Meth. Appl. Sci., 5 (1995), 725-738.
- [3] J. Bass, Thermoelasticity McGraw-Hill Encyclopedia of Physics, Parker SP (ed). McGraw-Hill, New York, 1982.
- [4] J. Karam F., Thermally Coupled Flow of Non-Newtonian Fluids in Conjugated Problems Doctoral Thesis, COPPE, Federal University of Rio de Janeiro, 1996 (in Portuguese).
- [5] S. Wardi, A convergence result for an iterative method for the equations of a stationary quasi-newtonian flow with temperature dependent viscosity. M2AN: Modél. Math. Anal. Numér., 32 (1998), 391-404.
- [6] G. Cimatti, A bound for the temperature in the thermistor problem. IMA J. Appl. Math., 40 (1988), 15-22.
- [7] G. Cimatti, Remarks on existence and uniqueness for the thermistor problem under mixed boundary conditions. Quart. Appl. Math., 47 (1989), 117-121.
- [8] G. Cimatti and G. Prodi, Existence results for a nonlinear elliptic system modelling a temperature dependent electrical resistor. Ann. Mat. Pura Appl., 152 (1989), 227-236.
- [9] C.M. Elliott and S. Larsson, A finite element model for the time-dependent Joule heating problem. Math. Comp., 64 (1995), 1433-1453.
- [10] S.D. Howison, J.F. Rodrigues and M. Shillor, Stationary solutions to the thermistor problem. J. Math. Anal. Appl., 174 (1993), 573-588.
- [11] A.F.D. Loula and J. Zhu, Finite element analysis of a coupled nonlinear system. Comput. Appl. Math., 20 (2001), 321-339.
- [12] J. Zhu and A.F.D. Loula, Mixed finite element analysis of a thermally nonlinear coupled problem. Numer. Methods Partial Differential Eq., 22 (2006), 180-196.
- [13] V. Girault and P.A. Raviart, Finite Element Methods for Navier-Stokes Equations Springer-Verlag: Berlin, 1986.
- [14] G.P. Galdi, C.G. Simader and H. Sohr, On the Stokes problem in Lipschitz domains. Ann. Mat. pura appl., CLXVII (1994), 147-163.
- [15] K. Groger, A W1,p - estimate for solutions to mixed boundary value problems for second order elliptic differential equations. Math. Ann., 283 (1989), 679-687.
- [16] N.G. Meyers, An Lp-estimate for the gradient of solutions of second order elliptic divergence equations. Ann. Sc. Norm. Sup. Pisa, 17 (1963), 189-206.
- [17] R.K. Shah and A.L. London, Laminar Flow Forced Convection in Ducts. Academic Press, NY, 1978.
- [18] L. Graetz, Uber die Wärmeleitungsfähigkeit von Flüssigkeiten (On Thermal Conductivity of Liquids), Part 1: Ann. Phys. Chem., 18 (1883), 79-94;
- Part 2: Ann. Phys. Chem., 25 (1885), 337-357.
Publication Dates
-
Publication in this collection
10 May 2007 -
Date of issue
2007
History
-
Accepted
13 Dec 2005 -
Received
01 Nov 2005