Computational & Applied Mathematics
versión On-line ISSN 1807-0302
Comput. Appl. Math. v.22 n.1 São Carlos 2003
Finite element approximation of bipolar viscous fluids
H. Manouzi^{I}; A. Brahmi^{II}; M. Farhloul^{II}
^{I}Département de Mathématiques et de Statistique Université Laval, Québec, G1K7P4, Canada
^{II}Département de Mathématiques et de Statistique Université de Moncton, Moncton, NB, E1A 3E9, Canada
E-mail: hm@mat.ulaval.ca
ABSTRACT
A bipolar viscous fluid model is assumed to regularise the solution of Newtonian and quasi-Newtonian flows. In this article, a mixed finite element approximation of the bipolar viscous fluids is proposed. In this approximation the velocity of the fluid together with its laplacian are the most relevant unknowns. An existence and uniqueness results are proved. A mixed finite element approximation is derived and numerical results are presented.
Mathematical subject classification: 65N30, 76M10.
Key words: finite element, bipolar fluids, Navier-Stokes, mixed finite element, PDE.
1 Introduction
Mathematical model for fluid motion play an important role in theoretical and computational studies in the aeronautics, meteorological, plasma physics, etc. The Navier-Stokes equations are generally accepted as providing an accurate model for the incompressible motion of viscous fluids in many practical situations. In the classical Navier-Stokes equations theory, the viscosity of the fluid is modeled by the dependence of the stress on the first gradient of velocity. It is well known that the corresponding mathematical theory contains a number of unsolved problems. One of the fundamental open problems for the Navier-Stokes equations in 3-dimensional space is that of the uniqueness of the solution of the Cauchy-Dirichlet problem, which is not guaranteed in the functional classes in which global existence holds. These problems have led many authors that a stronger mechanism of dissipation and viscosity, namely the dependence of the stress on the higher gradients of velocity, must occur in the flows of viscous fluids. In materials in which the higher gradients of velocity and the higher gradients of deformation influence the response, the rate of work of internal forces cannot be expected to be only the product of the usual second order stress tensor with the first gradient of velocity. Instead of this, a more general expression must be assumed containing additionally the sum of products of higher order multipolar stress tensors with the higher gradients of velocity. Otherwise such materials cannot be compatible with the Clausius-Duhem inequality. In [13] a thermodynamics theory of constitutive equations of multipolar viscous fluids within the framework of the theory of Green and Rivlin [9] has been developed. This paper is organized as follows. In section 2 we formulate the problem, we introduce some notations and we give some preliminary results. In section 3 we formulate the bipolar Navier-Stokes problem and we prove the existence and the unicity of the solution. In section 4 we study a mixed variational formulation and its approximation. Finally in section 5 we give some numerical results.
2 Problem formulation and preliminary results
2.1 Problem formulation
We will use the standard notation for gradient, divergence and curl operators, i.e. for a scalar function p and a vector function v, we write Ñq, div v and curl v for the gradient, divergence and the curl respectively. We further recall that the vector-valued Laplace operator, D, is defined by
Dv = (Dv_{1}, Dv_{2}), v = (v_{1}, v_{2})
and we have the vector identity
Next we let W be an open bounded domain of ^{2}, with smooth boundary G. We consider the stationary flow of an incompressible nonlinear bipolar fluid. This problem is modeled by the following system of equations:
where u is the velocity field, s the stress tensor and f represents the external forces. The stress tensor s = (s_{ij}) has the form (see [4], [8], [9],[10], [11], [13])
where
and r, n, l and m are positive constants with r > 1. Since div u = 0, the divergence of the stress tensor takes the form
divs = -Ñp + 2n div (F _{r}((u)) (u) ) - mD^{2}u.
Correspondingly, the bipolar Navier-Stokes equations become
By choosing the appropriate constants r, m, n and l we find these classical models of fluids:
- Carreau's law model (see [15]): m = 0.
- Power law model (see [15]): m = 0, l = 0.
- Classical Navier-Stokes: r = 2, m = 0.
- Regularized Navier-Stokes equations (see [14]): r = 2 , m ¹ 0.
These laws occur in many mathematical models of physical processes of polymer, of visco-elastic and visco-plastic fluids or of fluids with large stresses (See [10]). We will study the system of equations (4) with the boundary conditions
where t is the unit tangent to G.
2.2 Preliminary results
Let H^{m}(W) and H^{m}(G) denote the usual Sobolev spaces of order m defined on W and G, respectively; also H^{m}(W) and H^{m}(G) denote their vector-valued countreparts. The inner product on H^{m}(W) or H^{m}(W) is defined by (., .)_{m}. The space H^{m}(W) is equipped with the norm
where
For vector valued functions the above norm extends naturally as follows: for v = (v_{1}, v_{2})
We also define the Hilbert spaces
H (div, W) = {v | v Î L^{2}(W), div v Î L^{2}(W)}
H (curl, W) = {v | v Î L^{2}(W), curl v Î L^{2}(W)}
These spaces are equipped with their natural hilbertian norms.
Let (W) denote the space of infinitely differentiable functions having compact support and let H(W) be the closure of (W) with respect to the norm . _{m,W} and H^{–}^{m}(W) the dual space of H(W). The duality pairing between H(G), H(G) and their dual spaces H^{–}(G), H^{–}(G) is denoted by á., .ñ _{–}_{, }_{, G}.
We also have the following trace results (see, e. g., [6], [16])
Theorem 2.1. Let W be a bounded open set of ^{2} with a Lipschitz continuous boundary G. We have the following properties:
1. There exists a unique linear continuous map g_{0} : H^{1}(W) H(G) such that
g_{0}v = v_{|G}, "v Î H^{1}(W).
2. There exists a unique linear continuous map g_{n} : H(div, W) H^{–}(G) such that
g_{n} v = (v · n)_{|G}, "v Î H(div, W)
where n denotes the unit outward normal to G.
3. There exists a unique linear continuous map g_{t}: H(curl, W) H^{–}(G) such that
g_{t} v = (v · t)_{|G}, "v Î H(curl, W).
We use the convention H^{0}(W) = L^{2}(W), H^{0}(W) = L^{2}(W), (., .)_{0} = (., .).
Also we introduce some function spaces which are related to the study of our problem:
For details concerning Sobolev spaces consult [1], [5], [6].
Lemma 2.1. Let u Î L^{2}(W) such that Du Î , then we have the following Green Formula:
Proof. The proof of this equality follows from the classical Green formula:
and from the fact that curl curl v = –Dv and curl curl Du = –D^{2}u for all u and v satisfying the conditions of lemma 2.1.
The following preliminary results are classical and will be assumed ( see, e.g., [5], [7]).
Lemma 2.2. Let W be bounded with a smooth boundary G. Then the mapping
v Dv_{0,W}
defines a norm on X which is equivalent to the norm induced by H^{2}(W), i.e.
$a > 0, Dv_{0,W} > av_{X}.
From lemma 2.2 we see that in X the norm u_{2,W} induced by H^{2}(W) is equivalent to the norm Du_{0,W}. We then denote
as a norm in X derived from the inner product (u,v)_{X} = (Du, Dv)_{0,W}.
Lemma 2.3. Let W be bounded with a smooth boundary G. Then there exists a positive constant C = C(W) such that
3 The bipolar Navier-Stokes problem
3.1 Variational formulation
A variational formulation of the boundary value problem (4) and (5) reads as follows:
Remark 3.1. The term
is well defined for every (u,v) Î H^{2}(W) × H^{2}(W).
We associate with the problem (7) the following problem
Theorem 3.1. We have the following results:
1. If (u, p) is a smooth solution of (4)-(5) then u is a solution of (8).
2. If u is a smooth solution of (8) then there exists a function p Î M such that (u, p) is a solution of (4)-(5).
Proof. 1. Let v Î V, if (u, p) Î X × M is a solution of (4)-(5) with curl Du Î H^{1}(W) then we have
We derive from the Green formula (lemma 2.1) that
Since v Î V, div u = 0 and Du · t = 0, we have u Î V and
This proves the first part of the Theorem.
2. Now, let u Î V be a solution of (8) with curl Du Î H^{1}(W) then, we have
Let us introduce the operator
Note that < L, v > = 0 "v Î and, since -(Du, Dv) = (Du, curl curl v) = (curl Du, curl v), we get L Î H^{–}^{1}(W). Then, since W is simply-connected and from De Rham's lemma, there exists one and only one function p Î M such that
L = –Ñ p.
Then, we readly derive that
Now, Owing to Green's formula (lemma 2.1), we have
Hence, it follows that
Since f Î L^{2}(W), we have
Moreover, as u Î V, we have
Now, by using the Green's formula (lemma 2.1), we have on the one hand:
and, by (8), we get
Finally, from [2], we have
g_{t}(Du) = 0 on G.
3.2 Existence and unicity of the solution
In order to show the existence and the unicity of the the solution of problem (8) we need to state some intermediate results.
Definition 3.1. Let X be a reflexive Banach space and let A be an operator from X to X¢. We say that
- The operator A is monotone if and only if
< Au – Av, u – v > > 0 "(u, v) Î X × X
- The operator A is pseudo-monotone if A satisfies (u_{n} u) and
(limp sup áAu_{n}, u_{n} – uñ < 0)
áAu, u – wñ < limp sup áAu_{n}, u_{n} – wñ
for all w Î X.
- The operator A is strongly continuous if for every sequence u_{n}with u_{n} u we have Au_{n} Au.
Lemma 3.1. Let A_{1} and A_{2} be two operators from X to X¢ with
- A_{1} is monotone
- A_{2} is strongly continuous
Then, A_{1} + A_{2} is pseudo-monotone.
Proof. See [12], [17].
Lemma 3.2. Let X be a reflexive and separable Banach space and let A be an operator from X to X¢ with
- A is continuous
- A is coercive i.e.
- A is pseudo-monotone.
Then, for every b Î X¢ the equation Au = b has at least one solution u Î X.
Proof. See [12], [17].
In the following lemma the notation (s, t) means .
Lemma 3.3. For any s, t Î ^{2×2}, there exist constants C_{1}, C_{2}, C_{3}, C_{4} such that:
1. If 1 < r < 2, then
2. If r > 2, then
Proof. See [15].
Lemma 3.4. Let u, v Î H^{1}(W) and let w Î H^{1}(W) with div w = 0 and (w · n)|_{G} = 0. Then, we have
- ((w · Ñ) \u, v) + ((w · Ñ)v, u) = 0
- ((w · Ñ)v, v) = 0
Proof. See [6], [16].
Now, we introduce the operators A_{0}, A_{1} and A_{2} defined, for all u, v Î V, by
- (A_{0}u, v) = m(Du, Dv)
- (A_{1}u, v) = ((u · Ñ)u, v)
- (A_{2}u, v) = 2n(F_{r}(e(u))e(u), e(v)).
Let the operator A = A_{0} + A_{1} + A_{2}. Then, the problem (8) can be written as follows
Proposition 3.1. The operator A is pseudo-monotone.
Proof. According to lemma 3.1, we first prove that A_{1} is strongly continuous. To this end, let u_{n} be a sequence of elements of V such that u_{n} u in V. Owing to lemma 3.4, we have, for all v Î V,
(A_{1}u_{n}, v) – (A_{1}u, v) = (((u_{n }– u) · Ñ)u_{n}, v) + ((u · Ñ)v, u – u_{n})
and then
|(A_{1}u_{n}, v) – (A_{1}u, v)|
< C(u_{n }– u_{1,W}u_{n}_{1,W}v_{1,W} + u_{1,W}v_{1,W}u – u_{n}_{1,W}).
Now, since the embedding of V in H^{1}(W) is compact, we have u_{n} u in H^{1}(W). Therefore
(A_{1}u_{n}, v) (A_{1}u, v), "v Î V.
On the other hand, owing to the lemma 3.3, we prove that A_{0} + A_{2} is monotone. Indeed, for all u, v Î X we have
(A_{2}u – A_{2}v, u – v) = (F_{r}(e(u))e(u) – F_{r}(e(v))e(v), e(u) – e(v)) > 0
and
(A_{0}u – A_{0}v, u – v) = m(D(u – v), D(u – v)) > 0.
Therefore, owing to lemma 3.1, the operator (A_{0 }+ A_{2}) + A_{1} i.e. A is pseudo-monotone.
Proposition 3.2. The operator A is continuous and coercive.
Proof. It is obvious, since A_{0}, A_{1} and A_{2} are continuous, that A is continuous. Beside that, since ((v · Ñ)v, v) = 0, "v Î V, we have, for all v Î V,
In virtue of lemma 3.3, we have (F_{r}(e(v))e(v), e(v)) > 0. Then, we derive that
Finally, by using lemma 2.2 we conclude that A is coercive.
Theorem 3.2. The problem (8) has at least one solution u.
Proof. From propositions 3.1 and 3.2 we have
- A is continuous
- A is pseudo-monotone
- A is coercive.
Then, in virtue of lemma 3.2 and since V is reflexive and separable, we have the existence, of u Î V such that Au = f, or equivalently the problem(8) has at least one solution.
Theorem 3.3. If m is sufficiently large or f ''sufficiently small" so that
m^{2} > C(W)f_{X¢}
then problem (8) has exactly one and only one solution.
Proof. We assume that there are two solutions u_{1} and u_{2} in V of (8) then we have
(Au_{1} – Au_{2}, v) = 0 "v Î V.
Then, for all v Î V, we have
m(Du_{1} – Du_{2}, Dv) + 2n(F_{r}(e(u_{1}))e(u_{1}) – (F_{r}(e(u_{2}))e(u_{2}),
e(v)) + ((u_{1} · Ñ)u_{1}, v) – ((u_{2} · Ñ)u_{2}, v) = 0.
By taking v = u_{1} – u_{2}, we get
Since, owing to lemma 3.3,
2n(F_{r}(e(u_{1}))e(u_{1}) – F_{r}(e(u_{2}))e(u_{2}), e(u_{1}) – e(u_{2})) > 0,
we readily derive that
Now, for some positive constant C_{1} depending on W, we have
By using lemma 2.2, we get
On the other hand, since u_{1} is a solution of (8) and by taking v = u_{1} in (8), we get
Therefore,
We easily derive that u_{1} satisfies
and then
As a consequence of (25) and (26), we obtain
Therefore
Finally, we readily derive the uniqueness of the solution when the condition
holds.
4 The mixed finite element formulation
In order to approximate numerically problem (7) we introduce the mixed formulation which allows us to reduce the order of problem (7). That makes the numerical approximation easier. Let us take w = –Du as a new variable. Then, since div (Du) = 0 and owing to the identity (1), we may recast problem (4) as second order system for u and w:
In order to study the last problem we introduce the following space
Q = H_{0}(curl, W) = {v Î H(curl, W) and (v · t)|_{G} = 0}.
A variational formulation of the mixed problem reads as follows
Theorem 4.1. 1) If (u, p) Î X × M is a solution of (4)-(5) and w = –Du with curl w Î H^{1}(W) then (u, p, w) is a solution of problem (28).
2) If (u, p, w) is a solution of problem (28) then (u, p) is a solution of (4)-(5).
Proof. 1. Let (u, p) a solution of (4)-(5) and w = –Du with curl w Î H^{1}(W) . Then, as in the proof of theorem 3.1, we have Du · t = 0 i.e. w · t = 0. This yields w Î Q.
Beside that, for all v Î , we have
and
Now, by substitution in (29) we get
Owing to De Rham's lemma there exists a unique function p Î M such that
On the other hand by using the new variable w = –Du and the fact that div u = 0, we get for all z Î Q
Then, taking into account that z · t = 0 and u Î X, we find
(w, z)-(curl u, curl z) = 0 "z Î Q.
Finally, from div u = 0, we get
(div u, q) = 0 "q Î M.
2. Conversely, suppose that (u, p, w) is a solution of (28) such that u Î H^{2}(W) and curl Du Î H^{1}(W). First, since (div u, q) = 0, "q Î M and u Î H(W), we have div u = 0. Now, since div u = 0, we have
Thus w = –Du, and Du Î Q since w Î Q. Therefore,
Thus, by the fact that div Du = Ddiv u = 0 , Du Î Q, curl Du Î H^{1}(W) and by using the lemma 2.1, we have
(curl w, curl v) = (D^{2}u, v) "v Î .
Therefore,
Owing to De Rham's lemma and like in the proof of theorem 3.1, there exists one and only one function p Î M such that:
Finally, by gathering, we get
Problem (28) is easier to approximate numerically than problem (7) since the finite element approximation of (28) involves only the construction of finite-dimensional sub-spaces of the spaces H_{0}(curl, W) , H(W) and L^{2}(W). In this section we will only consider the discrete approximations to the mixed variational formulation (28). We assume for simplicity that the domain W is a two-dimensional polygon. Let _{h} be a given triangulation of into closed triangles, with h = max{diam (K)} where the maximumis taken over all K Î _{h}. Let P_{k}(D) denote the space of polynomials of degree less than or equal to k with respect to the set D Ì ^{2} and P_{k}(D) denotes the space of 2-vector valued functions each of whose components belongs to P_{k}(D). To approximate the velocity and the pressure, we use the Crouzeix-Raviart finite element:
where, for any K Î _{h},
Here, l _{i}(x), i = 1,2,3, denote the barycentric coordinates of the point x Î K with respect to the vertices of K. Finally we introduce the laplacian's velocity finite element space:
Now, let
Then (28) is discretized by the following problem
5 Numerical results
The performance of the numerical method described above has been tested on three examples of a comparative nature: the Poiseuille flow, the flow over a backward-facing step and the driven cavity flow.
- Poiseuille flow
At inlet of the channel, a parabolic profile for the velocity is used. Figures 1, 2 give the velocity fields for the Reynolds number Re = 1000 and for r = 1.7 and r = 1 respectively. The velocity distribution in the middle of the channel is shown in Figure 3 for different values of the parameter r = 1.8, 1.6, 1.4, 1.2, 1.1, 1.
- Backward-facing step flow
With flow in channels of constant geometry such as the Poiseuille flow, a predominant flow direction prevails. Instead, where the geometry changes abruptly, the flow separates and develops a recirculation region. Inlet value for the velocity is the same as that of Poiseuille flow. The Reynolds number is based on step height and a ratio outlet height to step height of 2. The streamline patterns are presented in Figures 4-11. The computations were performed for differents Reynolds numbers and for differents parameters r and m.
- The driven cavity flow
The depicted domain is the square cavity [0,1] × [0,1], and the boundary conditions are, u · t = 1 and u · n = 0 at y = 1, and u = 0 on the other parts of the boundary. In Figures 12-17 the streamlines contours are given for different Reynolds number and different parameters r and m.
For all the above numerical tests the value of the parameter l has been fixed to the value l = 10^{-11}.
REFERENCES
[1] R.A. Adams, Sobolev spaces, Academic Press, New York, (1975). [ Links ]
[2] C. Amrouche and V. Girault, Problèmes généralisés de Stokes, Publication du Laboratoire d'Analyse Numérique R91001, Université Pierre et Marie Curie, Paris. [ Links ]
[3] H. Bellout, R. Bloom and J. Necas, Phenomenological behavior of multipolar viscous fluids, Quarterly of Applied Mathematics, Volume IX, Number 4 (1991), 687-728. [ Links ]
[4] A. Brahmi, Simulation numérique des fluides bipolaires, Mémoire de maîtrise, Université Laval, Québec, (1994). [ Links ]
[5] P.G. Ciarlet, The Finite Element Method for Elliptic Problems, North-Holland, Amsterdam, 1978. [ Links ]
[6] V. Girault and P.A. Raviart, Finite element for Navier-Stokes equations, theory and algorihms, Springer Ser. Comput. Math., 5, Springer-Verlag, Berlin, New Tork, (1986). [ Links ]
[7] R. Glowinski and O. Pironneau, Numerical methods for the first biharmonic equation and for the two-dimensional Stokes problem, SIAM Rev., 21, 2 (1979), 167-212. [ Links ]
[8] K.K. Golovkin, New Equations Modeling the Motion of a Viscous Fluid, and their Unique Solvability, Proc. Steklov Inst. Math., 102 (1967). [ Links ]
[9] A.E. Green and R.S. Rivlin, Simple force and stress multipoles, Arch. Rat. Mech. Anal., 17 (1964), 325-353. [ Links ]
[10] O.A. Ladyzhenskaya, Modification of the Navier Stokes Equations for Large Velocity Gradients, in Boundary Value Problems of Mathematical Physics and Related Aspects of Function Theory II, Consultants Bureau, New York, (1970). [ Links ]
[11] O.A. Ladyzhenskaya, Nonstationary Navier-Stokes equations, Amer. Math. Soc. Transl. (2) 25 (1962), 151-160. [ Links ]
[12] J.L. Lions, Quelques méthodes de résolution des problèmes aux limites non linéaires, Dunod, Paris, (1969). [ Links ]
[13] J. Necas and M. Silhavy, Multipolar viscous fluids, Quart. Appl. Math., 49 (1991), 247-265. [ Links ]
[14] Y.U. Ou and S.S. Sritharan, Analysis of regularized Navier-Stokes equations I, Quart. Appl. Math., 49 (1991), 651-685. [ Links ]
[15] D. Sandri, Sur l'approximation numérique des écoulements quasi-newtoniens dont la viscosité suit la loi puissance ou la loi de Carreau, M^{2}AN 27 (1993), 131-155. [ Links ]
[16] R. Temam, Navier-Stokes equations, North-Holland, Amsterdam, (1977). [ Links ]
[17] E. Zeidler, Nonlinear functional analysis and its applications II, Springer-Verlag, (1985). [ Links ]
#551/02.