Mechanics Research Communications 35 (2008) 595–602
Contents lists available at ScienceDirect
Mechanics Research Communications journal homepage: www.elsevier.com/locate/mechrescom
Modeling particle segregation in dilute particle/fluid circular pipe flows R. Yu a, J. Peddieson b,*, S. Munukutla b a b
Unilux Boilers, 30 Commerce Park Drive, Niskayuna, New York 12309, USA Department of Mechanical Engineering and Center for Energy Systems Research, Tennessee Technological University, Cookeville, Tennessee 38505, USA
a r t i c l e
i n f o
Article history: Received 8 November 2007 Received in revised form 12 June 2008 Available online 21 June 2008
a b s t r a c t Simulation of the particle segregation phenomenon observed in pipe flows of solid particle/ liquid dispersions is carried out using a two fluid model of dilute fully developed flow. Both analytical and numerical results are reported and a new shear lift model is proposed which is capable of predicting a variety of observed phenomena. Ó 2008 Elsevier Ltd. All rights reserved.
Keywords: Multiphase flow Pipe flow Particle segregation
1. Introduction Chen et al. (2002) recently discussed the use of a dilute two fluid model to simulate near wall segregation of particles in upward turbulent flows of negatively buoyant particle/fluid dispersions in circular pipes. This behavior (also known as the core/annulus phenomenon) has been observed by Lee and Durst (1982) and Tsuji et al. (1984). The present note extends the work of Chen et al. (2002) in several directions using the dilute version of equations describing particle/fluid fully developed circular pipe flows reported by Jean and Peddieson (2001). Issues related to alternate segregation mechanisms, positively buoyant dispersions, laminar flows, and horizontal flows are addressed through both analytical and numerical means. 2. Governing equations For a dilute particle/fluid dispersion, the fluid behavior is unaffected by the presence of particles and, according to Jean and Peddieson (2001), the particle behavior for fully developed flow in a circular pipe of diameter D and average fluid velocity V is governed by axial and radial linear momentum balances which can be written in the respective dimensionless forms
aðr/v0d Þ0 =r þ b/ðvc vd vs Þ ¼ 0 0 0 ðc þ dðvd vc Þ2 Þ/0 þ eðr/v02 d Þ =r þ f/ðvd vc Þvc ¼ 0
ð1Þ
In (1) 2r/D is a radial polar coordinate in the plane of the pipe’s cross section with origin at the axis of symmetry; / is the particle volume fraction; Vvc, Vvd, and Vvs are the respective fluid phase, particle phase, and (constant) settling velocities; a, b, c, d, e, and f are material constants; and a superposed prime denotes differentiation with respect to r. The material constants a–f are respective dimensionless measures of particle phase viscous shear stresses, interphase drag, particle phase pressure, difference between particle average surface pressure and fluid pressure, particle phase viscous normal stresses, and shear induced lift. The first of (1) is the result of an axial linear momentum balance involving particle phase shear * Corresponding author. Tel.: +1 931 372 3274; fax: +1 931 372 6340. E-mail address:
[email protected] (J. Peddieson). 0093-6413/$ - see front matter Ó 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.mechrescom.2008.06.009
596
R. Yu et al. / Mechanics Research Communications 35 (2008) 595–602
stresses (idealized as linearly viscous), interphase drag (idealized as being linear in the phasic relative velocity), and gravity. The second of (1) is the result of a radial linear momentum balance involving a particle phase pressure gradient (with the particle phase pressure idealized as being directly proportional to the volume fraction), particle phase normal stresses (idealized as being quadratically viscous), and a shear induced lift force (idealized in the manner reviewed by Sokolichin et al. (2004)). Complete details of the modeling can be found in Jean and Peddieson (2001). Many competing two fluid models exist for the simulation of particle/fluid two phase flows and a complete review of the literature is not appropriate in this short communication. The reader is, therefore, referred to Jean and Peddieson (2001), Celik and Gel (2002), Zhang and Reese (2003), Sokolichin et al. (2004), Liu and Glasser (2006), and the numerous references contained therein. It is to be observed that no attempt is being made herein to relate the material constants introduced above to the microscopic properties (particle size, shape, etc.) of the dispersion. Many microscopic models exist but all (see above cited references) involve numerous modeling approximations or derive from limited experimental data. The philosophy adopted herein is that the best use for such models is to treat them as guides to the forms of terms appearing in the macroscopic equations but treat the associated coefficients as macroscopic quantities to be measured from macroscopic experiments performed on each particle/fluid dispersion. To make such measurements experimenters must both know what predicted effects to look for and when the inclusion of different terms in a model can produce predictions of the same effect. Parametric studies carried out using macroscopic models can both provide experimenters with information of this kind and indirectly account for the multiplicity of microscopic models available for the estimation of macroscopic coefficients. Measured values for coefficients will enhance the ability of macroscopic models to predict macroscopic phenomena, which the present authors believe to be the proper test for such models. In addition, once measured values of a coefficient are available for several particle/fluid combinations attempts to correlate these values with microscopic properties can be made. Such correlations, if they can be found, will be more reliable than formulas based on highly simplified microscopic models. Some concrete examples of the issues raised in this paragraph will be discussed subsequently. Eq. (1b) can be rearranged to read 0 2 0 02 /0 =/ ¼ ðeðrv02 d Þ =r þ fðvd vc Þvc Þ=ðc þ dðvd vc Þ þ evd Þ
ð2Þ
Substituting (2) into (1a) yields 2 00 02 0 02 0 2 ð1 2ev0d =ðc þ dðvd vc Þ2 þ ev02 d ÞÞvd þ ð1 ðevd þ frðvd vc Þvc Þ=ðc þ dðvd vc Þ þ evd ÞÞvd =r þ s ðvc vd vs Þ ¼ 0
ð3Þ where
s ¼ ðb=aÞ1=2
ð4Þ
With vc and vs known, (3) can be solved for vd. With vd (and, of course, vc and vs) then known, (2) can be solved for /. This uncoupling is an example of analytical work facilitated by the idealized nature of (1) mentioned earlier. All subsequent work is based on (2) and (3). Even without explicit solutions, some information about predicted particle segregation patterns can be inferred from (2). In this (and all subsequent) discussion it will be assumed that c > 0, d > 0. For upward flow it is expected that v0c < 0, v0d < 0. For a negatively buoyant dispersion it is expected that vc vd > 0 over most of the cross section. Then, the tendency of the lift term will be to create central segregation (/0 < 0) for f > 0 and wall segregation (/0 > 0) for f < 0. Similarly, the normal stress term will tend to create central segregation for e > 0 and wall segregation for e < 0. Chen et al. (2002) obtained an expression for vd by fitting to experimental measurements (thus obviating the need to solve an equation corresponding to (3)) and then used the equivalent of (2) with d = e = 0 to predict particle segregation patterns. They found, by fitting their predicted volume fraction expression to experimental data, that the equivalent of a negative f was needed to produce wall segregation. This exercise, incidentally, provides a concrete example of the process of combining macroscopic predictions with experimental observations to estimate a macroscopic coefficient. It is generally believed, however, (see, for instance, Sokolichin et al. (2004)) that f is positive for particles of compact shapes (such as spheres). A positive f would predict central segregation as observed by Young (1960) for laminar flows rather than the wall segregation observed by Lee and Durst (1982) and Tsuji et al. (1984) for turbulent flows. In addition, it can be shown that a positive value of f is needed to predict the wall segregation of positively buoyant solid particles in laminar flows observed by Young (1960), as well as the wall segregation of positively buoyant bubbles in upward gas bubble/liquid pipe flows and central segregation in corresponding downward flows discussed by Lahey (1990). Sokolichin et al. (2004) were so frustrated by the need for f’s of different signs and orders of magnitude to predict different observed phenomena that they refused to recommend the inclusion of lift forces in their review article on two phase gas bubble/ liquid flows. This is a concrete illustration of the ability of comparisons of macroscopic predictions with measurements to identify a significant uncertainty that currently exists in particle/fluid two phase flow modeling. In the next two sections, explicit solutions will be presented which will provide additional insight into the issues discussed heretofore.
R. Yu et al. / Mechanics Research Communications 35 (2008) 595–602
597
3. Analytical solutions Since a major contributor to particle phase stresses is thought to be particle/particle interactions, the assumption of an inviscid particle phase (a = e = 0) would seem plausible for the dilute dispersions under consideration herein. Under these circumstances (since s = 1), (3) has the solution
vd ¼ vc vs
ð5Þ
and (2) has the corresponding solution
/ ¼ /w expðfvs vc =ðc þ dv2s ÞÞ
ð6Þ
where the boundary condition
/ð1Þ ¼ /w
ð7Þ
has been used. As mentioned earlier, Chen et al. (2002) determined an approximate expression for vd by fitting an equation to experimental data. The equation employed expressed vd as a linear function of vc. Eq. (5) is a special case of this and, consequently, (aside from notational differences) (6) is a special case of the particle segregation equation given therein if d = 0. Eq. (6) makes explicit the qualitative trend discussed previously for upward flow (vc > 0) of a negatively buoyant suspension (vs > 0). Central segregation occurs for f > 0 and wall segregation occurs for f < 0. These trends hold for any vc, and thus for both laminar and turbulent fluid flows. Chen et al. (2002) used the equivalent of f < 0 to predict wall segregation. It is generally thought, however, that f > 0 for compact particles, at least in laminar flows. It is possible, of course, that turbulent flow conditions could reverse the sign of f. Intermediate segregation has been observed for neutrally buoyant particles by Young (1960), Segre and Silberberg (1961, 1962). Joseph and Ocando (2002) and Ko et al. (2006) reported on numerical simulations of single particle motion in a plane Poiseuille flow. They concluded that intermediate peaking could also occur for finite buoyancy and that shear lift expressions similar to the last term in (1b) cannot predict this phenomenon because they do not account for the spatial variation of the fluid shear rate. One of the most simplistic ways to allow the prediction of intermediate segregation and related phenomena is to replace the v0c appearing in the last term of (1b) by v0c þ X where X is a constant. Doing this modifies (6) to
/ ¼ /w expðfvs ðvc Xð1 rÞÞ=ðc þ dv2s ÞÞ
ð8Þ
First, consider the form of (8) corresponding to the parabolic velocity profile
vc ¼ 2ð1 r2 Þ
ð9Þ
representative of laminar fluid flow. Substituting (9) into (8) yields
/ ¼ /w expðfvs ð1 rÞð2ð1 þ rÞ XÞ=ðc þ dv2s ÞÞ
ð10Þ
If X = 0, (10) reduces to (9) which exhibits a maximum at r = 0 (central segregation) for f > 0. If 0 < X < 4, (10) exhibits a maximum in the region 0 < r < 1 (intermediate segregation) for f > 0. In particular, X = 2.4 produces the Segre/Silberberg radius r = 0.6. If X P 4, (9) exhibits a maximum at r = 1 (wall segregation) for f > 0. Thus the value of X controls the segregation pattern. Second, consider the form of (8) corresponding to the uniform velocity profile
vc ¼
1 for 0 6 r < 1 0
ð11Þ
for r ¼ 1
which can be regarded as a crude representation of very high Reynolds number turbulent fluid flow. Substituting (11) into (8) produces
( / ¼ /w
expðfvs ð1 Xð1 rÞÞ=ðc þ dv2s ÞÞ for 0 6 r < 1 1
for r ¼ 1
ð12Þ
For any X > 0, (12) exhibits a maximum infinitesimally close to r = 1 (essentially wall segregation) for f > 0. This maximum occurs at r = 1 for X = 1. For any X > 1, the centerline volume fraction is lower than the wall volume for f > 0. For more realistic turbulent velocity profiles it is expected that the maximum volume fraction would occur at a small finite distance from the wall, still approximating wall segregation. Thus armed with the information presented in this paragraph, experimenters would now be in a position to evaluate the proposed shear lift model by (as a first step) determining whether measurements confirm the qualitative features of the predictions. The discussion surrounding (8)–(12) exhibits the possibilities inherent in the crude modification of the shear induced lift representation proposed herein. In this preliminary work the parameter X was treated as a constant. It should probably be regarded as an internal variable which can evolve to different constant values in different situations. It is believed that this approach is deserving for further investigation.
598
R. Yu et al. / Mechanics Research Communications 35 (2008) 595–602
It is apparent from the form of the last term in (1b) that shear induced lift requires both a finite fluid shear rate and a finite phasic relative velocity. Particle segregation is known to occur in horizontal flows, but according to (5) the phasic relative velocity vanishes for horizontal flows (vs = 0). The single particle simulations of Joseph and Ocando (2002) (as well as numerous references contained therein) indicate that particles lag behind the fluid in horizontal flows due to wall effects. In a two fluid continuum model, the only mechanism capable of accounting for particle/wall interactions is particle phase stresses. Thus, even in dilute dispersions where particle/particle interactions are rare, particle phase stresses may be necessary. It is of interest to determine whether this will produce the required phasic relative velocity. A closed form solution to (2) and (3) is possible for the special case of e = 0 if the idealized fluid turbulent velocity profile (11) is employed. The corresponding forms of (2) and (3) are
/0 ¼ 0;
v00d þ v0d =r s2 vd ¼ s2 ðvs 1Þ
ð13Þ
which have the solutions
/ ¼ /w ;
vd ¼ ð1 vs Þð1 Io ðsrÞ=ðIo ðsÞ þ gsI1 ðsÞÞÞ
ð14Þ
subject to (7) and the boundary condition
vd ð1Þ ¼ gv0d ð1Þ
ð15Þ
(g being a dimensionless positive particle phase wall slip parameter). In (14) I denotes a modified Bessel function of the first kind with the subscript indicating the order. Using (11) and (14) it can be shown that
vr ¼ vd vc ¼
ðF þ vs ð1 FÞÞ
for 0 6 r < 1
ð1 vs Þð1 Fð1ÞÞ for r ¼ 1
ð16Þ
where
F ¼ Io ðsrÞ=ðIo ðsÞ þ gsI1 ðsÞÞ
ð17Þ
From (17) it is clear that
0 6 F 6 1;
061F 61
ð18Þ
Thus, for vs > 1 one finds vr 6 0 for 0 6 r 6 1 while for vs < 1 one finds vr 6 0 for 0 6 r < 1 and vr P 0 for r = 1 (with the transition from negative to positive relative velocity occurring infinitesimally close to the wall). More realistic turbulent fluid velocity profiles are expected to produce qualitatively similar results (with the transition from negative to positive relative velocity occurring at a small finite distance from the wall). This will be demonstrated in the next section. The above discussion shows that particle phase viscosity is an appropriate mechanism for the production of the phasic relative velocity needed to produce shear induced lift. One would expect wall effects to be localized which could call into question the use of a constant particle phase viscosity. However the use of a small particle phase viscosity produces
: F¼ expðsð1 rÞÞ=ð1 þ gsÞ
ð19Þ
which is localized in the vicinity of the wall. Thus, at least in this case, the use of a constant particle phase viscosity does not produce unphysical behavior. Recognizing that particle phase viscosity (wall effects) will produce a relative velocity in horizontal flows, it is to be expected that the discussion of particulate segregation carried out previously in connection with (12) would remain qualitatively valid for horizontal flows (with vs now playing the role of an average relative velocity). To investigate this issue and the effect of normal stresses requires numerical computation. This will be discussed in the next section.
4. Numerical solutions This section reports numerical solutions of (2) and (3). A standard iterative finite difference methodology based on three point central difference representation of derivatives and solution of tridiagonal linear algebraic equations by Potter’s method at each iteration was used to determine vd from (3). Then with vd known, (2) was solved by the trapezoidal rule to find /. For simplicity attention is restricted to horizontal flows (vs = 0). As discussed in the previous section upward flows with negative buoyancy will be qualitatively similar. Following Chen et al. (2002), the power law representation of the fluid phase velocity profile
vc ¼ ð2m þ 1Þðm þ 1Þðl rÞ1=m =ð2m2 Þ
ð20Þ
is used as the basis for the numerical computations. The form of (20) produces an infinite slope of the fluid velocity profile at the pipe wall. This causes both inaccuracy in (20) as a representation of the fluid velocity profile near the wall and numerical difficulties in the solution of (2) and (3). It was decided to deal with the numerical problem by replacing (20) with
vc ¼ ð2m þ 1Þðm þ 1Þð1 ð1 dt ÞrÞ1=m =ð2m2 Þ
ð21Þ
R. Yu et al. / Mechanics Research Communications 35 (2008) 595–602
599
where
dt ¼ 168ðlogðRÞÞ1:25 =R
ð22Þ
is the sum of the dimensionless thicknesses of the viscous sublayer and the adjacent buffer layer that exist in the immediate vicinity of the wall. In (22)
R ¼ VD=mc
ð23Þ
(mc being the fluid phase kinematic viscosity) is the Reynolds number. Daily and Harleman’s (1966) expression for dt in terms of the Reynolds number and the friction factor and a standard relation between the friction factor and the Reynolds number (see, for instance, White, 1974) have been combined to produce (22) which determines dt if R is regarded as known. Since the R also determines m (see, for instance, Daily and Harleman, 1966), specification of R renders (21) determinate and thus gives an analytical representation of the fluid velocity which exhibits a fictitious wall slip velocity and a finite wall slope. The inequality dt 1 is satisfied for turbulent pipe flow and, therefore, the use of (21) will not seriously impair the accuracy of the predicted fluid velocity in the core region (everywhere but the immediate vicinity of the wall). Several preliminary simulations demonstrated that the qualitative behavior of the particulate phase of a dilute suspension with turbulent fluid flow was not strongly dependent on the Reynolds number. The simulations to be discussed will, therefore, focus on the typical case of R = 105 (m = 7) used by Chen et al. (2002). Figs. 1 and 2 report predictions which illustrate the effect the wall slip modulus g. In all cases the volume fraction variation is largest near the wall. It is reasonable to expect this in turbulent flows because the fluid shear rate is extremely high in the area adjacent to the wall. In these simulations particle segregation is produced by shear induced lift only. For no wall slip (g = 0) the relative velocity is negative everywhere and central segregation is predicted for positive f (in qualitative agreement with (6) for an inviscid particle phase). For perfect wall slip (g = 1), on the other hand, the relative velocity is positive near the wall (as suggested by (16)) and this produces wall segregation in upward flow. It can be seen that both the particulate velocity and volume fraction are quite sensitive to the values of g (especially for g 6 10). For f > 0, a transition from central segregation to wall segregation occurs as g increases. Figs. 3 and 4 illustrate the capability of the particle phase normal stresses to create wall segregation for e < 0. Both the particulate velocity and volume fraction exhibit significant sensitivity to the value of the slip parameter g with a transition from significant wall segregation to slight intermediate segregation occurring as g increases. Numerical computation is needed to determine the effect of particle phase normal stresses since the analytical work discussed in Section 3 did not yield useful results for e 6¼ 0. The fact that both shear lift and particle phase normal stress effects can predict wall segregation is useful in the interpretation of measurements. It tells experimenters that an ambiguity exists as to what is actually being measured and that this must be addressed.
Fig. 1. Velocity profiles (a = 1, b = 5, c = 1, d = 1, e = 0, f = 1, /w = 0.01, R = 105).
600
R. Yu et al. / Mechanics Research Communications 35 (2008) 595–602
Fig. 2. Volume fraction profiles (a = 1, b = 5, c = 1, d = 1, e = 0, f = 1, /w = 0.01, R = 105).
Fig. 3. Velocity profiles (a = 0.5, b = 1, c = 1, d = 1, e = 8, f = 0.5, /w = 0.01, R = 105).
Additional simulations (not presented graphically for the sake of brevity) were performed to determine the influence of the drag modulus b on predictions. As expected, increasing b caused velocity equilibrium to be approached and produced a particulate velocity boundary layer. Outside the boundary layer a uniform volume fraction was approached. The influence of
R. Yu et al. / Mechanics Research Communications 35 (2008) 595–602
601
Fig. 4. Volume fraction profiles (a = 0.5, b = 1, c = 1, d = 1, e = 8, f = 0.5, /w = 0.01, R = 105).
b was more significant for zero particulate wall slip than for infinite particulate wall slip. These trends are completely consistent with those exhibited by the idealized closed form solution (14). 5. Conclusion In the foregoing, a mathematical model of steady fully developed pipe flow of a dilute two phase particulate suspension was employed to investigate the prediction of particle segregation in fully developed circular pipe flows. Approximate closed form solutions were obtained based on simplified forms of the governing equations. In addition, the full governing equations were solved numerically using a finite difference scheme. Simulations were carried out for a variety of parametric combinations. A representative sample of the predicted results was presented graphically and used to demonstrate interesting aspects of the model predictions. Some significant conclusions reached on the basis of this work are as follows. First, the standard shear induced lift force produces central segregation in vertically upward and horizontal flows of negatively buoyant suspensions. Second, a modification of the shear induced lift force expression suggested by the work of Joseph and Ocando (2002) appears capable of predicting a variety of segregation patterns including central, wall, and intermediate segregation. Third boundary layers are exhibited for small values of the particle phase viscosity. The exact solutions for this situation are close to corresponding solutions for an inviscid particulate phase except that a small deviation in particulate velocity appeared in the boundary layers. An inviscid particulate phase model is, therefore, of quantitative utility for low particulate phase viscosity suspensions. Fourth, particulate wall slip appears to have an important effect on predictions. The slip model used herein is tentative. Thus, future research is needed to deal with this important issue in a more quantitative way. Fifth, the viscous particulate normal stress effect is capable of producing wall segregation. Sixth, predictions of dilute suspension flows exhibiting fluid turbulence based on realistic fluid velocity profiles were qualitatively similar to those based on a slug representation of the fluid flow. References Celik, I., Gel, A., 2002. Flow Turb. Combs. 68, 289. Chen, X., Hsu, C., Chiew, Y., 2002. ASCE J. Eng. Mech. 128, 487. Daily, J., Harleman, D., 1966. Fluid Dynamics. Adison-Wesley, Reading, MA. Jean, T., Peddieson, J., 2001. Int. J. Eng. Sci. 39, 1167. Joseph, D., Ocando, D., 2002. J. Fluid Mech. 454, 263. Ko, T., Patankar, N., Joseph, D., 2006. Compt. Fluids 35, 121. Lahey, R., 1990. Nucl. Eng. Design 122, 17. Lee, S., Durst, F., 1982. Int. J. Multiphase Flow 15, 735. Liu, X., Glasser, B., 2006. AIChE J. 52, 940.
602
R. Yu et al. / Mechanics Research Communications 35 (2008) 595–602
Segre, G., Silberberg, A., 1961. Nature 189, 209. Segre, G., Silberberg, A., 1962. J. Fluid Mech. 14, 115. Sokolichin, A., Eigenberger, G., Lapin, A., 2004. AIChE J. 50, 24. Tsuji, Y., Morikawa, Y., Shiomi, H., 1984. J. Fluid Mech. 139, 417. White, F., 1974. Viscous Fluid Flows. McGraw-Hill, New York. Young, D., 1960. The Coring Phenomena in the Flow of Suspensions in Vertical Tubes, ASME Preprint 60-HYD-12. Zhang, Y., Reese, J., 2003. AIChE J. 49, 3048.