Finite element analysis of a coupled thermally dependent viscosity flow problem

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: [e-mail] / [e-mail] / [e-mail] / [e-mail]

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 + Ñu^T) is the strain rate tensor, and |D(u)|² 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 K₂> K₁ > 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 W^m,r(W) denote the standard Sobolev space with its norm , for m > 0 and 1 < r < ¥ . We write H^m( W ) = W^m,2( W ) when r = 2, with the norm , and L^r( W ) = W^0,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 W^m,r( W ) = [W^m,r( W )]^d, ( W ) = [( W )]^d, H^m( W ) = [H^m( W )]^d, ( W ) = [( W )]^d, and L^r( W ) = [L^r( W )]^d.

Throughout this work, we assume that f Î L²( W ), then the variational formulation of problem (1.1) can be defined as:

where

(· ,·) denotes the duality between L^r( W )^d and L^r' ( W )^d, d = 1, 2, 3, r' is the dual number of r. ( W ) = {q Î L^r( 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 S_r 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 C¹ domain or a bounded Lipschitz domain with sufficiently small Lipschitz constant depending on d and r, is of class S_r [14].

From now on, we assume W is of class S_r for some r > 2. For 1 < s < r, we define M_s > 1 by

Remark 2. Similarly to [15] and [16], we can see that M_s < ¥ and especially M₂ = 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 C_s > 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 Î L²( W ), we denote by u_x Î V the solution of

and define by q_x Î ( W ) the solution of

By (3.10), (1.2) and (3.6), we have

where the constant C₀ depends only on W , K₁ and K₂.

Let B_R be the ball in L²( W ) defined by

Then, the map T defined by

is compact since ( W ) is compactly imbedded in L²( W ), and satisfies that T(L²( W )) Ì B_R. 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 x_j® x in L²( 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 , K₁, K₂ 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 ₁, x ₂Î , such that

Suppose (3.8) has two solutions (u₁, q₁) and (u₂, q₂), and let = u₁ - u₂ and = q ₁ - q₂. 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 , K₁, K₂, 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 K₂ and L, inversely to K₁, 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 , K₁, K₂, 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 q⁰, and n = 1, 2, ... , we can get an iterative solution of problem (3.7) {(uⁿ,pⁿ, qⁿ)} by:

And an iterative solution of problem (3.8) {(uⁿ, qⁿ)} can be obtained by:

Theorem 4.The solutions {(uⁿ,pⁿ, qⁿ)} of (4.1) and {(uⁿ, qⁿ)} 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 {(uⁿ,pⁿ, qⁿ)} defined by (4.1) converges in ( W )×( W ) ×( W ) to (u,p, q); and the sequence {(uⁿ, qⁿ)} 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 e₀ > 0 and an infinite subsequence of {( uⁿ, qⁿ )}, 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 {(uⁿ,pⁿ, qⁿ)} 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- uⁿ, p - pⁿ and q - qⁿ, 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 - pⁿ 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 S^h be the Lagrangian finite element space of C⁰( W ) piecewise linear polynomials, and = S^h Ç ( W ). Let also X^h Ì H¹( W ) and Q^h Ì L²( W ) be two finite element spaces, and = X Ç ( W ) and = Q^h Ç ( W ) such that the following hypotheses hold:

Hypothesis H1 (Approximation property of ). There exists an operator P_XÎ (H²( W ); X^h) Ç (H²( W ) Ç ( W ); ) such that:

Hypothesis H2 (Approximation property of Q^h). There exists an operator P _QÎ (L²( W );Q^h) such that:

Hypothesis H3 (Uniform inf-sup condition). For each q_h Î there exists a v_h Î such that:

with a constant C > 0 independent of h, q_h and v_h.

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 ⁰, 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 , K₁, K₂and L such that

Proof. By (5.6.i) and (4.1.i), we have, " v_h Î

Let v_h = . By noticing that b(, q_h) = 0, " q_h Î 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 , K₁, K₂, L and b^* such that

Proof. By (5.6.ii) and (4.1.ii), we have, " v_h Î

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 , K₁, K₂and 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 a²-b² = (a-b)² + 2b(a-b), then R₂ can be split into:

By the Hölder inequality and the inverse inequality (5.14),

To estimate R₂₂, 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 ⁰, 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 , K₁, K₂ 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 (uⁿ,pⁿ, q ⁿ) 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 W^1,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 L^p-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 W^1,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 L^p-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