06.0 PP 30 59 Single-Particle Motion

2
Single-Particle Motion
2.1 Introduction
This chapter briefly reviews the issues and problems involved in constructing the
equations of motion for individual particles, drops, or bubbles moving through a
fluid. For convenience we use the generic name particle to refer to the finite pieces
of the disperse phase or component. The analyses are implicitly confined to those
circumstances in which the interactions between neighboring particles are negligible.
In very dilute multiphase flows in which the particles are very small compared with
the global dimensions of the flow and are very far apart compared with the particle
size, it is often sufficient to solve for the velocity and pressure, u i (xi , t) and p(xi , t), of
the continuous suspending fluid while ignoring the particles or disperse phase. Given
this solution one could then solve an equation of motion for the particle to determine
its trajectory. This chapter focuses on the construction of such a particle or bubble
equation of motion.
The body of fluid mechanical literature on the subject of flows around particles or
bodies is very large indeed. Here we present a summary that focuses on a spherical
particle of radius R and employs the following common notation. The components of
the translational velocity of the center of the particle is denoted by Vi (t). The velocity
that the fluid would have had at the location of the particle center in the absence of
the particle is denoted by Ui (t). Note that such a concept is difficult to extend to the
case of interactive multiphase flows. Finally, the velocity of the particle relative to the
fluid is denoted by Wi (t) = Vi − Ui .
Frequently the approach used to construct equations for Vi (t) (or Wi (t)) given
Ui (xi , t) is to individually estimate all the fluid forces acting on the particle and to
equate the total fluid force, Fi , to m p d Vi /dt (where m p is the particle mass, assumed
constant). These fluid forces may include forces due to buoyancy, added mass, drag,
and so on. In the absence of fluid acceleration (dUi /dt = 0) such an approach can
be made unambiguously; however, in the presence of fluid acceleration, this kind of
heuristic approach can be misleading. Hence we concentrate in the next few sections
on a fundamental fluid mechanical approach that minimizes possible ambiguities. The
30
https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

2.2 Flows Around a Sphere 31
Figure 2.1. Notation for a spherical particle.
classical results for a spherical particle or bubble are reviewed first. The analysis is
confined to a suspending fluid that is incompressible and Newtonian so that the basic
equations to be solved are the continuity equation.
∂u j
=0 (2.1)
∂x j
and the Navier–Stokes equations

∂u i ∂u i ∂p ∂ 2ui
ρC + uj =− + ρC νC (2.2)
∂t ∂x j ∂ xi ∂x j∂x j
where ρC and νC are the density and kinematic viscosity of the suspending fluid. It is
assumed that the only external force is that due to gravity, g. Then the actual pressure
is p = p − ρC gz, where z is a coordinate measured vertically upward. Furthermore,
to maintain clarity we confine our attention to rectilinear relative motion in a direction
conveniently chosen to be the x1 direction.
2.2 Flows Around a Sphere
2.2.1 At High Reynolds Number

For steady flows about a sphere in which dUi /dt = d Vi /dt = dWi /dt = 0, it is con-
venient to use a coordinate system, xi , fixed in the particle as well as polar coordinates
(r, θ ) and velocities u r , u θ as defined in Figure 2.1.
Then Eqs. (2.1) and (2.2) become the following:
1 ∂ 2 1 ∂
(r u r ) + (u θ sin θ) = 0 (2.3)
r ∂r
2 r sin θ ∂θ
and

∂u r ∂u r u θ ∂u r u 2θ
ρC + ur + − (2.4)
∂t ∂r r ∂θ r
∂p 1 ∂ ∂u r 1 ∂ ∂u r 2u r 2 ∂u θ
=− + ρC νC 2 r2 + 2 sin θ − 2 − 2
∂r r ∂r ∂r r sin θ ∂θ ∂θ r r ∂θ

32 Single-Particle Motion

∂u θ ∂u θ u θ ∂u θ ur u θ 1 ∂p
ρC + ur + + =− (2.5)
∂t ∂r r ∂θ r r ∂θ
1 ∂ ∂u θ 1 ∂ ∂u θ 2 ∂u r uθ
+ ρC νC 2 r2 + 2 sin θ + 2 − .
r ∂r ∂r r sin θ ∂θ ∂θ r ∂θ r sin2 θ
2
The Stokes streamfunction, ψ, is defined to satisfy continuity automatically:

1 ∂ψ 1 ∂ψ
ur = ; uθ = − (2.6)
r2 sin θ ∂θ r sin θ ∂r
and the inviscid potential flow solution is as follows:
Wr2 2 D
ψ =− sin θ − sin2 θ (2.7)
2 r
2D
u r = −W cos θ − 3 cos θ (2.8)
r
D
u θ = +W sin θ − 3 sin θ (2.9)
r
D
φ = −W r cos θ + 2 cos θ, (2.10)
r
where, because of the boundary condition (u r )r =R = 0, it follows that D = −W R 3 /2.
In potential flow one may also define a velocity potential, φ, such that u i = ∂φ/∂ xi .
The classic problem with such solutions is the fact that the drag is zero, a circumstance
termed D’Alembert’s paradox. The flow is symmetric about the x2 x3 plane through
the origin and there is no wake.
The real viscous flows around a sphere at large Reynolds numbers, Re =
2W R/νC > 1, are well documented. In the range from ∼103 to 3 × 105 , laminar
boundary layer separation occurs at θ ∼ = 84◦ and a large wake is formed behind the
sphere (see Figure 2.2). Close to the sphere the near-wake is laminar; further down-
stream transition and turbulence occurring in the shear layers spreads to generate a
turbulent far-wake. As the Reynolds number increases the shear layer transition moves
forward until, quite abruptly, the turbulent shear layer reattaches to the body, resulting
in a major change in the final position of separation (θ ∼= 120◦ ) and in the form of the
turbulent wake (Figure 2.2). Associated with this change in flow pattern is a dramatic
decrease in the drag coefficient, CD (defined as the drag force on the body in the
negative x1 direction divided by 12 ρC W 2 π R 2 ), from a value of ∼0.5 in the laminar
separation regime to a value of ∼0.2 in the turbulent separation regime (Figure 2.3). At
values of Re less than ∼103 the flow becomes quite unsteady with periodic shedding
of vortices from the sphere.
2.2.2 At Low Reynolds Number
At the other end of the Reynolds number spectrum is the classic Stokes solution for
flow around a sphere. In this limit the terms on the left-hand side of Eq. (2.2) are

Figure 2.2. Smoke visualization of the nom-

inally steady flows (from left to right) past
a sphere showing, at the top, laminar sepa-
ration at Re = 2.8 × 105 and, on the bottom,
turbulent separation at Re = 3.9 × 105 . Pho-
tographs by F.N.M. Brown, reproduced with
the permission of the University of Notre
Dame.
neglected and the viscous term is retained. This solution has the following form:

Wr2 A
ψ = sin θ −
2
+ + Br (2.11)
2 r
2A 2B
u r = cos θ −W + 3 + (2.12)
r r
Figure 2.3. Drag coefficient on a sphere

as a function of Reynolds number.
Dashed curves indicate the drag crisis
regime in which the drag is very sen-
sitive to other factors such as the free
stream turbulence.


A B
u θ = − sin θ −W − 3 + , (2.13)
r r
where A and B are constants to be determined from the boundary conditions on the
surface of the sphere. The force, F, on the particle in the x1 direction is as follows:

4 4W 8A 2B
F1 = π R ρC νC −
2
+ 4 + 2 . (2.14)
3 R R R
Several subcases of this solution are of interest in the present context. The first is
the classic Stokes (1851) solution for a solid sphere in which the no-slip boundary
condition, (u θ )r =R = 0, is applied [in addition to the kinematic condition (u r )r =R = 0].
This set of boundary conditions, referred to as the Stokes boundary conditions, leads
to the following:
W R3 3W R
A=− , B=+ , and F1 = −6πρC νC W R. (2.15)
4 4
The second case originates with Hadamard (1911) and Rybczynski (1911) who sug-
gested that, in the case of a bubble, a condition of zero shear stress on the sphere
surface would be more appropriate than a condition of zero tangential velocity, u θ .
Then it transpires that
WR
A = 0, B=+ , and F1 = −4πρC νC W R. (2.16)
2
Real bubbles may conform to either the Stokes or Hadamard–Rybczynski solutions
depending on the degree of contamination of the bubble surface, as discussed in more
detail in section 3.3. Finally, it is of interest to observe that the potential flow solution
given in Eqs. (2.7) to (2.10) is also a subcase with
W R3
A=+ , B = 0, and F1 = 0. (2.17)
2
However, another paradox, known as the Whitehead paradox, arises when the valid-
ity of these Stokes flow solutions at small (rather than zero) Reynolds numbers is
considered. The nature of this paradox can be demonstrated by examining the mag-
nitude of the neglected term, u j ∂u i /∂ x j , in the Navier–Stokes equations relative to
the magnitude of the retained term νC ∂ 2 u i /∂ x j ∂ x j . As is evident from Eq. (2.11),
far from the sphere the former is proportional to W 2 R/r 2 , whereas the latter be-
haves like νC W R/r 3 . It follows that although the retained term dominates close to
the body (provided Re = 2W R/νC 1), there will always be a radial position, rc ,
given by R/rc = Re beyond which the neglected term will exceed the retained vis-
cous term. Hence, even if Re 1, the Stokes solution is not uniformly valid. Rec-
ognizing this limitation, Oseen (1910) attempted to correct the Stokes solution by
retaining in the basic equation an approximation to u j ∂u i /∂ x j that would be valid in

the far field, −W ∂u i /∂ x1 . Thus the Navier–Stokes equations are approximated by the
following:
∂u i 1 ∂p ∂ 2ui
−W =− + νC . (2.18)
∂ x1 ρC ∂ x i ∂x j∂x j
Oseen was able to find a closed form solution to this equation that satisfies the Stokes
boundary conditions approximately as follows:

2 r sin θ
2 2
R sin2 θ 3νC (1 + cos θ) Wr

ψ = −W R + + 1 − e 2νC (1−cos θ ) , (2.19)
2R 2 4r 2W R
which yields a drag force as follows:

3
F1 = −6πρC νC W R 1 + Re . (2.20)
16
It is readily shown that Eq. (2.19) reduces to Eq. (2.11) as Re → 0. The correspond-
ing solution for the Hadamard–Rybczynski boundary conditions is not known to the
author; its validity would be more questionable because, unlike the case of Stokes’
boundary conditions, the inertial terms u j ∂u i /∂ x j are not identically zero on the
surface of the bubble.
Proudman and Pearson (1957) and Kaplun and Lagerstrom (1957) showed that
Oseen’s solution is, in fact, the first term obtained when the method of matched
asymptotic expansions is used in an attempt to patch together consistent asymptotic
solutions of the full Navier–Stokes equations for both the near field close to the sphere
and the far field. They also obtained the next term in the expression for the drag
force:

3 9 Re
F1 = −6πρC νC W R 1 + Re + 2
Re ln + 0(Re ) . (2.21)
2
16 160 2
The additional term leads to an error of 1% at Re = 0.3 and does not, therefore, have
much practical consequence.
The most notable feature of the Oseen solution is that the geometry of the stream-
lines depends on the Reynolds number. The downstream flow is not a mirror image
of the upstream flow as in the Stokes or potential flow solutions. Indeed, closer exam-
ination of the Oseen solution reveals that, downstream of the sphere, the streamlines
are further apart and the flow is slower than in the equivalent upstream location.
Furthermore, this effect increases with Reynolds number. These features of the Oseen
solution are entirely consistent with experimental observations and represent the initial
development of a wake behind the body.
The flow past a sphere at Reynolds numbers between ∼ 0.5 and several thousand
has proven intractable to analytical methods though numerical solutions are numerous.
Experimentally, it is found that a recirculating zone (or vortex ring) develops close to
the rear stagnation point at about Re = 30 (see Taneda 1956 and Figure 2.4). With

Figure 2.4. Streamlines of steady flow (from left to right) past a sphere at various Reynolds numbers
(from Taneda 1956, reproduced by permission of the author).
further increase in the Reynolds number this recirculating zone or wake expands.
Defining locations on the surface by the angle from the front stagnation point, the
separation point moves forward from ∼130◦ at Re = 100 to ∼115◦ at Re = 300.
In the process the wake reaches a diameter comparable to that of the sphere when
Re ≈ 130. At this point the flow becomes unstable and the ring vortex that makes up

the wake begins to oscillate (Taneda 1956). However, it continues to be attached to

the sphere until ∼Re = 500 (Torobin and Gauvin 1959).
At Reynolds numbers above about 500, vortices begin to be shed and then convected
downstream. The frequency of vortex shedding has not been studied as extensively
as in the case of a circular cylinder and seems to vary more with Reynolds number.
In terms of the conventional Strouhal number, Str, defined as follows:
Str = 2 f R/W (2.22)
the vortex shedding frequencies, f , that Moller (1938) observed correspond to a

range of Str varying from 0.3 at Re = 1000 to ∼ 1.8 at Re = 5000. Furthermore, as
Re increases above 500 the flow develops a fairly steady near-wake behind which
vortex shedding forms an unsteady and increasingly turbulent far-wake. This process
continues until, at a value of Re of the order of 1000, the flow around the sphere
and in the near-wake again becomes quite steady. A recognizable boundary layer has
developed on the front of the sphere and separation settles down to a position about 84◦
from the front stagnation point. Transition to turbulence occurs on the free shear layer
(which defines the boundary of the near-wake) and moves progressively forward as the
Reynolds number increases. The flow is similar to that of the top picture in Figure 2.2.
Then the events described in the previous section occur with further increase in the
Reynolds number.
Because the Reynolds number range between 0.5 and several hundred can often
pertain in multiphase flows, one must resort to an empirical formula for the drag force
in this regime. A number of empirical results are available; for example, Klyachko
(1934) recommends the following:
2

Re 3
F1 = −6πρC νC W R 1 + , (2.23)
6
which fits the data fairly well up to Re ≈ 1000. At Re = 1 the factor in the square
brackets is 1.167, whereas the same factor in Eq. (2.20) is 1.187. Conversely, at
Re = 1000, the two factors are respectively 17.7 and 188.5.
2.2.3 Molecular Effects
When the mean free path of the molecules in the surrounding fluid, λ, becomes
comparable with the size of the particles, the flow will clearly deviate from the
continuum models, which are relevant only when λ R. The Knudsen number,
Kn = λ/2R, is used to characterize these circumstances, and Cunningham (1910)
showed that the first-order correction for small but finite Knudsen number leads to
an additional factor, (1 + 2AKn), in the Stokes drag for a spherical particle. The nu-
merical factor, A, is roughly a constant of order unity (see, for example, Green and
Lane 1964).

When the impulse generated by the collision of a single fluid molecule with the
particle is large enough to cause significant change in the particle velocity, the resulting
random motions of the particle are called Brownian motion (Einstein 1956). This leads
to diffusion of solid particles suspended in a fluid. Einstein showed that the diffusivity,
D, of this process is given by the following:
D = kT /6πµC R, (2.24)
where k is Boltzmann’s constant. It follows that the typical rms displacement of the
1
particle in a time, t, is given by (kT t/3πµC R) 2 . Brownian motion is usually significant
only for micron- and sub-micron-sized particles. The example quoted by Einstein is
that of a 1-µm-diameter particle in water at 17◦ C for which the typical displacement
during 1s is 0.8 µm.
A third, related, phenomenon is the response of a particle to the collisions of
molecules when there is a significant temperature gradient in the fluid. Then the
impulses imparted to the particle by molecular collisions on the hot side of the particle
will be larger than the impulses on the cold side. The particle will therefore experience
a net force diving it in the direction of the colder fluid. This phenomenon is known
as thermophoresis (see, for example, Davies 1966). A similar phenomenon known as
photophoresis occurs when a particle is subjected to nonuniform radiation. One could
include in this list the Bjerknes forces described in section 3.4 because they constitute
sonophoresis, namely forces acting on a particle in a sound field.
2.3 Unsteady Effects
2.3.1 Unsteady Particle Motions

Having reviewed the steady motion of a particle relative to a fluid, we must now
consider the consequences of unsteady relative motion in which the particle, the fluid,
or both are accelerating. The complexities of fluid acceleration are delayed until the
next section. First we consider the simpler circumstance in which the fluid either
is at rest or has a steady uniform streaming motion (U = constant) far from the
particle. Clearly the second case is readily reduced to the first by a simple Galilean
transformation and it is assumed that this has been accomplished.
In the ideal case of unsteady inviscid potential flow, it can then be shown by using
the concept of the total kinetic energy of the fluid that the force on a rigid particle in
an incompressible flow is given by Fi as follows:
dVj
Fi = −Mi j , (2.25)
dt
where Mi j is called the added mass matrix (or tensor), although the name induced
inertia tensor used by Batchelor (1967) is, perhaps, more descriptive. The reader is
referred to Sarpkaya and Isaacson (1981), Yih (1969), or Batchelor (1967) for detailed

2.3 Unsteady Effects 39
descriptions of such analyses. The above-mentioned methods also show that Mi j for
any finite particle can be obtained from knowledge of several steady potential flows.
In fact,
ρC
Mi j = u ik u jk d(volume), (2.26)
2 volume of fluid
where the integration is performed over the entire volume of the fluid. The velocity
field, u i j , is the fluid velocity in the i direction caused by the steady translation of the
particle with unit velocity in the j direction. Note that this means that Mi j is necessarily
a symmetric matrix. Furthermore, it is clear that particles with planes of symmetry will
not experience a force perpendicular to that plane when the direction of acceleration
is parallel to that plane. Hence if there is a plane of symmetry perpendicular to the k
direction, then for i = k, Mki = Mik = 0, and the only off-diagonal matrix elements
that can be nonzero are Mi j , j = k, i = k. In the special case of the sphere all the
off-diagonal terms will be zero.
Tables of some available values of the diagonal components of Mi j are given by
Sarpkaya and Isaacson (1981), who also summarize the experimental results, partic-
ularly for planar flows past cylinders. Other compilations of added mass results can
be found in Kennard (1967), Patton (1965), and Brennen (1982). Some typical values
for three-dimensional particles are listed in Table 2.1. The uniform diagonal value
for a sphere (often referred to simply as the added mass of a sphere) is 2ρC π R 3 /3 or
one-half the displaced mass of fluid. This value can readily be obtained from Eq. (2.26)
using the steady flow results given in Eqs. (2.7) to (2.10). In general, of course, there
is no special relation between the added mass and the displaced mass. Consider, for
example, the case of the infinitely thin plate or disk with zero displaced mass that has
a finite added mass in the direction normal to the surface. Finally, it should be noted
that the literature contains little, if any, information on off-diagonal components of
added mass matrices.
Now consider the application of these potential flow results to real viscous flows at
high Reynolds numbers (the case of low Reynolds number flows is discussed in section
2.3.4). Significant doubts about the applicability of the added masses calculated from
potential flow analysis would be justified because of the experience of D’Alembert’s
paradox for steady potential flows and the substantial difference between the stream-
lines of the potential and actual flows. Furthermore, analyses of experimental results
requires the separation of the added mass forces from the viscous drag forces. Usually
this is accomplished by heuristic summation of the two forces so that
dVj 1
Fi = −Mi j
− ρC ACi j |V j |V j , (2.27)
dt 2
where Ci j is a lift and drag coefficient matrix and A is a typical cross-sectional area
for the body. This is known as Morison’s equation (see Morison et al. 1950).
Actual unsteady high Reynolds number flows are more complicated and not nec-
essarily compatible with such simple superposition. This is reflected in the fact that

Table 2.1. Added masses (diagonal terms in Mi j ) for some three-dimensional bodies
(particles)
Particle Matrix element Value

Sphere (T) Mii 2
ρ π R3
3 C
Disk(T) M11 8
ρ R3
3 C
Ellipsoids (T) Mii = K ii 34 ρC πab2 a/b K 11 K 22 (K 33 )

2 0.209 0.702
5 0.059 0.895
10 0.021 0.960

3 R3
Sphere near plane Mii = K ii 34 ρC π R 3 K 11 = 1
2
1+ 8 H3
+ ···

wall (T) 3 R3
K 22 = 1
1+ + ···
(R/H 1) 2 16 H 3
K 33 = K 22
Sphere near free Mii = K ii 34 ρC π R 3 H/R K 11

surface (E)
8.0 0.52
(R/H 1)
4.0 0.59
2.0 0.54
1.0 0.44
0.0 0.25
Note: (T), Potential flow calculations; (E), experimental data from Patton (1965).
the coefficients Mi j and Ci j appear from the experimental results to be functions not
only of Re but also of the reduced time or frequency of the unsteady motion. Typically,
experiments involve either oscillation of a body in a fluid or acceleration from rest.
The most extensively studied case involves planar flow past a cylinder (for exam-
ple, Keulegan and Carpenter 1958), and a detailed review of this data is included in
Sarpkaya and Isaacson (1981). For oscillatory motion of the cylinder with velocity
amplitude UM and period t ∗ , the coefficients are functions of both the Reynolds
number, Re = 2UM R/νC , and the reduced period or Keulegan–Carpenter number,
Kc = UM t ∗ /2R. When the amplitude, UM t ∗ , is less than about 10R(Kc < 5), the in-
ertial effects dominate and Mii is only a little less than its potential flow value over
a wide range of Reynolds numbers (104 < Re < 106 ). However, for larger values of

Kc, Mii can be substantially smaller than this and, in some range of Re and Kc, may
actually be negative. The values of Cii (the drag coefficient) that are deduced from
experiments are also a complicated function of Re and Kc. The behavior of the co-
efficients is particularly pathological when the reduced period, Kc, is close to that of
vortex shedding (Kc of the order of 10). Large transverse or lift forces can be gen-
erated under these circumstances. To the author’s knowledge, detailed investigations
of this kind have not been made for a spherical body, but one might expect the same
qualitative phenomena to occur.
2.3.2 Effect of Concentration on Added Mass

Although most multiphase flow effects are delayed until later chapters it is conve-
nient at this point to address the issue of the effect on the added mass of the par-
ticles in the surrounding mixture. It is to be expected that the added mass coeffi-
cient for an individual particle would depend on the void fraction of the surrounding
medium. Zuber (1964) first addressed this issue using a cell method and found that
the added mass, Mii , for spherical bubbles increased with volume fraction, α, as
follows:
Mii (α) (1 + 2α)
= = 1 + 3α + O(α 2 ). (2.28)
Mii (0) (1 − α)
The simplistic geometry assumed in the cell method (a concentric spherical shell of
fluid surrounding each spherical particle) caused later researchers to attempt improve-
ments to Zuber’s analysis; for example, van Wijngaarden (1976) used an improved
geometry (and the assumption of potential flow) to study the O(α) term and found
that
Mii (α)
= 1 + 2.76α + O(α 2 ), (2.29)
Mii (0)
which is close to Zuber’s result. However, even more accurate and more recent analyses
by Sangani et al. (1991) have shown that Zuber’s original result is, in fact, remarkably
accurate even up to volume fractions as large as 50% (see also Zhang and Prosperetti
1994).
2.3.3 Unsteady Potential Flow
In general, a particle moving in any flow other than a steady uniform stream will
experience fluid accelerations, and it is therefore necessary to consider the structure
of the equation governing the particle motion under these circumstances. Of course,
this will include the special case of acceleration of a particle in a fluid at rest (or with
a steady streaming motion). As in the earlier sections we confine the detailed solu-
tions to those for a spherical particle or bubble. Furthermore, we consider only those
circumstances in which both the particle and fluid acceleration are in one direction,

chosen for convenience to be the x1 direction. The effect of an external force field such
as gravity is omitted; it can readily be inserted into any of the solutions that follow by
the addition of the conventional buoyancy force.
All the solutions discussed are obtained in an accelerating frame of reference fixed
in the center of the fluid particle. Therefore, if the velocity of the particle in some
original, noninertial coordinate system, xi∗ , was V (t) in the x1∗ direction, the Navier–
Stokes equations in the new frame, xi , fixed in the particle center are as follows:
∂u i ∂u i 1 ∂P ∂ 2ui
+ uj =− + νC , (2.30)
∂t ∂x j ρC ∂ xi ∂x j∂x j
where the pseudo-pressure, P, is related to the actual pressure, p, by the following:
dV
P = p + ρC x 1
. (2.31)
dt
Here the conventional time derivative of V (t) is denoted by d/dt, but it should be
noted that in the original xi∗ frame it implies a Lagrangian derivative following the
particle. As before, the fluid is assumed incompressible (so that continuity requires
∂u i /∂ xi = 0) and Newtonian. The velocity that the fluid would have at the xi origin
in the absence of the particle is then W (t) in the x1 direction. It is also convenient to
define the quantities r , θ, u r , u θ as shown in Figure 2.1 and the Stokes streamfunction
as in Eq. (2.6). In some cases we shall also be able to consider the unsteady effects
due to growth of the bubble so the radius is denoted by R(t).
First consider inviscid potential flow for which Eq. (2.30) may be integrated to
obtain the Bernoulli equation as follows:
∂φ P 1 2
+ + u + u r2 = constant (2.32)
∂t ρC 2 θ
where φ is a velocity potential (u i = ∂φ/∂ xi ) and ψ must satisfy the following equa-
tion:

∂2 sin θ ∂ 1 ∂
Lψ = 0 where L ≡ 2 + 2 . (2.33)
∂r r ∂θ sin θ ∂θ
This is of course the same equation as in steady flow and has harmonic solutions, only
five of which are necessary for present purposes:

Wr2 D 2Ar 3 B
ψ = sin θ −
2
+ + cos θ sin θ
2
− 2 + E cos θ (2.34)
2 r 3 r

D 1 B E
φ = cos θ −W r + 2 + cos θ − 2
Ar + 3 +
2
(2.35)
r 3 r r

2D 1 3B E
u r = cos θ −W − 3 + cos θ − 2
2Ar − 4 − 2 (2.36)
r 3 r r

D B
u θ = − sin θ −W + 3 − 2 cos θ sin θ Ar + 4 . (2.37)
r r

The first part, which involves W and D, is identical to that for steady translation. The
second, involving A and B, will provide the fluid velocity gradient in the x 1 direction,
and the third, involving E, permits a time-dependent particle (bubble) radius. The W
and A terms represent the fluid flow in the absence of the particle, and the D, B, and
E terms allow the following boundary condition:
dR
(u r )r =R = (2.38)
dt
to be satisfied provided the following hold:
W R3 2A R 5 dR
D=− , B= , E = −R 2 . (2.39)
2 3 dt
In the absence of the particle the velocity of the fluid at the origin, r = 0, is simply −W
in the x1 direction and the gradient of the velocity ∂u 1 /∂ x1 = 4A/3. Hence A is
determined from the fluid velocity gradient in the original frame as the following:
3 ∂U
A= . (2.40)
4 ∂ x1∗
Now the force, F1 , on the bubble in the x1 direction is given by the following:
π
F1 = −2π R 2
p sin θ cos θdθ, (2.41)
0
which upon using Eqs. (2.31), (2.32), and (2.35) to (2.37) can be integrated to yield
the following:
F1 D 4 2 dV
=− (W R) − RW A + R . (2.42)
2π R ρC
2 Dt 3 3 dt
Reverting to the original coordinate system and using v as the sphere volume for
convenience (v = 4π R 3 /3), one obtains the following:
1 dV 3 DU 1 dv
F1 = − ρC v ∗ + ρC v ∗ + ρC (U − V ) ∗ , (2.43)
2 dt 2 Dt 2 dt
where the two Lagrangian time derivatives are defined by the following:
D ∂ ∂
∗
≡ ∗ +U ∗ (2.44)
Dt ∂t ∂ x1
d ∂ ∂
∗
≡ ∗ + V ∗. (2.45)
dt ∂t ∂ x1
Equation (2.43) is an important result, and care must be taken not to confuse the
different time derivatives contained in it. Note that in the absence of bubble growth,
viscous drag, and body forces, the equation of motion that results from setting
F1 = m p d V /dt ∗ is as follows:

2m p d V DU
1+ = 3 , (2.46)
ρC v dt ∗ Dt ∗

where m p is the mass of the particle. Thus for a massless bubble the acceleration of
the bubble is three times the fluid acceleration.
In a more comprehensive study of unsteady potential flows Symington (1978) has
shown that the result for more general (i.e., noncolinear) accelerations of the fluid and
particle is merely the vector equivalent of Eq. (2.43):
1 d Vi 3 DUi 1 dv
Fi = − ρC v ∗ + ρC v ∗ + ρC (Ui − Vi ) ∗ , (2.47)
2 dt 2 Dt 2 dt
where
d ∂ ∂ D ∂ ∂
∗
= ∗ + Vj ∗ ; ∗
= ∗ + Uj ∗ . (2.48)
dt ∂t ∂x j Dt ∂t ∂x j
The first term in Eq. (2.47) represents the conventional added mass effect due to the
particle acceleration. The factor 3/2 in the second term due to the fluid acceleration may
initially seem surprising. However, it is made up of two components: (1) 12 ρC d Vi /dt ∗ ,
which is the added mass effect of the fluid acceleration, and (2) ρC v DUi /Dt ∗ , which
is a buoyancy-like force due to the pressure gradient associated with the fluid accel-
eration. The last term in Eq. (2.47) is caused by particle (bubble) volumetric growth,
dv/dt ∗ , and is similar in form to the force on a source in a uniform stream.
Now it is necessary to ask how this force given by Eq. (2.47) should be used in the
practical construction of an equation of motion for a particle. Frequently, a viscous
drag force FiD , is quite arbitrarily added to Fi to obtain some total effective force on
the particle. Drag forces, FiD , with the following conventional forms:
CD
FiD = ρC |Ui − Vi | (Ui − Vi ) π R 2 (Re 1) (2.49)
2
FiD = 6πµC (Ui − Vi ) R (Re 1) (2.50)
have both been employed in the literature. It is, however, important to recognize that
there is no fundamental analytical justification for such superposition of these forces.
At high Reynolds numbers, we noted in the last section that experimentally observed
added masses are indeed quite close to those predicted by potential flow within certain
parametric regimes, and hence the superposition has some experimental justification.
At low Reynolds numbers, it is improper to use the results of the potential flow
analysis. The appropriate analysis under these circumstances is examined in the next
section.
2.3.4 Unsteady Stokes Flow
To elucidate some of the issues raised in the last section, it is instructive to examine
solutions for the unsteady flow past a sphere in low Reynolds number Stokes flow.
In the asymptotic case of zero Reynolds number, the solution of Section 2.2.2 is
unchanged by unsteadiness, and hence the solution at any instant in time is identical
to the steady-flow solution for the same particle velocity. In other words, because the

fluid has no inertia, it is always in static equilibrium. Thus the instantaneous force is
identical to that for the steady flow with the same Vi (t).
The next step is therefore to investigate the effects of small but nonzero inertial
contributions. The Oseen solution provides some indication of the effect of the con-
vective inertial terms, u j ∂u i /∂ x j , in steady flow. Here we investigate the effects of the
unsteady inertial term, ∂u i /∂t. Ideally it would be best to include both the ∂u i /∂t term
and the Oseen approximation to the convective term, U ∂u i /∂ x. However, the resulting
unsteady Oseen flow is sufficiently difficult that only small-time expansions for the
impulsively started motions of droplets and bubbles exist in the literature (Pearcey and
Hill 1956).
Consider, therefore, the unsteady Stokes equations in the absence of the convective
inertial terms:
∂u i ∂P ∂ 2ui
ρC =− + µC . (2.51)
∂t ∂ xi ∂x j∂x j
Because both the equations and the boundary conditions used below are linear in u i ,
we need only consider colinear particle and fluid velocities in one direction, say x1 .
The solution to the general case of noncolinear particle and fluid velocities and ac-
celerations may then be obtained by superposition. As in Section 2.3.3 the colinear
problem is solved by first transforming to an accelerating coordinate frame, xi , fixed
in the center of the particle so that P = p + ρC x1 d V /dt. Elimination of P by taking
the curl of Eq. (2.51) leads to the following:

1 ∂
L− Lψ = 0, (2.52)
νC ∂t
where L is the same operator as defined in Eq. (2.33). Guided by both the steady
Stokes flow and the unsteady potential flow solution, one can anticipate a solution of
the following form:
ψ = sin2 θ f (r, t) + cos θ sin2 θ g(r, t) + cos θ h(t), (2.53)
plus other spherical harmonic functions. The first term has the form of the steady
Stokes flow solution; the last term would be required if the particle were a growing
spherical bubble. After substituting Eq. (2.53) into Eq. (2.52), the equations for f , g,
h are as follows:

1 ∂ ∂2 2
L1 − L 1 f = 0 where L 1 ≡ 2 − 2 (2.54)
νC ∂t ∂r r

1 ∂ ∂2 6
L2 − L 2 g = 0 where L 2 ≡ 2 − 2 (2.55)
νC ∂t ∂r r

1 ∂ ∂2
L0 − L 0 h = 0 where L 0 ≡ 2 . (2.56)
νC ∂t ∂r

Moreover, the form of the expression for the force, F1 , on the spherical particle (or
bubble) obtained by evaluating the stresses on the surface and integrating is as follows:

F1 dV 1 ∂2 f νC 2 ∂ f 2 ∂2 f ∂3 f
= + + + − 3 . (2.57)
4
ρ π R3
3 C
dt r ∂r ∂t r r 2 ∂r r ∂r 2 ∂r r =R
It transpires that this is independent of g or h. Hence only the solution to Eq. (2.54) for
f (r, t) need be sought to find the force on a spherical particle, and the other spherical
harmonics that might have been included in Eq. (2.53) are now seen to be unnecessary.
Fourier or Laplace transform methods may be used to solve Eq. (2.54) for f (r, t),
and we choose Laplace transforms. The Laplace transforms for the relative velocity,
W (t), and the function f (r, t) are denoted by Ŵ (s) and fˆ(r, s) as follows:
∞ ∞
−st
Ŵ (s) = e W (t)dt; fˆ(r, s) = e−st f (r, t)dt. (2.58)
0 0
Then Eq. (2.54) becomes the following:
L 1 − ξ 2 L 1 fˆ = 0, (2.59)
1
where ξ = (s/νC ) 2 , and the solution after application of the condition that û 1 (s, t) far
from the particle be equal to Ŵ (s) is as follows:
2
Ŵ r A(s) 1
fˆ = − + + B(s) + ξ e−ξr , (2.60)
2 r r
where A and B are functions of s whose determination requires application of the
boundary conditions on r = R. In terms of A and B the Laplace transform of the
force F̂ 1 (s) is as follows:

F̂ 1 dV s ∂ fˆ νC 4Ŵ 8A
= + + − + 4 + C Be−ξr , (2.61)
4
ρ
3 C
π R 3 dt r ∂r R r r
r =R
where
3ξ 3 3ξ 2 8ξ 8
C = ξ4 + + 2 + 3 + 4. (2.62)
r r r r
The classical solution (see Landau and Lifshitz 1959) is for a solid sphere (i.e., constant
R) using the no-slip (Stokes) boundary condition for which

∂ f
f (R, t) = =0 (2.63)
∂r r =R
and hence
Wˆ R 3 3Wˆ RνC 3Wˆ RνC ξ R
A=+ + {1 + ξ R}; B=− e (2.64)
2 2s 2s

so that
1
Fˆ1 dV 3 ˆ 9νC Wˆ 9νC2 1 ˆ
= − s W − − s 2 W. (2.65)
4
3
ρC π R 3 dt 2 2R 2 2R
For a motion starting at rest at t = 0 the inverse Laplace transform of this yields the
following:
t
F1 dV 3 dW 9νC 9 νC 12 dW (t̃) d t̃
= − − W − , (2.66)
4
ρ π R3
3 C
dt 2 dt 2R 2 2R π d t̃ (t − t̃) 12
0
where t̃ is a dummy time variable. This result must then be written in the original
coordinate framework with W = V − U and can be generalized to the noncolinear
case by superposition so that
1 d Vi 3 dUi 9vµC
Fi = − vρC ∗ + vρC ∗ + (Ui − Vi )
2 dt 2 dt 2R 2
t∗
9vρC νC 12 d(Ui − Vi ) d t̃
+ 1 , (2.67)
2R π d t̃ ∗
(t − t̃) 2
0
∗
where d/dt is the Lagrangian time derivative following the particle. This is then
the general force on the particle or bubble in unsteady Stokes flow when the Stokes
boundary conditions are applied.
Compare this result with that obtained from the potential flow analysis, Eq. (2.47)
with v taken as constant. It is striking to observe that the coefficients of the added
mass terms involving d Vi /dt ∗ and dUi /dt ∗ are identical to those of the potential
flow solution. On superficial examination it might be noted that dUi /dt ∗ appears in
Eq. (2.67), whereas DUi /Dt ∗ appears in Eq. (2.47); the difference is, however, of order
W j ∂Ui /d x j and terms of this order have already been dropped from the equation of
motion on the basis that they were negligible compared with the temporal derivatives
like ∂ Wi /∂t. Hence it is inconsistent with the initial assumption to distinguish between
d/dt ∗ and D/Dt ∗ in the present unsteady Stokes flow solution.
The term 9νC W/2R 2 in Eq. (2.67) is, of course, the steady Stokes drag. The new
phenomenon introduced by this analysis is contained in the last term of Eq. (2.67).
This is a fading memory term that is often named the Basset term after one of its
identifiers (Basset 1888). It results from the fact that additional vorticity created at the
solid particle surface due to relative acceleration diffuses into the flow and creates a
1
temporary perturbation in the flow field. Like all diffusive effects it produces an ω 2
term in the equation for oscillatory motion.
Before we conclude this section, comment should be included on three other ana-
lytical results. Morrison and Stewart (1976) considered the case of a spherical bubble
for which the Hadamard–Rybczynski boundary conditions rather than the Stokes con-
ditions are applied. Then, instead of the conditions of Eq. (2.63), the conditions for

zero normal velocity and zero shear stress on the surface require that
2
∂ f 2 ∂f
f (R, t) = − =0 (2.68)
∂r 2 r ∂r r =R
and hence in this case (see Morrison and Stewart 1976)
ˆ 3 3ŴR (1 + ξ R)
WR 3Ŵ Re+ξ R
A(s) = + + 2 ; B(s) = − (2.69)
2 ξ (3 + ξ R) ξ 2 (3 + ξ R)
so that
F̂ 1 dV 9Ŵ νC 3 6νC Ŵ
= − − Ŵ s + 1 . (2.70)
4
3
πρC R 3 dt R 2 2 1
R 2 1 + s 2 R/3νC2
The inverse Laplace transform of this for motion starting at rest at t = 0 is as follows:
t
F1 dV 3 dW 3νC W 6νC dW (t̃)
= − − − 2
4
ρ πR
3 C
3 dt 2 dt R 2 R d t̃
0 1
9νC t − t̃ 9νC (t − t̃) 2
× exp erfc d t̃. (2.71)
R2 R2
Comparing this with the solution for the Stokes conditions, we note that the first two
terms are unchanged and the third term is the expected Hadamard–Rybczynski steady
drag term [see Eq. (2.16)]. The last term is significantly different from the Basset term
in Eq. (2.67) but still represents a fading memory.
More recently, Magnaudet and Legendre (1998) have extended these results further
by obtaining an expression for the force on a particle (bubble) whose radius is changing
with time.
Another interesting case is that for unsteady Oseen flow, which essentially consists
of attempting to solve the Navier–Stokes equations with the convective inertial terms
approximated by U j ∂u i /∂ x j . Pearcey and Hill (1956) have examined the small-time
behavior of droplets and bubbles started from rest when this term is included in the
equations.
2.4 Particle Equation of Motion
2.4.1 Equations of Motion

In a multiphase flow with a very dilute discrete phase the fluid forces discussed in
Sections 2.1 to 2.3.4 determine the motion of the particles that constitute that discrete
phase. In this section we discuss the implications of some of the fluid force terms. The
equation that determines the particle velocity, Vi , is generated by equating the total
force, FiT , on the particle to m p d Vi /dt ∗ . Consider the motion of a spherical particle

2.4 Particle Equation of Motion 49
(or bubble) of mass m p and volume v (radius R) in a uniformly accelerating fluid. The
simplest example of this is the vertical motion of a particle under gravity, g, in a pool of
otherwise quiescent fluid. Thus the results are written in terms of the buoyancy force.
However, the same results apply to motion generated by any uniform acceleration
of the fluid, and hence g can be interpreted as a general uniform fluid acceleration
(dU/dt). This will also allow some tentative conclusions to be drawn concerning
the relative motion of a particle in the nonuniformly accelerating fluid situations that
can occur in general multiphase flow. For the motion of a sphere at small relative
Reynolds number, Re 1 (where Re = 2W R/νC and W is the typical magnitude of
the relative velocity), only the forces due to buoyancy and the weight of the particle
need be added to Fi as given by Eqs. (2.67) or (2.71) to obtain FiT . This addition is
simply given by (ρC v − m p )gi where g is a vector in the vertically upward direction
with magnitude equal to the acceleration due to gravity. Conversely, at high relative
Reynolds numbers, Re 1, one must resort to a more heuristic approach in which
the fluid forces given by Eq. (2.47) are supplemented by drag (and lift) forces given by
1
ρ ACi j |W j |W j as in Eq. (2.27). In either case it is useful to nondimensionalize the
2 C
resulting equation of motion so that the pertinent nondimensional parameters can be
identified.
Examine first the case in which the relative velocity, W (defined as positive in
the direction of the acceleration, g, and therefore positive in the vertically upward
direction of the rising bubble or sedimenting particle), is sufficiently small so that the
relative Reynolds number is much less than unity. Then, using the Stokes boundary
conditions, the equation governing W may be obtained from Eq. (2.66) as follows:
12 t∗
dw 9 dw d t˜
w+ + =1 (2.72)
dt∗ π 1 + 2m p /ρC v d t˜ (t∗ − t̃) 12
0
where the dimensionless time, t∗ = t/tu and the relaxation time, tu , is given by the
following:
tu = R 2 1 + 2mp /ρC v /9νC (2.73)
and
w = W/W∞ ,
where W∞ is the steady terminal velocity given by the following:
W∞ = 2R 2 g 1 − m p /ρC v /9νC . (2.74)
In the absence of the Basset term the solution of Eq. (2.72) is simply the following:
w = 1 − e−t/tu (2.75)
and therefore the typical response time is given by the relaxation time, tu (see, for
example, Rudinger 1969 and Section 1.2.7). In the general case that includes the Basset

Figure 2.5. The velocity, W, of a particle released from rest at t∗ = 0 in a quiescent fluid and its approach
to terminal velocity, W∞ . Horizontal axis is a dimensionless time defined in text. Solid lines represent
the low-Reynolds-number solutions for various particle mass/displaced mass ratios, m p /ρC v, and the
Stokes boundary condition. The dashed line is for the Hadamard–Rybczynski boundary condition and
m p /ρC v = 0. The dash-dot line is the high-Reynolds-number result; note that t∗ is nondimensionalized
differently in that case.
term the dimensionless solution, w(t∗ ), of Eq. (2.72) depends only on the parameter
m p /ρC v (particle mass/displaced fluid mass) appearing in the Basset term. Indeed,
the dimensionless Eq. (2.72) clearly illustrates the fact that the Basset term is much
less important for solid particles in a gas where m p /ρC v 1 than it is for bubbles
in a liquid where m p /ρC v 1. Note also that for initial conditions of zero relative
velocity (w(0) = 0) the small-time solution of equation 2.72 takes the following form:
2 3
ω = t∗ − 1 12 t∗ + · · · .
2
(2.76)
π 2 1 + 2m p /ρC v
Hence the initial acceleration at t = 0 is given dimensionally by the following:
2g 1 − m p /ρC v / 1 + 2m p /ρC v
or 2g in the case of a massless bubble and −g in the case of a heavy solid particle
in a gas where m p ρC v. Note also that the effect of the Basset term is to reduce
the acceleration of the relative motion, thus increasing the time required to achieve
terminal velocity.
Numerical solutions of the form of w(t∗ ) for various m p /ρC v are shown in
Figure 2.5 where the delay caused by the Basset term can be clearly seen. In fact
in the later stages of approach to the terminal velocity the Basset term dominates
over the added mass term, (dw/dt∗ ). The integral in the Basset term becomes

1/2
approximately 2t∗ dw/dt∗ so that the final approach to w = 1 can be approximated
by the following:
12
1 9
w = 1 − C exp −t∗2 , (2.77)
π{1 + 2m p /ρC v}
where C is a constant. As can be seen in Figure 2.5, the result is a much slower
approach to W∞ for small m p /ρC v than for larger values of this quantity.
The case of a bubble with Hadamard–Rybczynski boundary conditions is very
similar except that
W∞ = R 2 g 1 − m p /ρC v /3νC (2.78)
and the equation w(t∗ ) is as follows:

t∗
3 dw dw
w+ +2 t∗ − t˜ d t˜ = 1, (2.79)
2 dt∗ d t̃
0
where the function, (ξ ), is given by the following:

12
2m p 2m p
(ξ ) = exp 1 + ξ erfc 1+ ξ . (2.80)
ρC v ρC v
For the purposes of comparison the form of w(t∗ ) for the Hadamard–Rybczynski
boundary condition with m p /ρC v = 0 is also shown in Figure 2.5. Though the altered
Basset term leads to a more rapid approach to terminal velocity than occurs for the
Stokes boundary condition, the difference is not qualitatively significant.
If the terminal Reynolds number is much greater than unity then, in the ab-
sence of particle growth, Eq. (2.47) heuristically supplemented with a drag force
of the form of Eq. (2.49) leads to the following equation of motion for unidirectional
motion:
dw
w2 + = 1, (2.81)
dt∗
where w = W/W∞ , t∗ = t/tu , and the relaxation time, tu , is now given by the
following:
1
tu = 1 + 2m p /ρC v 2R/3CD g(1 − m p /vρC ) 2
(2.82)
and
1
W∞ = 8Rg 1 − m p /ρC v /3CD 2 . (2.83)
The solution to Eq. (2.81) for w(0) = 0,
w = tanh t∗ , (2.84)

is also shown in Figure 2.5 though, of course, t∗ has a different definition in this
case.
The relaxation times given by the Eqs. (2.73) and (2.82) are particularly valuable
in assessing relative motion in disperse multiphase flows.When this time is short
compared with the typical time associated with the fluid motion, the particle will
essentially follow the fluid motion and the techniques of homogeneous flow (see
Chapter 9) are applicable. Otherwise the flow is more complex and special effort is
needed to evaluate the relative motion and its consequences.
For the purposes of reference in Section 3.2 note that, if we define a Reynolds
number, Re, and a Froude number, Fr, by the following:
2W∞ R W∞
Re = ; Fr = 12 , (2.85)
νC 2Rg(1 − m p /ρC v)
then the expressions for the terminal velocities, W∞ , given by Eqs. (2.74), (2.78), and
(2.83), can be written as follows:
1 1 1
Fr = (Re/18) 2 , Fr = (Re/12) 2 , and Fr = (4/3CD ) 2 (2.86)
respectively. Indeed, dimensional analysis of the governing Navier–Stokes equations
requires that the general expression for the terminal velocity can be written as follows:
F(Re, Fr) = 0 (2.87)
or, alternatively, if CD is defined as 4/3Fr 2 , then it could be written as follows:
F ∗ (Re, CD ) = 0. (2.88)
2.4.2 Magnitude of Relative Motion

Qualitative estimates of the magnitude of the relative motion in multiphase flows can
be made from the analyses of the last section. Consider a general steady fluid flow
characterized by a velocity, U , and a typical dimension, ; it may, for example, be useful
to visualize the flow in a converging nozzle of length, , and mean axial velocity, U .
A particle in this flow will experience a typical fluid acceleration (or effective g) of
U 2 / for a typical time given by /U and hence will develop a velocity, W , relative
to the fluid. In many practical flows it is necessary to determine the maximum value
of W (denoted by Wm ) that could develop under these circumstances. To do so, one
must first consider whether the available time, /U , is large or small compared with
the typical time, tu , required for the particle to reach its terminal velocity as given by
Eq. (2.73) or (2.82). If tu /U then Wm is given by Eqs. (2.74), (2.78), or (2.83)
for W∞ and qualitative estimates for Wm /U would be as follows:
1 12
mp UR R mp 2 1 R
1− and 1− , (2.89)
ρC v νC ρC v 1
CD2

When W R/νC 1 and W R/νC 1 respectively. We refer to this as the quasistatic

regime. Conversely, if tu /U, Wm can be estimated as W∞ /U tu so that Wm /U is
of the order of
2 1 − mp /ρC v
(2.90)
1 + 2mp /ρC v
for all W R/νC . This is termed the transient regime.
In practice, W R/νC will not be known in advance. The most meaningful quantities
that can be evaluated prior to any analysis are a Reynolds number, U R/νC , based on
flow velocity and particle size, a size parameter

R m p
X = 1 − (2.91)
ρC v
and the parameter

m p 2m p

Y = 1 − 1+ . (2.92)
ρC v ρC v
The resulting regimes of relative motion are displayed graphically in Figure 2.6. The
transient regime in the upper right-hand sector of the graph is characterized by large
relative motion, as suggested by Eq. (2.90). The quasistatic regimes for W R/νC 1
and W R/νC 1 are in the lower right- and left-hand sectors respectively. The shaded
boundaries between these regimes are, of course, approximate and are functions of the
parameter Y , that must have a value in the range 0 < Y < 1. As one proceeds deeper
into either of the quasistatic regimes, the magnitude of the relative velocity, Wm /U ,
becomes smaller and smaller. Thus, homogeneous flows (see Chapter 9) in which
the relative motion is neglected require that either X Y 2 or X Y/(U R/νC ).
Conversely, if either of these conditions is violated, relative motion must be included
in the analysis.
2.4.3 Effect of Concentration on Particle Equation of Motion

When the concentration of the disperse phase in a multiphase flow is small (less than,
say, 0.01% by volume) the particles have little effect on the motion of the continuous
phase and analytical or computational methods are much simpler. Quite accurate
solutions are then obtained by solving a single-phase flow for the continuous phase
(perhaps with some slightly modified density) and inputting those fluid velocities into
equations of motion for the particles. This is known as one-way coupling.
As the concentration of the disperse phase is increased a whole spectrum of com-
plications can arise. These may effect both the continuous phase flow and the disperse
phase motions and flows with this two-way coupling pose many modeling challenges.
A few examples are appropriate. The particle motions may initiate or alter the tur-
bulence in the continuous phase flow; this particularly challenging issue is briefly

Figure 2.6. Schematic of the various regimes of relative motion between a particle and the surrounding
flow.
addressed in Section 1.3. Moreover, particles may begin to collide with one another,
altering their effective equation of motion and introducing random particle motions that
may need to be accounted for; Chapter 13 is devoted to flows dominated by such colli-
sions. These collisions and random motions may generate additional turbulent motions
in the continuous phase. Often the interactions of particles become important even if
they do not actually collide. Fortes et al. (1987) have shown that in flows with high
relative Reynolds numbers there are several important mechanisms of particle/particle
interactions that occur when a particle encounters the wake of another particle. The
following particle drafts the leading particle and impacts it when it catches up with it
and the pair then begin tumbling. In packed beds these interactions result in the
development of lateral bands of higher concentration separated by regions of low,
almost zero, volume fraction. How these complicated interactions could be incorpo-
rated into a two-fluid model (short of complete and direct numerical simulation) is
unclear.
At concentrations that are sufficiently small so that the complications of the pre-
ceding paragraph do not arise, there are still effects on the coefficients in the particle

equation of motion that may need to be accounted for. For example, the drag on a
particle or the added mass of a particle may be altered by the presence of neighboring
particles. These issues are somewhat simpler to deal with than those of the preceding
paragraph and we cover them in this chapter. The effect on the added mass was ad-
dressed earlier in Section 2.3.2. In the next section we address the issue of the effect
of concentration on the particle drag.
2.4.4 Effect of Concentration on Particle Drag
Section 2.2 reviewed the dependence of the drag coefficient on the Reynolds number
for a single particle in a fluid and the effect on the sedimentation of that single particle in
an otherwise quiescent fluid was examined as a particular example in Section 2.4. Such
results would be directly applicable to the evaluation of the relative velocity between
the disperse phase (the particles) and the continuous phase in a very dilute multiphase
flow. However, at higher concentrations, the interactions between the flow fields around
individual particles alter the force experienced by those particles and therefore change
the velocity of sedimentation. Furthermore, the volumetric flux of the disperse phase is
no longer negligible because of the finite concentration and, depending on the boundary
conditions in the particular problem, this may cause a nonnegligible volumetric flux of
the continuous phase. For example, particles sedimenting in a containing vessel with
a downward particle volume flux, − jS (upward is deemed the positive direction), at a
concentration, α, will have a mean velocity as follows:
−u S = − jS /α (2.93)
and will cause an equal and opposite upward flux of the suspending liquid, jL = − jS ,
so that the mean velocity of the liquid is as follows:
u L = jL /(1 − α) = − jS /(1 − α). (2.94)
Hence the relative velocity is as follows:
u SL = u S − u L = jS /α(1 − α) = u S /(1 − α). (2.95)
Thus care must be taken to define the terminal velocity and here we focus on the more
fundamental quantity, namely the relative velocity, u SL , rather than on quantities such
as the sedimentation velocity, u S , that are dependent on the boundary conditions.
Barnea and Mizrahi (1973) have reviewed the experimental, theoretical and em-
pirical data on the sedimentation of particles in otherwise quiescent fluids at various
concentrations, α. The experimental data of Mertes and Rhodes (1955) on the ratio
of the relative velocity, u SL , to the sedimentation velocity for a single particle, (u SL )0
(equal to the value of u SL as α → 0), are presented in Figure 2.7. As one might antici-
pate, the relative motion is hindered by the increasing concentration. It can also be seen
that u SL /(u SL )0 is not only a function of α but varies systematically with the Reynolds

Figure 2.7. Relative velocity of sedimenting particles, u SL [normalized by the velocity as α → 0, (u SL )0 ]

as a function of the volume fraction, α. Experimental data from Mertes and Rhodes (1955) are shown for
various Reynolds numbers, Re, as follows: Re = 0.003 (+), 0.019 (×), 0.155 (), 0.98 (), 1.45 (),
4.8 (∗), 16 (), 641, (), 020 (), and 2180 (). Also shown are the analytical results of Brinkman
[Eq. (2.97)] and Zick and Homsy and the empirical results of Wallis [Eq. (2.100)] and Barnea and Mizrahi
[Eq. (2.98)].
number, 2R(u SL )0 /νL , where νL is the kinematic viscosity of the suspending medium.
Specifically, u SL /(u SL )0 increases significantly with Re so that the rate of decrease
of u SL /(u SL )0 with increasing α is lessened as the Reynolds number increases. One
might intuitively expect this decrease in the interactions between the particles because
the far field effects of the flow around a single particle decline as the Reynolds number
increases.
We also note that complementary to the data of Figure 2.7 is extensive data on the
flow through packed beds of particles. The classical analyses of that data by Kozeny
(1927) and, independently, by Carman (1937) led to the widely used expression for
the pressure drop in the low Reynolds number flow of a fluid of viscosity, µC , and
superficial velocity, jCD , through a packed bed of spheres of diameter, D, and solids
volume fraction, α, namely
dp 180α 3 µC jCD
= , (2.96)
ds (1 − α)3 D 2
where the 180 and the powers on the functions of α were empirically determined. This
expression, known as the Carman–Kozeny equation, is used shortly.
Several curves that are representative of the analytical and empirical results are
also shown in Figure 2.7 (and in Figure 2.8). One of the first approximate, analytical

Figure 2.8. The drift flux, jSL (normalized by the velocity (u SL )0 ) corresponding to the relative velocities
of Figure 2.7 (see that caption for codes).
models to include the interactions between particles was that of Brinkman (1947)
for spherical particles at asymptotically small Reynolds numbers who obtained the
following:
u SL (2 − 3α)2
= 1 (2.97)
(u SL )0 4 + 3α + 3(8α − 3α 2 ) 2
and this result is included in Figures 2.7 and 2.8. Other researchers (see, for example,
Tam 1969 and Brady and Bossis 1988) have studied this low Reynolds number limit
quite closely. Exact solutions for the sedimentation velocity of a various regular arrays
of spheres at asymptotically low Reynolds number were obtained by Zick and Homsy
(1982) and the particular result for a simple cubic array is included in Figure 2.7.
Clearly, these results deviate significantly from the experimental data and it is currently
thought that the sedimentation process cannot be modeled by a regular array because
the fluid mechanical effects are dominated by the events that occur when particles
happen to come close to one another.
Switching attention to particle Reynolds numbers greater than unity, it was men-
tioned earlier that the work of Fortes et al. (1987) and others has illustrated that the
interactions between particles become very complex because they result, primarily,
from the interactions of particles with the wakes of the particles ahead of them. Fortes
et al. (1987) have shown this results in a variety of behaviors they term drafting,
kissing, and tumbling that can be recognized in fluidized beds. As yet, these behaviors
have not been amenable to theoretical analyses.

Figure 2.9. Data indicating the variation

in the bubble relative velocity, u GL , with
the void fraction, α, and the overall flow
rate (as represented by jL ) in a vertical,
10.2-cm diameter tube. The dashed line
is the correlation of Wallis [Eq. (2.100)].
Adapted from Bernier (1982).
The literature contains numerous empirical correlations but three will suffice for
present purposes. At small Reynolds numbers, Barnea and Mizrahi (1973) show that
the experimental data closely follow an expression of the following form:
u SL (1 − α)
≈ 1 (2.98)
(u SL )0 (1 + α 3 )e5α/3(1−α)
By way of comparison the Carman–Kozeny Eq. (2.96) implies that a sedimenting
packed bed would have a terminal velocity given by the following:
u SL 1 (1 − α)2
= , (2.99)
(u SL )0 80 α 2
which has magnitudes comparable to the Eq. (2.98) at the volume fractions of packed
beds.
At large rather than small Reynolds number, the ratio u SL /(u SL )0 seems to be better
approximated by the following empirical relation:
u SL
≈ (1 − α)b−1 , (2.100)
(u SL )0
where Wallis (1969) suggests a value of b = 3. Both of these empirical formulae are
included in Figure 2.7.
In later chapters discussing sedimentation phenomena, we use the drift flux, jSL ,
more frequently than the relative velocity, u SL . Recalling that jSL = α(1 − α)u SL , the
data from Figure 2.7 are replotted in Figure 2.8 to display jSL /(u SL )0 .

It is appropriate to end by expressing some reservations regarding the generality

of the experimental data presented in Figures 2.7 and 2.8. At the higher concentra-
tions, vertical flows of this type often develop instabilities that produce large-scale
mixing motions whose scale is of the same order as the horizontal extent of the
flow, usually the pipe or container diameter. In turn, these motions can have a sub-
stantial effect on the mean sedimentation velocity. Consequently, one might expect
a pipe size effect that would manifest itself nondimensionally as a dependence on a
parameter such as the ratio of the particle to pipe diameter, 2R/d, or, perhaps, in
1
a Froude number such as (u SL )0 /(gd) 2 . Another source of discrepancy could be a
dependence on the overall flow rate. Almost all of the data, including that of Mertes
and Rhodes (1955), have been obtained from relatively quiescent sedimentation or
fluidized bed experiments in which the overall flow rate is small and, therefore, the
level of turbulence is limited to that produced by the relative motion between the
particles and the suspending fluid. However, when the overall flow rate is increased
so that even a single-phase flow of the suspending fluid would be turbulent, the mean
sedimentation velocities may be significantly altered by the enhancement of the mix-
ing and turbulent motions. Figure 2.9 presents data from some experiments by Bernier
(1982) in which the relative velocity of bubbles of air in a vertical water flow were
measured for various total volumetric fluxes, j. Small j values cause little deviation
from the behavior at j = 0 and are consistent with the results of Figure 2.7. However,
at larger j values for which a single-phase flow would be turbulent, the decrease in
u GL with increasing α almost completely disappears. Bernier surmised that this dis-
appearance of the interaction effect is due to the increase in the turbulence level in the
flow that essentially overwhelms any particle/particle or bubble/bubble interaction.

06.0 PP 30 59 Single-Particle Motion

Uploaded by

Copyright:

Available Formats

06.0 PP 30 59 Single-Particle Motion

Uploaded by

Document Information

Original Title

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

06.0 PP 30 59 Single-Particle Motion

Uploaded by

Copyright:

Available Formats

2

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

Figure 2.1. Notation for a spherical particle.

2.2 Flows Around a Sphere

2.2.1 At High Reynolds Number

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

The Stokes streamfunction, ψ, is defined to satisfy continuity automatically:

2.2.2 At Low Reynolds Number

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

Figure 2.2. Smoke visualization of the nom-

Figure 2.3. Drag coefficient on a sphere

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

which yields a drag force as follows:

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

the wake begins to oscillate (Taneda 1956). However, it continues to be attached to

Str = 2 f R/W (2.22)

the vortex shedding frequencies, f , that Moller (1938) observed correspond to a

2.2.3 Molecular Effects

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

2.3 Unsteady Effects

2.3.1 Unsteady Particle Motions

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

Particle Matrix element Value

Ellipsoids (T) Mii = K ii 34 ρC πab2 a/b K 11 K 22 (K 33 )

Sphere near free Mii = K ii 34 ρC π R 3 H/R K 11

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

2.3.2 Effect of Concentration on Added Mass

2.3.3 Unsteady Potential Flow

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

2.3.4 Unsteady Stokes Flow

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

ψ = sin2 θ f (r, t) + cos θ sin2 θ g(r, t) + cos θ h(t), (2.53)

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

Then Eq. (2.54) becomes the following:

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

2.4 Particle Equation of Motion

2.4.1 Equations of Motion

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

W∞ = R 2 g 1 − m p /ρC v /3νC (2.78)

and the equation w(t∗ ) is as follows:

where the function, (ξ ), is given by the following:

The solution to Eq. (2.81) for w(0) = 0,

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

2.4.2 Magnitude of Relative Motion

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

When W R/νC 1 and W R/νC 1 respectively. We refer to this as the quasistatic

2.4.3 Effect of Concentration on Particle Equation of Motion

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

https://doi.org/10.1017/CBO9780511807169.003 Published online by Cambridge University Press

2.4.4 Effect of Concentration on Particle Drag

u L = jL /(1 − α) = − jS /(1 − α). (2.94)