Finite element approximation of bipolar viscous fluids

Finite element approximation of bipolar viscous fluids

H. Manouzi^I; A. Brahmi^II; M. Farhloul^II

^IDépartement de Mathématiques et de Statistique Université Laval, Québec, G1K7P4, Canada

^IIDépartement de Mathématiques et de Statistique Université de Moncton, Moncton, NB, E1A 3E9, Canada

E-mail: [e-mail]

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₁, Dv₂), v = (v₁, v₂)

and we have the vector identity

Next we let W be an open bounded domain of

², 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²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₁, v₂)

We also define the Hilbert spaces

H (div, W) = {v | v Î L²(W), div v Î L²(W)}

H (curl, W) = {v | v Î L²(W), curl v Î L²(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

² with a Lipschitz continuous boundary G. We have the following properties:

1. There exists a unique linear continuous map g₀ : H¹(W) H(G) such that

g₀v = v_|G, "v Î H¹(W).

2. There exists a unique linear continuous map g_n : H(div, W) H^–(G) such that

g_nv = (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⁰(W) = L²(W), H⁰(W) = L²(W), (., .)₀ = (., .).

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²(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²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 byH²(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²(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²(W) × H²(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¹(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¹(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²(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

  _nwith u

  _n

  u we have Au

  _n

  Au.

Lemma 3.1. Let A₁ and A₂ be two operators from X to X¢ with

• A

  ₁

  is monotone

• A

  ₂

  is strongly continuous

Then, A₁ + A₂ 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₁, C₂, C₃, C₄ such that:

1. If 1 < r < 2, then

2. If r > 2, then

Proof. See [15].

Lemma 3.4. Let u, v Î H¹(W) and let w Î H¹(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₀, A₁ and A₂ defined, for all u, v Î V, by

• (A

  ₀

  u,

  v) = m(D

  u, D

  v)

• (A

  ₁

  u,

  v) = ((

  u · Ñ)

  u,

  v)

• (A

  ₂

  u,

  v) = 2n(F

  _r(e(

  u))e(

  u), e(

  v)).

Let the operator A = A₀ + A₁ + A₂. 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₁ 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₁u_n, v) – (A₁u, v) = (((u_n – u) · Ñ)u_n, v) + ((u · Ñ)v, u – u_n)

and then

|(A₁u_n, v) – (A₁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¹(W) is compact, we have u_n u in H¹(W). Therefore

(A₁u_n, v) (A₁u, v), "v Î V.

On the other hand, owing to the lemma 3.3, we prove that A₀ + A₂ is monotone. Indeed, for all u, v Î X we have

(A₂u – A₂v, u – v) = (F_r(e(u))e(u) – F_r(e(v))e(v), e(u) – e(v)) > 0

and

(A₀u – A₀v, u – v) = m(D(u – v), D(u – v)) > 0.

Therefore, owing to lemma 3.1, the operator (A₀ + A₂) + A₁ i.e. A is pseudo-monotone.

Proposition 3.2. The operator A is continuous and coercive.

Proof. It is obvious, since A₀, A₁ and A₂ 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² > C(W)f

_X¢

then problem (8) has exactly one and only one solution.

Proof. We assume that there are two solutions u₁ and u₂ in V of (8) then we have

(Au₁ – Au₂, v) = 0 "v Î V.

Then, for all v Î V, we have

m(Du₁ – Du₂, Dv) + 2n(F_r(e(u₁))e(u₁) – (F_r(e(u₂))e(u₂),

e(v)) + ((u₁ · Ñ)u₁, v) – ((u₂ · Ñ)u₂, v) = 0.

By taking v = u₁ – u₂, we get

Since, owing to lemma 3.3,

2n(F_r(e(u₁))e(u₁) – F_r(e(u₂))e(u₂), e(u₁) – e(u₂)) > 0,

we readily derive that

Now, for some positive constant C₁ depending on W, we have

By using lemma 2.2, we get

On the other hand, since u₁ is a solution of (8) and by taking v = u₁ in (8), we get

Therefore,

We easily derive that u₁ 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₀(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¹(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¹(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²(W) and curl Du Î H¹(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¹(W) and by using the lemma 2.1, we have

(curl w, curl v) = (D²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₀(curl, W) , H(W) and L²(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 Ì ² 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.

#551/02.

• [1] R.A. Adams, Sobolev spaces, Academic Press, New York, (1975).
• [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.
• [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.
• [4] A. Brahmi, Simulation numérique des fluides bipolaires, Mémoire de maîtrise, Université Laval, Québec, (1994).
• [5] P.G. Ciarlet, The Finite Element Method for Elliptic Problems, North-Holland, Amsterdam, 1978.
• [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).
• [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.
• [8] K.K. Golovkin, New Equations Modeling the Motion of a Viscous Fluid, and their Unique Solvability, Proc. Steklov Inst. Math., 102 (1967).
• [9] A.E. Green and R.S. Rivlin, Simple force and stress multipoles, Arch. Rat. Mech. Anal., 17 (1964), 325-353.
• [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).
• [11] O.A. Ladyzhenskaya, Nonstationary Navier-Stokes equations, Amer. Math. Soc. Transl. (2) 25 (1962), 151-160.
• [12] J.L. Lions, Quelques méthodes de résolution des problčmes aux limites non linéaires, Dunod, Paris, (1969).
• [13] J. Necas and M. Silhavy, Multipolar viscous fluids, Quart. Appl. Math., 49 (1991), 247-265.
• [14] Y.U. Ou and S.S. Sritharan, Analysis of regularized Navier-Stokes equations I, Quart. Appl. Math., 49 (1991), 651-685.
• [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²AN 27 (1993), 131-155.
• [16] R. Temam, Navier-Stokes equations, North-Holland, Amsterdam, (1977).
• [17] E. Zeidler, Nonlinear functional analysis and its applications II, Springer-Verlag, (1985).

Publication Dates

• Publication in this collection
  19 July 2004
• Date of issue
  2003