Influence of Shear Thinning on Viscoelastic Fluid–Structure Interaction

Apr 13, 2011 - We compute numerically the steady flow of a viscoelastic FENE-P fluid in a two-dimensional collapsible channel, in which a zero-thickne...
0 downloads 0 Views 2MB Size
ARTICLE pubs.acs.org/IECR

Influence of Shear Thinning on Viscoelastic FluidStructure Interaction in a Two-Dimensional Collapsible Channel Debadi Chakraborty and J. Ravi Prakash* Department of Chemical Engineering, Monash University, Melbourne, VIC 3800, Australia ABSTRACT: We compute numerically the steady flow of a viscoelastic FENE-P fluid in a two-dimensional collapsible channel, in which a zero-thickness membrane under constant tension has replaced part of the top wall. This geometry has been studied extensively as a benchmark to examine fluid-structure interactions between a flowing Newtonian fluid and a deformable wall. The rheological behavior of the viscoelastic FENE-P fluid is described in terms of a conformation tensor model. The mesh equation and transport equations are discretized using the discrete elastic viscous stress split-traceless gradient (DEVSS-TG/SUPG) mixed finiteelement method. The influence of viscoelasticity and shear thinning on flow patterns and stress profiles is examined, and comparisons with Newtonian predictions are reported. The existence of a limiting Weissenberg number beyond which computations fail is demonstrated. The extent of shear thinning is shown to be a key factor in determining the nature of the fluidstructure interaction.

’ INTRODUCTION Numerical simulation of blood transportation through the human cardiovascular system is an intense area of research. Blood shows complex rheological properties because of its various suspended components. The flow in small blood vessels, such as arterioles, venules, and capillaries is identified with the microcirculation, while flow in larger arteries is identified with the macrocirculation. It has been well-known that in mediumto-large arteries, the NavierStokes equations for an incompressible viscous fluid are considered to be a good model for blood flow. Most studies in literature have so far been focused on the macrocirculation, where the Reynolds numbers are high, tube diameters are large, and blood can be considered a Newtonian fluid. Because of the particulate nature of blood, it behaves as a non-Newtonian fluid at low shear rates, and in small blood vessels.13 It is clear that, to model blood flow in the microcirculation, where the flow rates are low, blood must be treated as a shear-thinning, viscoelastic, and thixotropic fluid. The objective of this work is to examine the flow of a FENE-P viscoelastic fluid in a simple two-dimensional (2D) collapsible channel, as a preliminary study of the types of behavior that might arise when more-realistic models of blood and blood vessels are simulated under conditions in which the viscoelastic character of blood and the elastic nature of blood vessel walls become important. Laboratory experiments performed on flow through collapsible tubes have demonstrated the complex and nonlinear nature of the dynamics of such a system, with a multiplicity of selfexcited oscillations.4 As discussed in the literature, the simplest numerical model that captures some of the rich behavior observed in fluidstructure interaction problems is a 2D model, with fluid flowing in a rigid parallel sided channel, where part of one wall is replaced by a membrane under tension, as shown in Figure 1. This geometry has been studied extensively by Pedley and co-workers57 in the case of Newtonian fluids. The flexible wall has been treated earlier as an elastic membrane of zero thickness, with stretching along the flow direction and bending r 2011 American Chemical Society

stiffness of the membrane neglected.68 While initial studies seemed to suggest that, within the range of Reynolds numbers and transmural pressures considered, a converged numerical solution could only be achieved for relatively large values of membrane tension.6,8 Luo and Pedley7 subsequently found that, with the help of an algorithm capable of time-dependent simulations, converged steady-state solutions could be obtained at arbitrary values of membrane tension, for sufficiently small values of the Reynolds number. Interestingly, Luo and Pedley7 showed that self-excited oscillations developed when the membrane tension was reduced below a critical value for sufficiently large Reynolds numbers. More recently, several improvements to the basic model of Pedley and co-workers have been made, with more realistic treatments of the elastic wall9,10 and an extension of the geometry to enable description of flows in three-dimensional (3D) compliant tubes.11 Many of the complex dynamical features observed experimentally have been reproduced in these simulations. Attempts to understand the origin of self-excited oscillations have also been made.9 However, as mentioned earlier, in all these studies, the fluid has always been treated as Newtonian. Recently, Chakraborty et al.12 studied the effect of viscoelasticity on flow in a 2D collapsible channel, considering three different viscoelastic fluid models—the Oldroyd-B, the FENE-P, and Owens’ model for blood13—along with a zero-thickness membrane model with constant tension for the collapsible wall.6 It was shown that the predictions for the different viscoelastic fluids differ significantly from each other, with the key factor being the extent of shear thinning predicted by the individual models. In particular, it was shown that viscoelastic fluids behave identically to Newtonian fluids, provided that the viscosity of the Special Issue: Ananth Issue Received: January 25, 2011 Accepted: April 13, 2011 Revised: April 7, 2011 Published: April 13, 2011 13161

dx.doi.org/10.1021/ie200173b | Ind. Eng. Chem. Res. 2011, 50, 13161–13168

Industrial & Engineering Chemistry Research

ARTICLE

Figure 1. Geometry of the two-dimensional (2D) collapsible channel; the segment BC is an elastic membrane.

Figure 2. Contribution of the microstructure to the total viscosity (ηp) for the FENE-P fluid in steady shear flow as a function of the local fi, at a constant value of the relaxation time λ0. Weissenberg number W

two fluids at the location of the maximum shear rate in the channel is the same. However, the influence of the degree of shear thinning of the viscoelastic fluid was not examined systematically. The finite extensibility parameter bM, in the FENE-P model controls the extent of shear thinning experienced by the fluid, and is consequently a convenient parameter for examining the influence of shear thinning. In this work, we compute the flow of a FENE-P model fluid for various values of bM and delineate the role of shear thinning. A discrete elastic viscous stress splittraceless gradient (DEVSS-TG/SUPG) finite-element method is used to solve for the fluid velocity and stress field, and the shape of the collapsible membrane, which is an unknown free boundary.14 The conformation tensor is used to represent the microstructural state of the FENE-P fluid. The plan of this paper is as follows. We first discuss the formulation of the problem and the FENE-P fluid model. The various predictions of the FENE-P fluid model, along with a comparison of predictions for a Newtonian fluid, then are presented. Finally, concluding remarks are drawn in the last section.

’ PROBLEM FORMULATION To explore the nature of fluidstructure interaction between a steadily flowing FENE-P fluid and a zero-thickness membrane, we have adopted the geometry studied earlier by Luo and Pedley.6 As illustrated in Figure 1, in units of W (the width of the channel), the dimensions of the channel are Lu = 7W, L = 5W, and Ld = 7W. Here, Q is the flow rate, pe the external pressure on the membrane, pd the pressure on the wall at the downstream boundary (denoted by location E), L the length of the deformable membrane, and h the minimum height of the gap between the bottom wall of the channel and the deformable membrane. We have also assumed, following Luo and Pedley,6 that the tension in the flexible wall is constant and

Figure 3. Minimum value of the smallest eigenvalue (m1) and maximum value of the largest eigenvalue (m3) in the entire flow domain, for a FENE-P model, as a function of Wi at R = 45 and bM = 100.

the shape of the flexible part is governed by the normal force acting upon it. Governing Equations. The equations of motion for steady, incompressible flow in the absence of body forces are =3v ¼ 0

ð1Þ

Fv 3 =v ¼ = 3 T

ð2Þ

where F is the density of the liquid, v is the velocity, and = denotes the gradient. The Cauchy stress tensor is given as T = pI þ τ þ σ (where p is the pressure, I the identity tensor, τ the viscous stress tensor, and σ the elastic stress tensor). The viscous stress tensor is τ = 2ηsD, where ηs is the solvent viscosity and D = 1/2(=v þ =vT) is the rate-of-strain tensor. The elastic stress depends on the microstructural state of the fluid, represented here by the dimensionless conformation tensor M.14 For the FENE-P model, this relation takes the form    ηp, 0 bM  1 σ ¼ MI ð3Þ bM  ððtr MÞ=3Þ λ0 where ηp,0 is the contribution of the microstructural elements to the zero-shear-rate viscosity, λ0 is the constant characteristic relaxation time of the microstructure, and bM is the finite extensibility parameter, defined as the ratio of maximum length squared of the microstructural element to its average length squared at equilibrium. The conformation tensor obeys the following evolution equation: DM þ v 3 =M  =vT 3 M  M 3 =v Dt    1 bM  1 ¼  MI λ0 bM  ððtr MÞ=3Þ 13162

ð4Þ

dx.doi.org/10.1021/ie200173b |Ind. Eng. Chem. Res. 2011, 50, 13161–13168

Industrial & Engineering Chemistry Research

ARTICLE

Figure 4. Minimum value of the smallest eigenvalue (m1) and maximum value of the largest eigenvalue (m3) in the entire flow domain, for the FENE-P model at different finite extensibility bM.

(3) On the rigid walls, we apply no-slip boundary conditions (v = 0). (4) The fully developed flow boundary condition is imposed, n 3 =v = 0, at the downstream boundary. (5) On the bottom wall of the downstream boundary (at location E in Figure 1), the pressure of the fluid (pd) is set equal to zero. (6) The conformation tensor does not change along the streamlines, because the flow is fully developed at the upstream inflow.14,18 Thus, v 3 =M ¼ 0 Figure 5. Comparison of contours of axial velocity (vx) in the domain at tension ratio of 45 for Newtonian (black line) and different values of bM for the FENE-P model, bM = 100 (red dashed line), bM = 50 (blue dots), and bM = 2 (thin green line) at Re = 1.0 and β = 0.0071.

To handle free boundaries, we use a boundary-fitted finite-elementbased elliptic mesh generation method,1417 which involves solving the following elliptic differential equation for the mapping x = x(ξ), ~ 3 =ξ ¼ 0 =3D

ð5Þ

where ξ is a vector of positions in the computational domain and the ~ , is a function of ξ, analogous to a diffusion coefficient, dyadic, D which controls the spacing of the coordinate lines.17 Boundary Conditions and Discretization. The following boundary conditions are prescribed: (1) A fully developed velocity profile is specified at the upstream boundary in the form vx = f(y), and vy = 0. (2) At the flexible wall, (a) We impose on the momentum equation, (i) t 3 v = 0, where t is the unit tangent to the flexible wall, and (ii) a force balance in the normal direction through the traction boundary condition: nn : T ¼  pe þ χ=II 3 n

ð6Þ

where n is the unit normal to the flexible wall, =II denotes the surface gradient operator, and χ is the fixed tension in the flexible wall. (b) We impose on the mapping equation, (i) n 3 v = 0 in the normal direction, and (ii) a uniform node distribution in the tangential direction.

ð7Þ

Equations 15 are converted to a set of algebraic equations by the DEVSS-TG finite-element method,14,19 which introduces the traceless interpolated velocity gradient L:14 L  =v þ

1 ð= vÞI ¼ 0 tr I 3

ð8Þ

In the transport equations, the rate-of-strain tensor D is calculated from the interpolated velocity gradient L. The weighted residual form of eqs 1, 2, 4, 5, and 8 yields a large set of coupled nonlinear algebraic equations, solved by Newton’s method with analytical Jacobian, frontal solver, and first-order arc length continuation in parameters.14,20 Dimensionless Numbers and Choice of Parameter Values. We have nondimensionalized the various physical quantities, by scaling lengths with W, velocities with U0, and pressure and stresses with η0U0/W. Here, U0 is the average inlet velocity, η0 the zero-shear-rate solution viscosity (η0 = ηs þ ηp,0), ηs is the solution viscosity, and ηp,0 the contribution of the microstructural elements to the zero-shear-rate viscosity. (For a Newtonian fluid, η0 is just the constant Newtonian viscosity.) Nondimensionalization of the governing equations and boundary conditions yields the following dimensionless numbers: Re ¼

FWU0 ; η0

β¼

ηs ; η0

Wi ¼

λ0 U 0 ; W

Ca ¼

η0 U0 ð9Þ χ

where Re is the Reynolds number, Ca is analogous to a capillary number, β is the viscosity ratio, and Wi is the inlet Weissenberg number. It is convenient to define a local Weissenberg number fi, which measures the nondimensional shear rate anywhere in W fi = λ0γ). _ the flow (W 13163

dx.doi.org/10.1021/ie200173b |Ind. Eng. Chem. Res. 2011, 50, 13161–13168

Industrial & Engineering Chemistry Research

ARTICLE

Figure 6. (a) The deformed shape of the flexible wall for the steady flow of the FENE-P model fluid in a 2D collapsible channel, compared with the profile for a Newtonian fluid, with Re = 1.0 and β = 0.0071, at various values of bM for R = 45 and Wi = 0.1. (b) Comparison of the narrowest channel gap between the flexible wall and the bottom wall at various values of bM for varying Wi values.

Figure 7. Dependence of the pressure profile along the flexible membrane on Wi for the FENE-P model at various values of bM. Here, the vertical lines with the same color as the pressure profile indicate the position of the narrowest channel gap for the corresponding cases.

Luo and Pedley6 introduced a dimensionless ratio R to represent the influence of membrane tension (R = Ca/Ca*, where Ca* is a reference dimensionless tension (defined using χ = 1.610245 N/m). We keep R fixed at R = 45, and we use the same value of the dimensionless transmural pressure difference, Pd = (pe  pd)W/(η0U0) = 9.3  104, as used by Luo and Pedley6 and in our earlier work.12 To attain a value of Re = 1, which is characteristic of the microcirculation, we used values of F = 1054 kg/m3, U0 = 1.338  102 m/s, W = 102 m, ηp,0 = 0.14 Pa s, ηs = 0.001 Pa s, and η0 = 0.141 Pa s. This yields a viscosity ratio β = 0.0071 (which signifies that the fluid is predominantly elastic). The FENE-P parameter bM is set in the range of 2100.

’ RESULTS AND DISCUSSION A comparison of the predictions of the current formulation for a Newtonian fluid with those of Luo and Pedley6 has been carried out by Chakraborty et al.12 Predictions of the profile of the membrane at different values of R were shown to be in excellent agreement with the reported results of Luo and Pedley.6 Before we present results for the various relevant variables, it is appropriate here to show the dependence of the FENE-P model’s predictions of the microstructure’s contribution to viscosity ηp

Figure 8. Dependence of the pressure drop in the channel on Wi for the FENE-P model at various values of bM. Note that, for a Newtonian fluid, ΔP = 7474.0.

fi, in a steady shear flow, at a on the local Weissenberg number W constant value of λ0 for different values of bM. As is evident from Figure 2, the FENE-P model undergoes increasing shear thinning with decreasing values of the finite extensibility parameter bM. Note that, for the values of U0 and W chosen here, this value of λ0 corresponds to an inlet Weissenberg number of Wi = 0.1. The dependence of ηp on γ_ for a FENE-P model in steady shear flow has been derived previously.21 In this section, we examine the effect of varying bM on other properties such as the limiting Weissenberg number Wi, the shape of the membrane, and pressure profile in the channel. Mesh Convergence and the Limiting Weissenberg Number. Since most viscoelastic computations break down numerically at some limiting value of Wi, because of the development of large stresses, it is essential to ensure mesh convergence over a range of parameters. Three different meshes (M1, M2, and M3), with 2145, 3008, and 3800 elements, respectively, have been considered here for examining the convergence of the numerical computations. A convenient set of variables with which to examine the convergence of numerical computations of viscoelastic fluids are the invariants of the conformation dyadic, M. Previous studies have revealed that, invariably, the breakdown of viscoelastic computations coincides with the smallest eigenvalue becoming negative, and the largest eigenvalue assuming significantly large values in some regions of the flow domain.14,20 The magnitude of 13164

dx.doi.org/10.1021/ie200173b |Ind. Eng. Chem. Res. 2011, 50, 13161–13168

Industrial & Engineering Chemistry Research

ARTICLE

Figure 9. Pressure profile along the flexible membrane for Newtonian fluids with a range of viscosities, and the pressure profile for a FENE-P fluid with Wi = 0.1 and bM = 2, 10, and 100. The values η = 0.0523, 0.0205, and 0.0156 Pa s are the calculated values of the reduced viscosity, respectively, at bM = 100, 10, and 2, at the location of the maximum fi under the flexible membrane. W

Figure 10. Dependence of the axial component of the conformation tensor Mxx on bM for a fixed value of Wi = 0.1.

the maximum in m3 and the minimum in m1, in the entire flow domain, is displayed in Figure 3, as a function of Wi, for the FENEP model, at R = 45 and bM = 100. The value of Wi at which mesh convergence no longer exists can be clearly located in Figure 3a for the minimum eigenvalue m1. As has been observed in previous studies,14,20 the maximum attainable value of Wi increases with mesh refinement; for example, in the FENE-P model at R = 45 and bM = 100, the limiting value of Wi is 0.26 on M2 and 0.27 on M3. Unless otherwise specifically stated, we have used the M3 mesh in all the remaining computations carried out in this work. Figure 4 examines the dependence of the limiting Weissenberg number (which is the value at which the curves terminate), on the finite extensibility parameter bM. Both figures indicate that Wi increases with decreasing bM. Therefore, shear thinning enhances the computability of the viscoelastic fluid by decreasing the value of m3, delaying the value of Wi at which m1 becomes negative. Flow Field and Flexible Membrane Shape. To examine the behavior of a viscoelastic liquid in this flow domain, we compare the velocity contours predicted by a Newtonian fluid with the predictions of the FENE-P model. Figure 5 compares the velocity profile prediction at Wi = 0.1 for different values of bM. It is clear that the Newtonian and FENE-P fluids do not show any

Figure 11. Dependence of the axial component of the conformation tensor Mxx on Wi, for (a) bM = 100, (b) bM = 10, and (c) bM = 2.

significant differences between them for bM = 100 at Wi = 0.1. The prediction of the FENE-P model for lower values of bM, on the other hand, displays some difference from the Newtonian 13165

dx.doi.org/10.1021/ie200173b |Ind. Eng. Chem. Res. 2011, 50, 13161–13168

Industrial & Engineering Chemistry Research

ARTICLE

Figure 12. Dependence of the total tangential shear stress on the membrane (τt þ σt) on bM for a fixed value of Wi = 0.1.

velocity profile, which must arise since different values of bM imply different degrees of shear thinning. In the remainder of this section, we focus our attention on examining the influence of the Weissenberg number Wi and the finite extensibility bM on the shape of the membrane, for the fixed values of the viscosity ratio β, the transmural pressure Pd, and the Reynolds number Re that have been adopted here. The dependence of the shape of the channel on Wi and bM is shown in Figure 6a, which is a zoomed-in view of the membrane shape close to the center of the membrane, for the FENE-P model fluid, at various values of bM, for a fixed value of Wi = 0.1. Shear thinning leads to a departure from the Newtonian shape. The influence of shear thinning and viscoelasticity on the extent of deformation of the membrane can be seen more clearly in Figure 6b, which reveals the effect of Weissenberg number (Wi) on the narrowest channel gap scaled by the value for a Newtonian fluid for various values of bM. The value of the narrowest channel gap is seen to be decreasing with increase in Wi at fixed values of bM, while at a fixed Wi, h/hNewtonian decreases with decreasing bM. Both of these parameters clearly have a similar effect on changes in the membrane shape, since they both have a similar effect on the viscosity of the fluid, as evident from Figure 2. Pressure and Stresses. From the nature of the boundary condition on the flexible membrane adopted here (eq 6), it is clear that the only reason the membrane can change shape as a result of changing either bM or Wi is if there is a change in the normal stress acting on the membrane. Patankar et al.22 analytically showed that the contribution from the elastic and viscous stresses to the component of stress normal to the surface of a rigid body, in an incompressible Oldroyd-B fluid, is zero. A similar result has also been observed by Chakraborty et al.12 in a 2D collapsible channel for the flow of incompressible OldroydB, FENE-P, and Owens model fluids, where they showed that the normal component of stress on the membrane wall is solely due to pressure. This is also confirmed in the present simulations, where we find that the normal component of stress on the wall is solely due to pressure for all values of bM. Figure 7 displays the pressure profile along the flexible membrane for the FENE-P model at various values of bM. It is clear that both increasing the value of Wi and decreasing the value of bM decreases the pressure on the membrane upstream of the narrowest section of the channel. Figure 8 shows the overall

Figure 13. Dependence of the total tangential shear stress on the membrane (τt þ σt) on Wi, for (a) bM = 100, (b) bM = 10, and (c) bM = 2.

pressure drop (ΔP) across the channel for the FENE-P model for different values of bM with varying Wi. The decrease in pressure drop is significant with increase in Wi and with decrease in bM. This decrease in pressure and pressure drop is clearly attributable to the decrease in the viscosity of the fluid in the channel, because of shear thinning. In turn, the decrease 13166

dx.doi.org/10.1021/ie200173b |Ind. Eng. Chem. Res. 2011, 50, 13161–13168

Industrial & Engineering Chemistry Research in pressure is the source of change in membrane shape seen earlier in Figure 6. If shear thinning is the primary cause for the reduced prediction of pressure on the membrane in a channel flow, then a Newtonian fluid with viscosity equal to the effective viscosity in the FENE-P model, would predict a similar value of pressure as fi = λ0γ_ the FENE-P model. We find the maximum value of W (which occurs just below the collapsible membrane) by scanning the flow field in the channel for different values of bM for an inlet Weissenberg number Wi = 0.1. Subsequently, utilizing Figure 2, the corresponding reduced viscosity is obtained at that particular value of bM. In Figure 9, the pressure profile along the flexible membrane predicted by the FENE-P model at Wi = 0.1 for different values of bM has been compared with pressure profiles predicted for Newtonian fluids with a range of viscosity values obtained from the above-mentioned procefi, at Wi = 0.1 and dure. We find that the maximum value of W bM = 100, for a FENE-P fluid flowing in the channel is 44.4. This corresponds to ηp = 0.0513. As can be seen in Figure 9, the pressure profile on the flexible membrane for a Newtonian fluid with η0 = 0.0523 is fairly similar to that for a FENE-P fluid. In particular, the value for the maximum pressure on the membrane is almost identical. The pressure profiles for a FENE-P fluid at other values of bM also show the same trend as the pressure profile of a Newtonian fluid calculated at the respective reduced viscosity. We now explore the influence of viscoelasticity, by examining the behavior of the axial component of the conformation tensor Mxx along the membrane wall. Figure 10 displays the dependence on bM, while Figure 11 displays the dependence on Wi. The value of Mxx is maximum near the narrowest region of the channel and increases with both Wi and bM. This signifies that the FENE-P fluid is experiencing the maximum stretching near the narrowest region of the channel. The value of bM imposes a upper bound to the mean stretchability of the spring in the FENE-P model and restricts the maximum value for Mxx. Figure 12 shows the dependence of the shear stress on bM for the FENE-P fluid at a fixed value of Wi, while Figure 13 explores the effect of Wi at various values of bM. It is immediately apparent from these figures that the common feature is a decrease in tangential shear stress with an increase in Wi and a decrease in bM. This is clearly due to the reduction of viscosity in the channel that accompanies the change in these parameters. The significant differences in the shear stresses predicted by the FENE-P model at different values of bM has no effect on membrane shape under the present set of boundary conditions, since the influence of shear stresses exerted by the fluid on the flexible membrane cannot be taken into account. Future advancement of this work requires a more realistic implementation of boundary conditions.

’ CONCLUSIONS The flow of a FENE-P fluid in a two-dimensional (2D) channel partly bounded by a zero-thickness membrane under constant tension, has been studied numerically with the help of the DEVSS-TG/SUPG mixed finite-element method. Viscoelastic flow in this benchmark geometry for fluidstructure interactions has recently been studied by Chakraborty et al.,12 considering three different viscoelastic fluid models: the Oldroyd-B, the FENE-P, and Owens’ model for blood.13 The aim here has been to examine the influence of shear thinning on the nature of the fluidstructure

ARTICLE

interaction. The most significant conclusions of this work are the following: (1) There is a limiting Weissenberg number Wi at each value of bM for the FENE-P fluid beyond which computations fail. The limiting value of Wi increases as the value of bM decreases. (2) The shape of the deformable membrane, and its dependence on Wi and bM, is entirely determined by the pressure on the membrane surface. (3) The pressure drop, the molecular conformation tensor fields, and the stresses in the flow domain are significantly affected by the extent of shear thinning of the FENE-P fluid. There are many aspects of the fluidstructure interaction problem for viscoelastic flows in 2D collapsible channels that have not been examined in our work, and these are important to examine in order to obtain a more-complete understanding. An immediate requirement is to improve the zero-thickness approximation for the membrane. In the present study, we have found that even though the value of the shear stress on the membrane surface varies significantly with changes in value of finite extensibility parameter and Wi, since the boundary condition accounts only for the normal stress on the membrane, the influence of shear stress on the membrane shape has not been reflected. We are currently developing a model for the collapsible wall that includes the influence of shear stress. In terms of the fluid model, it is important to use a constitutive model that accounts for thixotropy, since this is an important rheological feature of blood flow in the microcirculation, caused by the aggregation of blood cells in regions of low shear rate. Thixotropy leads to rheological properties that depend locally on microstructural dynamics. In terms of improvements in the solution scheme, it would be interesting to see if the use of a log-conformation tensor formalism leads to much higher upper limits to the Wi value. Indeed, most log-conformation tensor simulations, so far, have been restricted to confined flows, and the extension of this technique to free surface flows is a considerable challenge. Flows of Newtonian fluids in collapsible channels at relatively high Reynolds numbers, are well-known for the rich behavior that occurs when the flow is unsteady. This aspect is very worthy of study in the context of viscoelastic fluids.

’ AUTHOR INFORMATION Corresponding Author

*Tel.: þ61 3 9905 3274. Fax: þ61 3 9905 5686. E-mail: ravi. [email protected].

’ ACKNOWLEDGMENT We thank Matheo Pasquali for providing us with his finite element code for simulating coating flows, which we have modified and adapted to this work. This work was supported by an award under the Merit Allocation Scheme on the NCI National Facility at the Australian National University (ANU). The authors also would like to thank the VPAC (Australia), and SUNGRID (Monash University, Australia) for the allocation of computing time on their supercomputing facilities. ’ REFERENCES (1) Sequeira, A.; Janela, J. An Overview of Some Mathematical Models of Blood Rheology. In A Portrait of State-of-the-Art Research at 13167

dx.doi.org/10.1021/ie200173b |Ind. Eng. Chem. Res. 2011, 50, 13161–13168

Industrial & Engineering Chemistry Research

ARTICLE

the Technical University of Lisbon; Springer: Dordrecht, The Netherlands, 2007; pp 6587. (2) Cristini, V.; Kassab, G. S. Computer Modeling of Red Blood Cell Rheology in the Microcirculation: A Brief Overview. Ann. Biomed. Eng. 2005, 33, 1724–1727. (3) Baskurt, O. K.; Meiselman, H. J. Blood Rheology and Hemodynamics. Semin. Thromb. Hemost. 2003, 29, 435–450. (4) Bertram, C. D. Unstable equilibrium behaviour in collapsible tubes. J. Biomech. 1986, 19, 61–69. (5) Lowe, T. W.; Pedley, T. J. Computation of Stokes flow in a channel with a collapsible segment. J. Fluids Struct. 1995, 9, 885–905. (6) Luo, X. Y.; Pedley, T. J. A numerical simulation of steady flow in a 2-D collapsible channel. J. Fluids Struct. 1995, 9, 149–174. (7) Luo, X. Y.; Pedley, T. J. A numerical simulation of unsteady flow in a 2-D collapsible channel. J. Fluid Mech. 1996, 314, 191–225. (8) Rast, M. P. Simultaneous solution of the Navier-Stokes and elastic membrane equations by a finite-element method. Int. J. Num. Methods Fluids 1994, 19, 1115–1135. (9) Jensen, O. E.; Heil, M. High-frequency self-excited oscillations in a collapsible-channel flow. J. Fluid Mech. 2003, 481, 235–268. (10) Liu, H. F.; Luo, X. Y.; Cai, Z. X.; Pedley, T. J. Sensitivity of unsteady collapsible channel flows to modelling assumptions. Commun. Num. Methods Eng. 2009, 25, 483–504. (11) Hazel, A. L.; Heil, M. Steady finite-Reynolds-number flows in three-dimensional collapsible tubes. J. Fluid Mech. 2003, 486, 79–103. (12) Chakraborty, D.; Bajaj, M.; Yeo, L.; Friend, J.; Pasquali, M.; Prakash, J. R. Viscoelastic flow in a two-dimensional collapsible channel. J. Non-Newton. Fluid Mech. 2010, 165, 1204–1218. (13) Owens, R. G. A new microstructure-based constitutive model for human blood. J. Non-Newton. Fluid Mech. 2006, 140, 57–70. (14) Pasquali, M.; Scriven, L. E. Free surface flows of polymer solutions with models based on the conformation tensor. J. Non-Newton. Fluid Mech. 2002, 108, 363–409. (15) Christodoulou, K. N.; Scriven, L. E. Discretization of free surface flows and other moving boundary problems. J. Comput. Phys. 1992, 99, 39–55. (16) de Santos, J. M. Two-phase Cocurrent Downflow through Constricted Passages, Ph.D. Thesis, University of Minnesota, Minneapolis, MN, 1991. (Available from UMI, Ann Arbor, MI, Order No. 9119386.) (17) Benjamin, D. F. Roll Coating Flows and Multiple Roll Systems, Ph.D. Thesis, University of Minnesota, Minneapolis, MN, 1994. (Available from UMI, Ann Arbor, MI, Order No. 9512679.) (18) Xie, X.; Pasquali, M. A new, convenient way of imposing openflow boundary conditions in two- and three-dimensional viscoelastic flows. J. Non-Newton. Fluid Mech. 2004, 122, 159–176. (19) Guenette, R.; Fortin, M. A new mixed finite element method for computing viscoelastic flows. J. Non-Newton. Fluid Mech. 1995, 60, 27–52. (20) Bajaj, M.; Prakash, J. R.; Pasquali, M. A computational study of the effect of viscoelasticity on slot coating flow of dilute polymer solutions. J. Non-Newton. Fluid Mech. 2008, 149, 104–123. (21) Bird, R. B.; Curtiss, C. F.; Armstrong, R. C.; Hassager, O. Dynamics of Polymeric Liquids—Vol. 2: Kinetic Theory, 2nd ed.; John Wiley: New York, 1987. (22) Patankar, N. A.; Huang, P. Y.; Joseph, D. D.; Hu, H. H. Normal Stresses on the Surface of a Rigid Body in an Oldroyd-B Fluid. ASME J. Fluids Eng. 2002, 124, 279–290.

13168

dx.doi.org/10.1021/ie200173b |Ind. Eng. Chem. Res. 2011, 50, 13161–13168