Numerical simulations of non-Newtonian flows between eccentric cylinders by domain decomposition and stream-tube method

J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

Numerical simulations of non-Newtonian flows between eccentric cylinders by domain decomposition and stream-tube method Dana Grecova , Jean-Robert Clermontb,∗ a b

Department of Chemical Engineering, MacGill University, 3610 University Street, Montreal, Que., Canada H3A 2B2 Laboratoire de Rh´eologie, UMR 5520 CNRS, Universit´e Joseph-Fourier, Institut National Polytechnique de Grenoble, Domaine Universitaire, B.P. No 53, 38041 Grenoble Cedex 9, France Received 5 January 2004; received in revised form 27 September 2004; accepted 4 October 2004

Abstract In this paper, we study the influence of rheological properties of inelastic and viscoelastic fluids in flows occurring in the annulus of rotating eccentric cylinders. To compute this flow where vortex regions are encountered, the physical domain is split up into a finite number of sub-domains that are mapped into domains where open and closed streamlines are parallel and straight by means of local mapping functions to be determined numerically. Such simplicity for the mapped streamlines makes it easy to handle time-dependent constitutive models in viscoelasticity. The governing non-linear equations are solved by means of an optimization algorithm. The results point out the influence of non-Newtonian properties of the fluids in annular flows between eccentric cylinders as also in particular cases of journal bearings, notably the role of elasticity. Results with Newtonian inertial flows with large gaps underline the possibility of the approach to investigate flows with multiple vortices. © 2004 Elsevier B.V. All rights reserved. Keywords: Streamlines; Viscoelastic models; Mapping functions; Domain decomposition; Eccentric cylinders; Journal bearings

1. Introduction In this study, flows between eccentrically rotating cylinders are computed by means of the stream-tube method (STM) [1,2], for fluids of different types. Besides the interest of considering such configuration as a benchmark problem [3], such flows are directly connected to journal bearing lubrication. The addition of polymers to Newtonian oils for industrial purposes has resulted in shear thinning and changes in elastic properties of the materials, thus leading to complex behaviour in such eccentric rotating geometries [4–6]. In the literature, numerous numerical studies have been devoted to problems occurring in journal bearings [7–13], with various methods for computing the flow characteristics. In particular, many spectral approaches have been developed in order to simulate flows in eccentric geometries. For example, Davies and Li [14] adopted a pseudo∗

Corresponding author. Tel.: +33 4 76 82 52 92; fax: +33 4 76 82 51 64. E-mail address: [email protected] (J.-R. Clermont).

0377-0257/$ – see front matter © 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.jnnfm.2004.10.004

spectral method with bipolar transformation to simulate nonisothermal flows between eccentrically rotating cylinders, using a White–Metzner viscoelastic differential model for the fluid. Beris et al. [7,8] have studied steady flows between slightly eccentric cylinders, for inelastic and viscoelastic differential fluids, using a Galerkin finite-element method and a spectral finite-element technique. Similarly, the analysis proposed by Chawda et al. [15] involved numerical pseudospectral methods to compute steady-state flow between eccentric rotating cylinders and to explore the stability of threedimensional disturbances, for an upper-convected Maxwell fluid. In order to investigate specific lubrication problems related to the separate and combined effects of shear thinning, temperature thinning and temperature thickening in two-dimensional journal bearings, a recent study involving spectral elements has been proposed by Li et al. [9]. However, as pointed out in our previous papers [16,17], we still are not aware of calculations on memory-integral constitutive equations concerning flows between eccentric cylinders; this is probably because of difficulties arising when the kine-

176

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

matic tensors and stresses on pathlines have to be evaluated. In this paper, the method previously defined in [16,17] is applied to the flow conditions of circulating regions. Our contribution concerns notably the adoption of fluids obeying integral equations. The local transformations to be determined numerically in STM are related to sub-domains involving open or closed streamlines. The spatial variables are defined in simple computational domains. The fluids are assumed to be incompressible and move under steady and isothermal conditions. The respective radii of inner and outer cylinders C0 and C1 are r0 and r1 . Denoting by e, the distance between the axes of the cylinders, we define the eccentricity by ε = e/(r1 − r0 ). The article is organized as follows. In Section 2, the basic equations of the stream-tube method to be applied to flows between rotating cylinders are summarized. The rheological equations for the inelastic and viscoelastic fluids are presented in Section 3. The computational procedure and the numerical results are presented in Section 4. Concluding remarks are given in Section 5.

component does not change its sign, can thus be divided into two curves L1 and L2 , such that L = L1 ∪ L2 , according to the positions of streamline points of zero velocity for the component of V normal to the reference section Sm . The curves L1 and L2 are mapped into rectilinear lines L∗1 and L∗2 . As for an open streamline (Fig. 1a), a mapped closed stream tube of normalized transformed streamlines defines a straight cylinder, the basis of which is included in the reference section ∗ of Ω . When using Carteof the mapped sub-domain Ωm m sian coordinates, local transformations for non-overlapping sub-domains may be defined by with

2. General mappings—domain decomposition and local transformations

• either an open streamline L originating at a point M0 of the reference section or • one of the branches L1 and L2 of a closed streamline originating at points M1 and M2 of the reference section Sm (Fig. 1b).

2.1. Transformations The basic elements on streamline transformations for open and closed streamlines have been detailed in previous papers [2,16,17]. Fig. 1(a and b) presents streamline mappings related to main or vortex regions in a sub-domain Ωm . In both cases, we define plane reference sections that are assumed to be identical in the physical and mapped domains. We now focus on transforming a closed streamline L of the physical sub-domain Ωm (Fig. 1b), by the means of a reference section Sm , where M1 and M2 are the intersecting points of the closed streamline and the reference section. The closed streamline L, which involves two separated regions where the main velocity

M ∗ (X, Y, s) → M(x, y, z)

(1)

and x = αm (X, Y, s);

y = βm (X, Y, s);

z = γm (X, Y, s). (2)

The local variables (X, Y) of Eq. (2) are defined in the plane of ∗ of Ω∗ . In steady flows, we denote by the reference section Sm m χM (χM < +∞) the maximum curvilinear abscissa that may correspond, in the physical domain, to:

∗ as the The variable s is defined in the mapped domain Ωm length of a segment of streamline, using the reference section [17] from the curvilinear abscissa χ given by
t χ= |V (τ)| dτ. (3) t0

Let smin and smax denote the respective minimum and maximum values of s for a given mapped streamline of non-zero velocity, and χmin and χmax are the respective minimum and maximum values of χ on the corresponding streamline in

Fig. 1. Transformation of streamlines: (a) open streamline and (b) closed streamline.

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

the local physical sub-domain. The local transformation Tm : M* (X, Y, s) → M(x, y, z) such that χ − χmin =

(s − smin )(χmax − χmin ) . smax − smin

(4)

This equation implies the following compatibility relation for the local mapping functions α, β and γ
s 2 2 2 χ= [(α s ) + (βs ) + (γs ) ] ds (5) smin

where χ0 denotes the abscissa of the position M0 at the reference section Sm of sub-domain Ωm . The subscript m is omitted in Eq. (5) for reasons of conciseness. It should also be pointed out that the length χM of an open branch of streamline in Ωm that corresponds to the length sM (arbitrarily cho∗ under consideration) is sen in the mapped sub-domain Ωm unknown [17]. Omitting then the subscript m for reasons of brevity, we can express the Jacobian  = |∂(xi )/∂(ξ j )|, assumed to be nonzero, by the following equation  = α X (βY γs − βs γY ) − βX (α Y γs − γY α s )

+ γX (α Y βs − βY α s )

(6)

The derivative operators as also the respective natural and reciprocal basis vectors related to the coordinates (X, Y, s) can be evaluated from Eq. (2) [16,17]. Fig. 2 illustrates local transformations for a bounded physical domain Ω involving non-overlapping sub-domains m=m Ωi such that Ω = ∪m=1 0 Ωm . The open and closed stream tubes of the sub-domains are a priori unknown and correspond to main flow and vortex regions of domain Ω. A reference cross-section Sm is selected in each sub-domain.

177

∗ are assumed to The mapped streamlines in sub-domains Ωm be perpendicular to the reference sections. The analysis enables the computations to be performed in the simple mapped ∗ . The transformed streamlines in these subsub-domains Ωm domains are rectilinear and parallel to the generants of the straight cylinders whose basis is identical to the reference section Sm of the sub-domain Ωm under consideration. For each sub-domain, the kinematics can be expressed in terms of a reference kinematic function φ(X, Y) = w(X, Y, s = s0 ) at the reference section S(X, Y, s0 ). We use the derivative operators ∂/∂x, ∂/∂y and ∂/∂z in terms of ∂/∂X, ∂/∂Y and ∂/∂s [16] to express the components (u, v, w) of the velocity vector V as follows

u=

α s φ(X, Y ) ; 

v=

βs φ(X, Y ) ; 

w=

γs φ(X, Y ) ,  (7)

In Eq. (7),  denotes the Jacobian, assumed to be nonsingular, given by a relation of the type of Eq. (6). The incompressibility condition
·V = 0 is automatically verified from Eq. (7) [2]. Concerning the kinematics, it can be shown [2] that the velocity vector V in a stream tube is given in the natural basis ei (e1 , e2 , e3 ) by the simple relation φ(X, Y ) (8) e3 .  Eq. (8) highlights the simplicity of particle-tracking operations with variables of the mapped sub-domains where the streamlines are rectilinear, for fluids obeying time-dependent constitutive equations. The time evolution of particles along the pathlines may be easily evaluated using the w-velocity component given by Eq. (8) such that
ξ=s
ξ=s   dξ  1 t − t0 = dξ = φ(X, Y ) ξ=0 γs (X,Y,ξ) ξ=0 w(X, Y, ξ) (9)

V =

In Eq. (9), t0 , corresponding to s0 , denotes the reference time associated with the position of the material point at the reference section (X, Y, s0 ) in the sub-domain under consideration. We can assume the following relations at the reference section S for (X, Y, s = 0) [17] α X (X, Y, 0) = 1, (X, Y, 0) = 0, βX γX (X, Y, 0) = 0,

α Y (X, Y, 0) = 0 βY (X, Y, 0) = 1 γY (X, Y, 0) = 0.

(10)

For a sub-domain Ωm , the primary unknowns to be considered are the three mapping functions αm , βm , γ m related to ∗ and Ω and the pressure p , given a constitutive domains Ωm m m equation for the moving fluid and the sub-domains Ωm . Under isothermal conditions, we write the governing equations involving [17]: Fig. 2. Decomposition of domain Ω into sub-domains Ωm by local mappings.

• the momentum conservation law
·σ˜ = F, where σ˜ denotes the stress tensor expressed in terms of the extra-stress tensor T˜ by the relation σ˜ = −p·I˜ + T˜ and F are the external

178

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

Fig. 3. Local transformation function for the flow in the eccentric annulus geometry. (a) Sub-domain Ωm of Ω with a reference section Sm and (b) corresponding ∗ for the local transformation. mapped domain Ωm

forces. The three dynamic equations can thus be formally written as Emj (αm , βm , γ m , pm ) = 0 (j = 1, 2, 3), for m = 1, 2, . . ., M; • the compatibility Eq. (5). These equations are written with variables (X, Y, s) of the computational sub-domains. 2.2. Equations for flows between eccentric cylinders The two-dimensional flow domain Ω corresponding to the annulus of two cylinders is shown in Fig. 3(a). The trans∗ of sub-domain Ω , involving the refformed domain Ωm m ∗ erence section Sm and the parallel mapped streamlines, is shown in Fig. 3(b). The physical domain, which can be split up into two or more sub-domains, is referred to polar coordinates (x1 = r, x2 = θ). (x1 , x2 )m denote local polar coordinates in sub-domains Ωm of Ω, defined with respect to the refer∗ , identical ence flow section Sm . Using a reference section Sm ∗ of in shape to Sm , the corresponding mapped domain Ωm sub-domain Ωm is defined according to Fig. 3(b). 2.2.1. Case of two sub-domains We now consider the flow situation where the inner cylinder (of radius r0 ) rotates with an angular velocity ω and the outer cylinder (of radius r1 ) is at rest. According to results obtained for such a geometry with a Newtonian fluid and ignoring inertial effects [2,16], we choose to adopt two elementary sub-domains Ω1 and Ω2 of the half-domain Ω limited

by the azimuth angles θ = 0 and π, for reasons of symmetry. For m = 1 and 2, let (x1 , x2 )m be local polar coordinates in sub-domains Ωm of Ω, defined with respect to the reference ∗ , identical in flow section Sm . Using a reference section Sm ∗ shape to Sm , we denote by Ωm the transformed domain of sub-domain Ωm (m = 1 and 2) where the mapped streamlines are parallel to a given direction (Fig. 3b). Sub-domain Ω1 , of reference section at θ = 0, involves only open streamlines. Closed and open streamlines are to be taken into account in sub-domain Ω2 , of reference section at θ = π. The angle θ 1 , a priori unknown, corresponds to the limiting azimuthal section between the sub-domains Ω1 and Ω2 (Fig. 4). Sub-domain ∗ is referred to a coordinate system ξ j expressed by Ωm ξ 1 = s;

ξ 2 = R.

(11)

Two reference kinematic functions φ1 and φ2 are considered at sections θ = 0 and π, respectively. The variables are defined according to Eq. (6) and the variable R corresponds to the radius r at the reference sections of the sub-domains. The ∗ →Ω , two local transformations (for m = 1 and 2): Tm : Ωm m * where M (s, R) → M(r, θ) involve the unknown local functions α(s, R) and β(s, R) (or λa (s, R) and λb (s, R)) such that r = α(s, R) = λa (s, R) · R;

θ = β(s, R) =

λb (s, R) · s R (12)

The coordinate R corresponds to the radius r at the refer∗ (s , R) of sub-domain Ω∗ . The coordinate s ence section Sm 0 m is defined in relation to the curvilinear abscissa χ (Eq. (3)). It

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

179

The relevant unknowns of the problem are expected to be determined iteratively. The function ψ* (R), a priori unknown, must be updated at each step of the numerical procedure. We also put the following conditions at the reference section α R (s = 0, R) = 1;

βs (s = 0, R) = 0.

(16)

From Eq. (14), the velocity components become α s (s, R) · H ∗ (R) α(s, R) ·  β (s, R) vθ = − s · H ∗ (R).  vr = −

(17)

In the flow domain, the coordinate θ on a streamline of nonzero velocity can be expressed by the following equation:


θ = θ0 + H ∗ (R)

Fig. 4. (a) Physical domain Ω partitioned into sub-domains Ω1 and Ω2 . (b) Mapped sub-domains Ω1∗ and Ω2∗ with boundary conditions.

should be pointed out that the streamline lengths in mapped sub-domains Ω1∗ and Ω2∗ have been renormalized for computational purposes. − α · β is The Jacobian  = det |∂xi /∂ξ j | = α s · βR s R non-zero. Thus, derivative operators relating the coordinates (r, θ) to the new system of local coordinates (s, R) are expressed by ∂ 1 = ∂r [βR · ∂/∂s − βs · ∂/∂R] ∂ 1 = ∂θ [αR · ∂/∂s − α s · ∂/∂R]

(13)

vr = −

(14)

In Eq. (14) and the following relationships, we omit the subscript m indicating the local variables, for the purpose of simplicity. ψ* (R) denotes the transformed stream function of ∗ ψ defined from a reference kinematic function at section Sm of coordinates (s0 , R) such that dψ∗ (R) = vθ (s = 0, R) = H ∗ (R). dR

βs (ξ, R) dξ.  ξ=0

(18)

Concerning the particle tracking for a material point moving on its streamline in a sub-domain Ωm , the movement versus time can be expressed readily [16,17]. Ignoring inertia and body forces, the dynamic equations along the two axes (r, θ) of the physical domain Ω, with the local variables (s, R) of the sub-domain Ωm under consideration, are expressed as   ∂p ∂(αT rr ) ∂(αT rr ) ∂p α βs · · − βR · + βs · − βR ∂R ∂s ∂R ∂s − βR ·

∂(αT rr ) ∂(αT rθ ) ∂(T rθ ) − α s · + α R · =0 ∂s ∂R ∂s

∂p ∂p ∂(αT rθ ) ∂(αT rθ ) − α s · + βs · − βR · ∂R ∂s ∂R ∂s rr rr ∂(T ) ∂(T ) + α s · − α R · =0 ∂R ∂s

(19)

α s ·

Using variables (s, R), the velocity component equations written in terms of a stream function ψ(r, θ) as vr = (1/r)∂ψ(r, θ)/∂θ; vθ = −∂ψ(r, θ)/∂r in polar coordinates are given, in a ∗ , by the following relations sub-domain Ωm α s (s, R) dΨ ∗ (R) · α(s, R) ·  dR ∗ β (s, R) dΨ (R) vθ = − s · .  dR

ξ=s

(15)

(20)

2.2.2. Case of multiple sub-domains With inertial flows between eccentric cylinders, more than two sub-domains of the physical domain Ω are required for the flow analysis. In such cases, the sub-domains are considered in the total flow domain and, as soon as reference sections are defined, a similar approach to that proposed in the previous sub-sections for two sub-domains can be used. Inertial terms are thus written in the dynamic equations.

3. Constitutive equations Inelastic and viscoelastic models have been adopted in the flow computations: (1) Inelastic models: ˜ • a Newtonian equation (T˜ = 2 µ D);

180

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

• two generalized Newtonian fluids of different viscosity characteristics, where the stress tensor is given by ˜ D ˜ T˜ = 2µ(D)

Table 1 Parameters used in the K-BKZ memory-integral equation (after Chai and Yeow [19])

(21)

p

λp (s)

ap

˜ denotes the rate-of-deformation tensor and µ where D denotes the fluid viscosity which obeys Carreau Law [18] in Eq. (21). In the following, these fluids are denoted by CAR1 and CAR2. The viscosity curves for the two fluids, shown in Fig. 5, versus the velocity gradient are defined such that the CAR2 fluid has more variations in viscosity than the CAR1 fluid. (2) Viscoelastic models: • an upper-convected Maxwell model (UCM) (e.g. Bird et al. [18]) of equation:

1 2 1

1.04 8.9 × 10−2 5.2 × 10−4

0.55 5.46 3.06 × 10−3

δT˜ ˜ + T˜ = 2η0 D (22) δt where λ stands for the relaxation time and η0 is the constant viscosity of the material. In this equation, ␦/␦t denotes the Gordon/Schowalter convected derivative operator [19] such that

λ

δT˜ ∂T˜ ˜ − ξ D) ˜ T˜ − T˜ (L ˜ − ξ D) ˜ T = + V · ∇ T˜ − (L δt ∂t (23) • a multi-mode memory-integral model of the K-BKZ type that proved to fit the experimental data for a polymer solution (Chai and Yeow [20]), expressed by the following equation  
t  3 ap t − t T˜ (t) = exp − λp λp −∞

In Eq. (22), ap and λp are the moduli and corresponding relax˜ t the Finger tensor and B ˜ −1 ation times, respectively, B t is the Cauchy tensor related to the configuration at times t . Material constants of the fluid, α and β, determined experimentally, are found to be α = 500 and β = 0.001. The coefficients ap and λp of the model are given in Table 1 and the corresponding material functions are shown in Fig. 6, according to data presented in [20]. The Finger and Cauchy tensors to be considered in the constitutive Eq. (22) are related to the deformation gradient ∂xm tensor F˜ t (t) = tγ by the following relations [18]: ∂xt

˜ t (t ) = [T F˜ t (t )F˜ t (t )]−1 B ˜ t−1 (t ) = T F˜ t (t )F˜ t (t ) B

(25)

To evaluate the deformation gradient tensor F˜ t (t), we applied the analysis already defined in [16,17], starting from Adachi’s work [21,22]. The approach allows simple relations

p=1

×

˜ t (t) αB dt ˜ ˜ −1 α − 3 + β trBt + (1 − β) trB

(24)

t

Fig. 5. Viscosity curves for the two non-Newtonian, purely viscous fluids, versus the velocity gradient. (—) CAR1 fluid and (– · –) CAR2 fluid.

Fig. 6. Material functions for the K-BKZ viscoelastic M1 fluid, after Chai and Yeow [19].

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

181

to be obtained for the components of F˜ t (t). Thus, the Cauchy and Finger tensor components can be determined by using Eq. (25).

4. Numerical procedure and results In the present approach, where the stream-tube method is associated with domain decomposition techniques, the governing equations involve the dynamic equations, together with boundary conditions and compatibility equations at the interfaces of the sub-domains [16,17]. The primary unknowns are the pressure and the local mapping functions. Stress components Tij , related to the constitutive equations, are expressed in terms of the transformation functions and their derivatives. Owing to the rectangular mesh defined in the computational sub-domains, finite-difference formulae are adopted to approximate the derivatives involved in the equations and regular meshes were used to discretize the equations and unknowns. The functions Φm at the reference sections are unknown and evaluated iteratively by the Levenberg–Marquardt optimization algorithm [23] in mapped sub-domains Ωi∗ , given an initial estimate from the lubrication theory [16]. In the present application, the reference sections remain unchanged during the iterative process. New values of the separating angles between the domains are then obtained. Consequently, the sizes of the physical and mapped sub-domains Ωi are modified and new meshes are built in the computational sub-domains. The process is repeated up to verification of the convergence criteria. 4.1. Results without inertia When inertia is ignored, the calculations involved two subdomains with the half-domain of the total flow domain, for reasons of symmetry, as pointed out previously. The computations were run on a Pentium 2.4 GHz processor workstation. The convergence behaviour was tested by considering, with the K-BKZ fluid, three regular meshes M1, M2 and M3 in the computational sub-domains (see Fig. 4) with 3200, 2450 and 1800 unknowns, respectively, with double-precision variables, for a journal bearing configuration (r0 = 5.084 cm and r1 = 6 cm, ε = 0.5, ω = 150 rad/s) The evolution in norm of the objective function to 10−6 with the three meshes was found to require a number of iterations of 25, 35 and 42, respectively. For all the grids, the maximum CPU time for K-BKZ was found to be 180 s by iteration. Concerning the difference between the computed solutions with the three meshes with K-BKZ, the difference was not perceptible for the streamline plots. The flow calculations reported in this study have concerned the use of mesh M2. The dimensionless Weissenberg number, We, is given by We = λω,

(26)

Fig. 7. Computed streamlines in the developed flow domain for a Newtonian fluid at different eccentricities: (a) ε = 0.6, (b) ε = 0.7 and (c) ε = 0.8 (r0 = 5.084 cm and r1 = 6 cm).

where λ denotes the characteristic time for UCM or the mean characteristic time of the multi-mode viscoelastic model (Eq. (24)), which was found to be 0.6 s. A satisfactory convergence was obtained for the fluids considered. Concerning viscoelastic computations, it should be pointed out that our previous results [17] have been found to be in good agreement with available literature data [10] for a UCM differential equation and practically identical for the one-mode differential UCM constitutive equation and its integral equivalent. For the present calculations in a journal bearing geometry (r0 = 5.084 cm and r1 = 6 cm, ε = 0.5), the limit of convergence with the K-BKZ fluid was found to correspond to an angular speed of ω = 150 rad/s (We = 90). 4.1.1. Kinematic results Numerical streamline results are provided for four fluids investigated, namely Newtonian, generalized Newtonian (CAR1 and CAR2) and K-BKZ fluids, for the same values of eccentricity and angular speeds of the inner cylinder, when neglecting inertia. The geometrical characteristics, r0 = 5.084 cm and r1 = 6 cm, correspond to journal bearing specifications. Owing to the small gap, the plotted streamlines are presented on a developed domain with a dimensionless radius r* , for an angular velocity of ω = 5.24 rad/s. Though these results are presented for a small velocity speed, the streamlines of Figs. 7–10(a–c) clearly illustrate significant kinematic differences in behaviour for the four fluids,

182

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

Fig. 8. Computed streamlines in the developed flow domain for the CAR1 fluid at different eccentricities: (a) ε = 0.6, (b) ε = 0.7 and (c) ε = 0.8 (r0 = 5.084 cm and r1 = 6 cm).

Fig. 10. Computed streamlines in the developed flow domain for the KBKZ fluid at different eccentricities: (a) ε = 0.6, (b) ε = 0.7 and (c) ε = 0.8 (r0 = 5.084 cm and r1 = 6 cm).

Fig. 9. Computed streamlines in the developed flow domain for the CAR2 fluid at different eccentricities: (a) ε = 0.6, (b) ε = 0.7 and (c) ε = 0.8 (r0 = 5.084 cm and r1 = 6 cm).

with larger vortex zones for the viscoelastic fluid, for the eccentricities ε = 0.6, 0.7 and 0.8. With all the fluids tested, the recirculating regions increase in size, as expected, with the eccentricity. Concerning the generalized Newtonian CAR1 and CAR2 fluids, the computed streamlines indicate more important vortex zones than the Newtonian fluid. The more shear-thinning fluid (CAR2) presents larger vortex zones than the CAR1 fluid. For all the eccentricities, the recirculating regions are found to be greater in size with the viscoelastic K-BKZ fluid, which shows that the elastic effects are preponderant with respect to viscosity effects owing to the shear-thinning properties of the fluids considered. The differences in velocity, with four fluids investigated, are shown in Fig. 11, where the dimensionless velocities, v∗θ , for θ = 0, ε = 0.7 and ω = 20 rad/s (r0 = 5.084 cm and r1 = 6 cm) are plotted. In this case, r0 ω is used as a dimensionless factor. These velocity results show an increase in nonlinearity for the velocity curves, from the Newtonian case to the viscoelastic situation and also confirm the preponderant influence of the elastic properties, with respect to viscosity effects related to CAR1 and CAR2 fluids.

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

Fig. 11. Dimensionless velocity component v∗ θ for θ = 0, ε = 0.7, ω = 20 rad/s (r0 = 5.084 cm and r1 = 6 cm) for the four fluids: (· · ·) Newtonian fluid; (– – –) CAR1 fluid; (– · –) CAR2 fluid; (—) K-BKZ fluid.

183

Fig. 13. Comparison of our results (symbol ×) with those of Beris et al. [8]; () Huang et al. [10]; () Fan et al. [24]; (—) for the load amplitude, with a UCM fluid.

4.1.2. Stress results Fig. 12 presents an example of dimensionless stress results Tθθ* , with the same geometrical characteristics (r0 = 5.084 cm and r1 = 6 cm), with ε = 0.7 and θ = 0, for the four fluids investigated. The curves indicate a maximum, the intensity of which is greater for the viscoelastic material. These results confirm the previous comments on elastic effects versus the influence of viscosity. With regard to lubrication problems in journal bearings, we evaluate the torque C and the load F per unit length on the inner cylinder of radius r0 . The torque can be expressed by means of the stress component Trθ , according to the relation:
C = −r02

2π 0

T rθ |r=r0 dθ

(27)

Fig. 12. Dimensionless stress Tθθ* in the eccentric geometry (r0 = 5.084 cm and r1 = 6 cm), with ε = 0.7 and θ = 0 for the four fluids: (· · ·) Newtonian fluid; (– – –) CAR1 fluid; (– · –) CAR2 fluid; (—) K-BKZ fluid.

Fig. 14. Dimensionless load F* on the inner cylinder versus the angular speed for the Newtonian, Carreau1 and K-BKZ fluids (r0 = 5.96 cm, r1 = 6 cm and ε = 0.5). (– · –) Newtonian fluid; (– – –) CAR1 fluid; (—) K-BKZ fluid.

Fig. 15. Dimensionless torque C* on the inner cylinder versus the angular speed for the Newtonian, Carreau1 and K-BKZ fluids (r0 = 5.96 cm, r1 = 6 cm and ε = 0.5). (– · –) Newtonian fluid; (– – –) CAR1 fluid; (—) K-BKZ fluid.

184

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185

Fig. 16. Computed streamlines for inertial Newtonian flows: (a) CASE1, (b) CASE2 and (c) CASE3. Table 2 Geometrical and kinematic conditions for flow computations with inertia for r0 = 1.5 cm, r1 = 3.0 cm and eccentricity ε = 0.42 Configuration name

ω0 (rad/s)

ω1 (rad/s)

CASE1 CASE2 CASE3

5.2 −4.5 −5.2

5.2 −7.2 0.0

In order to provide comparisons with literature data, Fig. 13 presents our results on the computed load amplitude exerted by the fluid on the inner cylinder and numerical data from Beris et al. [8], Huang et al. [10] and Fan et al. [24] with a UCM fluid. It can be observed that our results are in agreement with those of these authors. The dimensionless curves of load F* and torque C* versus the angular speed, related to the inelastic Newtonian, Carreau1 and the K-BKZ viscoelastic fluids, in journal bearing conditions (r0 = 5.96 cm, r1 = 6 cm, eccentricity ε = 0.5) are presented in Figs. 14 and 15, respectively, for a large range of angular speeds. The angular velocity of 150 rad/s for the K-BKZ fluid corresponds to the limit of convergence with this fluid. These figures also point to the dominant influence of elastic properties on the curves, with respect to viscosity effects.

cylinders and investigate the influence of rheological properties on the flow characteristics, using domain decomposition and local transformation functions. The analysis enabled us to compute the main flow characteristics in journal bearing geometries at significant angular speeds close to journal lubrication situations as well as inertial flows with a large gap between the cylinders. It should also be underlined that elastic effects have been considered for a viscoelastic fluid obeying a K-BKZ memory-integral constitutive equation without particle–tracking problems. For a UCM fluid, good agreement was found with literature results. The study of inelastic and viscoelastic fluid flows in journal bearing conditions has led to significant differences in predictions for the kinematics and stresses, even at moderate angular speeds. In particular, the results highlight the preponderant influence of elastic properties on the flow characteristics, as compared to shearthinning effects. We also wish to point out the possibility of using the method to compute complex flow situations as those considered for inertial flows with a large gap. The numerical applications also help to highlight the efficiency and robustness of the Levenberg–Marquardt optimization algorithm used for solving the governing equations associated with the stream-tube method.

4.2. Results with inertia References When taking into account inertial effects, the computations were carried out with more than two sub-domains, as stated previously. The analysis was applied to a Newtonian fluid, under the conditions given in Table 2, allowing rotation of both cylinders. These flow configurations required four sub-domains for the simulations. The computed streamlines shown in Fig. 16(a–c) for CASE1, CASE2 and CASE3 illustrate the complexity of the flows from a kinematic point of view and also underline the possibility of using the approach to solve complex flow situations with several vortex zones.

5. Conclusions In this paper, we have applied an original approach to simulate two-dimensional complex flows between eccentric

[1] J.-R. Clermont, Analysis of incompressible three-dimensional flows using the concept of stream tubes in relation with a transformation of the physical domain, Rheol. Acta 27 (1988) 357. [2] J.-R. Clermont, M. Normandin, D. Radu, Some remarks on the concept of stream tubes for numerical simulations of complex fluid flows, Appl. Control Cybern. 29 (2) (2000) 535–554. [3] R.G. Owens, T.N. Phillips, Computational Rheology, Imperial College Press, London, 2002. [4] J. Frˆene, D. Nicolas, B. Degueurce, D. Berthe, M. Godet, Lubrification hydrodynamique—Paliers et But´ees, Editions Eyrolles, Paris, 1990. [5] B.P. Williamson, K. Walters, T.W. Bates, R.C. Coy, A.L. Milton, The viscoelastic properties of multigrade oils and their effect on journalbearing characteristics, J. Non-Newtonian Fluid Mech. 73 (1997) 115–126. [6] A. Berker, M.G. Bouldin, S.J. Kleis, W.E. Van Arsdale, Effects of polymer on flow in journal bearings, J. Non-Newtonian Fluid Mech. 56 (1995) 333–345.

D. Grecov, J.-R. Clermont / J. Non-Newtonian Fluid Mech. 126 (2005) 175–185 [7] A.N. Beris, R.C. Armstrong, R.A. Brown, Finite element calculation of viscoelastic flow in a journal bearing: I. Small eccentricites, J. Non-Newtonian Fluid Mech. 16 (1984) 141–172. [8] A.N. Beris, R.C. Armstrong, R.A. Brown, Spectral/finite element calculations of the flow of a Maxwell fluid between eccentric rotating cylinders, J. Non-Newtonian Fluid Mech. 22 (1987) 129–167. [9] X.K. Li, D. Rh. Gwynllyw, A.R. Davies, T.N. Phillips, On the influence of lubricant properties on the dynamics of two-dimensional journal bearings, J. Non-Newtonian Fluid Mech. 93 (2000) 29–59. [10] X. Huang, N. Phan-Thien, R.I. Tanner, Viscoelastic flow between eccentric rotating cylinders: unstructured control volume method, J. Non-Newtonian Fluid Mech. 64 (1996) 71–92. [11] M.A. Kelmanson, A boundary integral equation method for the study of slow flow in bearings with arbitrary geometries, Trans. ASME 108 (1984) 260–264. [12] P.S. Ramesh, M.H. Lean, A boundary-integral equation for Navier–Sokes equations—application to flow in annulus of eccentric cylinders, Int. J. Numer. Methods Fluids 13 (3) (1991) 355– 369. [13] R.-X. Dai, Q. Dong, A.Z. Szeri, Flow between eccentric rotating cylinders: bifurcation and stability, Int. J. Eng. Sci. 130 (1992) 1323–1340. [14] A.R. Davies, X.K. Li, Numerical modeling of pressure and temperature effects in viscoelastic flow between eccentrically rotating cylinders, J. Non-Newtonian Fluid Mech. 54 (1994) 331– 350.

185

[15] A. Chawda, M. Avgousti, Stability of viscoelastic flow between eccentric rotating cylinders, J. Non-Newtonian Fluid Mech. 63 (1996) 97–112. [16] D. Grecov Radu, J.-R. Clermont, Numerical study of flows of complex fluids between eccentric cylinders using transformation functions, Int. J. Numer. Methods Fluids 40 (5) (2002) 669–696. [17] D. Grecov Radu, M. Normandin, J.-R. Clermont, A numerical approach for computing flows by local transformations and domain decomposition using an optimization algorithm, Comput. Methods Appl. Mech. Eng. 191 (39–40) (2002) 4401–4419. [18] R.B. Bird, R.C. Armstrong, O. Hassager, Dynamics of Polymer Liquids, Wiley, 1976. [19] R.R. Huilgol, N. Phan-Thien, Fluid Mechanics of Viscoelasticity, Elsevier, Amsterdam, 1997. [20] M.S. Chai, Y.L. Yeow, Modelling of fluid M1 using multiplerelaxation-time constitutive equations, J. Non-Newtonian Fluid Mech. 35 (1990) 459–470. [21] K. Adachi, Calculation of strain histories with Protean coordinate systems, Rheol. Acta 22 (1983) 326–335. [22] K. Adachi, A note on the calculation of strain histories in orthogonal streamline coordinate systems, Rheol. Acta 25 (1986) 555–563. [23] A. Gourdin, M. Boumahrat, M´ethodes num´eriques appliqu´es, in: Technique et Documentation, Editions Lavoisier, Paris, 1989. [24] Y. Fan, R.I. Tanner, N. Phan-Thien, Galerkin/least-square finite element methods for steady viscoelastic flows, J. Non-Newtonian Fluid Mech. 84 (1999) 233–256.