Group Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Title: 6.2.1 - Influence of the history term in a Lagrangian method for oil-water separation
ALL VOLUMES CITATION THUMBNAILS PAGE IMAGE ZOOMABLE
Full Citation
STANDARD VIEW MARC VIEW
Permanent Link: http://ufdc.ufl.edu/UF00102023/00147
 Material Information
Title: 6.2.1 - Influence of the history term in a Lagrangian method for oil-water separation Particle Bubble and Drop Dynamics
Series Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Physical Description: Conference Papers
Creator: van Eijkeren, D.F.
Hoeijmakers, H.W.M.
Publisher: International Conference on Multiphase Flow (ICMF)
Publication Date: June 4, 2010
 Subjects
Subject: Lagrangian simulations
history term
swirling flow
oil-water separation
 Notes
Abstract: In the oil industry, both settling tanks and swirling flow separators are used to provide separation of phases. In the present paper the motion of oil droplets is investigated numerically for the case of swirling flow fields as well as for the case of settling in quiescent fluids. During phase separation the oil droplets are subjected to various forces, i.e. drag force, buoyancy force, pressure gradient force, virtual mass force, Saffman and Magnus lift force and the history or Basset force. In oil-water separation, the density of the continuous phase is of the same order of magnitude as that of the dispersed phase. This implies that, for a relatively dense continuous phase, the history term cannot be neglected. Several experimental and numerical data, e.g. Mordant and Pinton (2000); Mei et al. (1991); Loth and Dorgan (2009), show that using the classical Basset kernel for the history term overestimates this unsteady force significantly. In the present research, Lagrangian simulations have been performed using, for the history term, variants of the kernel proposed by Mei et al. (1991). Obtained results are compared with experimental and numerical data. The influence of the history force using various kernels has also been investigated for a generic swirling flow field, for the moment without accounting for effects of turbulence. The results show that inclusion of the history term leads to an increase of the time required for particles to reach the axis of rotation. This implies that in numerical methods for swirling flow separators the history term should be included. For the case of solid-body rotation, present results have been compared with experimental results by Candelier et al. (2005). This comparison shows that, in order to account for the forces adequately, additional effects have to be taken into account for solid-body rotation flow fields, as also suggested by Candelier (2008), who also accounts for non-Basset type of history effects.
General Note: The International Conference on Multiphase Flow (ICMF) first was held in Tsukuba, Japan in 1991 and the second ICMF took place in Kyoto, Japan in 1995. During this conference, it was decided to establish an International Governing Board which oversees the major aspects of the conference and makes decisions about future conference locations. Due to the great importance of the field, it was furthermore decided to hold the conference every three years successively in Asia including Australia, Europe including Africa, Russia and the Near East and America. Hence, ICMF 1998 was held in Lyon, France, ICMF 2001 in New Orleans, USA, ICMF 2004 in Yokohama, Japan, and ICMF 2007 in Leipzig, Germany. ICMF-2010 is devoted to all aspects of Multiphase Flow. Researchers from all over the world gathered in order to introduce their recent advances in the field and thereby promote the exchange of new ideas, results and techniques. The conference is a key event in Multiphase Flow and supports the advancement of science in this very important field. The major research topics relevant for the conference are as follows: Bio-Fluid Dynamics; Boiling; Bubbly Flows; Cavitation; Colloidal and Suspension Dynamics; Collision, Agglomeration and Breakup; Computational Techniques for Multiphase Flows; Droplet Flows; Environmental and Geophysical Flows; Experimental Methods for Multiphase Flows; Fluidized and Circulating Fluidized Beds; Fluid Structure Interactions; Granular Media; Industrial Applications; Instabilities; Interfacial Flows; Micro and Nano-Scale Multiphase Flows; Microgravity in Two-Phase Flow; Multiphase Flows with Heat and Mass Transfer; Non-Newtonian Multiphase Flows; Particle-Laden Flows; Particle, Bubble and Drop Dynamics; Reactive Multiphase Flows
 Record Information
Bibliographic ID: UF00102023
Volume ID: VID00147
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: 621-vanEijkeren-ICMF2010.pdf

Full Text



7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


Influence of the history term in a Lagrangian method for oil-water separation


D.F. van Eijkeren* and H.W.M. Hoeijmakers*

Faculty of Engineering Technology, University of Twente, PO Box 217, 7500 AE, Enschede, The Netherlands
d.f.vaneijkeren@utwente.nl and h.w.m.hoeijmakers@utwente.nl
Keywords: Lagrangian simulations, history term, swirling flow, oil-water separation




Abstract

In the oil industry, both settling tanks and swirling flow separators are used to provide separation of phases. In the
present paper the motion of oil droplets is investigated numerically for the case of swirling flow fields as well as for
the case of settling in quiescent fluids. During phase separation the oil droplets are subjected to various forces, i.e.
drag force, buoyancy force, pressure gradient force, virtual mass force, Saffman and Magnus lift force and the history
or Basset force. In oil-water separation, the density of the continuous phase is of the same order of magnitude as that
of the dispersed phase. This implies that, for a relatively dense continuous phase, the history term cannot be neglected.
Several experimental and numerical data, e.g. Mordant and Pinton (2000); Mei et al. (1991); Loth and Dorgan (2009),
show that using the classical Basset kernel for the history term overestimates this unsteady force significantly.
In the present research, Lagrangian simulations have been performed using, for the history term, variants of the kernel
proposed by Mei et al. (1991). Obtained results are compared with experimental and numerical data. The influence
of the history force using various kernels has also been investigated for a generic swirling flow field, for the moment
without accounting for effects of turbulence. The results show that inclusion of the history term leads to an increase
of the time required for particles to reach the axis of rotation. This implies that in numerical methods for swirling
flow separators the history term should be included. For the case of solid-body rotation, present results have been
compared with experimental results by Candelier et al. (2005). This comparison shows that, in order to account for
the forces adequately, additional effects have to be taken into account for solid-body rotation flow fields, as also
suggested by Candelier (2008), who also accounts for non-Basset type of history effects.


Nomenclature

Roman symbols
C coefficient for force abbreviated in subscript (-)
c empirical parameters in history term (-)
D particle or droplet diameter (m)
F force(kgms 2)
f force per unit mass (m s 2)
f drag factor (-)
fH history factor (-)
g gravitational constant (m s 2)
g gravitational acceleration vector (m s 2)
I moment of inertia (kg mn2)
K dimensionless kernel in integration (-)
m mass (kg)
R radius (m)
Re Reynolds number based on D(-)
t time (s)
T dimensionless time (-)
T Torque (kg n2s-2)
u flow velocity vector (m s 1)


v particle velocity vector (ms 1)
w v u (ms 1)
x coordinate in space (m)
Greek symbols
a dimensionless particle rotation (-)
3 dimensionless flow rotation (-)
F Circulation (m2 s 1)
p dynamic viscosity (kg m 1 s 1)
v kinematic viscosity (m2 s 1)
S arbitrary vector variable
p density (kg m 3)
p density ratio (-)
(a stress tensor (N m 2)
T characteristic and integration time (s)
4 arbitrary forcing vector
w rate of rotation (s1)
Subscripts
c continuous
d diffusive
D non-dimensional based on particle diameter











H history
v referring to viscous decay time
p particle, either solid or fluid
R non-dimensional based on initial radius
t terminal
V referring to relaxation time
w referring to rotation Reynolds number
0 at initial time to
Superscipts
T transposed
t at discrete time



Introduction

In oil industry, water is used to maintain pressure in oil
wells during the process of extraction of oil. Therefore, a
mixture of oil and water will be produced during the ex-
traction process. It is required to separate this oil-water
mixture, not only to obtain the desired oil, but also to
be able to re inject the water into the oil field, or to dis-
pense it in the environment. Separation of phases can
be divided in primary, or bulk, separation and secondary
separation. A mixture with a very small amount of oil,
as well as a mixture with a very small amount of wa-
ter are produced during the initial bulk separation. The
secondary separation of oil from produced water, is also
known as produced water treatment. Tiny oil droplets
that remain in the water after bulk separation will need
to be removed in order to obtain sufficiently clean water.
This requires more advanced separation methods.
Separation of oil and water is most often based on the
density difference between oil and water. A widely used
method based on the density difference is gravity set-
tling. A mixture is put in a settling tank and after some
time oil will float on top of the water and can be ex-
tracted from the original mixture. The method of sepa-
ration using gravity settling is most often implemented
as a batch process. Moreover, the settling times during
produced water treatment would be days or even longer.
Therefore, more advanced methods have been developed
to obtain separation of phases.
An important quantity to estimate the speed of set-
tling, and thus to identify improvements for separation,
is the terminal velocity vt of droplets. As a first approx-
imation the droplet will be approximated by a spherical
particle. Another approximation for an initial estimate
of vt is that the flow around the particle will be Stokes'
flow. Therefore the drag will be Stokes' drag. The ter-
minal velocity vt is then:


Dg (1-
18V,


g (1
P


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


with the particle relaxation time TV and density ratio p
defined as


=pD2
18V,


and p


respectively.
The most important parameter influencing the settling
time is the droplet diameter D. Coalescence of droplets
causes a huge improvement in settling time. However,
in dilute mixtures with very small droplets, coalescence
of droplets is not occurring very often. Another param-
eter in equation (1) that can be influenced is the gravita-
tional acceleration. Centrifugal forces create an artificial
gravitational field, stronger than the natural gravitational
field. This is used in swirling flow separators. Swirl is
generated by either an internal swirl element in the pipe
or a tangential inlet to the pipe. Due to the centrifugal
force the oil droplets will move towards the axis of rota-
tion. A major advantage with respect to settling tanks is
that a continuous, rather than a batch, separation process
is possible.
To control the efficiency of swirling flow separators,
it is not only necessary to be able to compute the flow
field, but also to be able to compute the droplet trajec-
tories. A small change in inward velocity of a particle
could imply, for example, that the separator should be
much longer than foreseen. Therefore, the equations of
motion for a particle have been investigated in order to
be able to predict particle behaviour adequately. These
studies involve both particle behaviour in settling tank
conditions and particle behaviour in a generic swirling
flow field. The equations of motion involve various flow
phenomena resulting in forces on particles.
One of the forces involved in particle tracking cap-
tures the effects of the history of the particle movement.
This history term is often neglected, not only because the
term is often deemed very small, but more importantly
because calculation of the force involves an integral over
the time from the start of the motion. Such an integral
is computationally very expensive, especially if the time
steps are very small. However, present research shows
that the history term cannot be neglected in the case of
oil droplets in water. This force will generally increase
the predicted settling time, both in conventional settling
tanks and in swirling flow separators. This will result in
either more, or larger, settling tanks, or longer swirling
flow separators.
First the equations of motion for a particle will be de-
scribed. Subsequently solutions for these equations of
motion will be compared with experimental data from
Mordant and Pinton (2000). The data will be used to
compare the effect of several kernels that are used in the
history term. Empirical coefficients that appear in the
kernels have been varied to obtain a proper behaviour of











the history term. Subsequently results will be presented
of a comparison of numerical results with numerical and
experimental data obtained by Candelier et al. (2005);
Candelier (2008). Finally results will be shown for a
generic swirling flow field and conclusions will be pre-
sented about the history kernel and its importance for
water treatment using swirling flow fields.


Equations of motion for a particle

The motion of a particle is dictated by the velocity due
to forces acting on that particle. Newton's second law
states that d = F. The change in time of the posi-
tion of a particle is then the velocity of the particle. The
torque on a particle is related to the rotation of the parti-
cle, coupled by its moment of inertia: dIP- T. Mak-
ing the assumption that no mass transfer is occurring and
the moment of inertia is constant, the differential equa-
tions for the position, velocity and rotation of a particle
can be written as:


dv
Spdt =m volumee + fVstress d+ rag + fSaffman
+fMagnus + faddedmass + history + unknown) (3b)


dwp T
dt LI


(3c)


where the forcing term per unit mass is split into several
additive contributions. The last forcing term captures
all forcing terms currently not taken into account. Us-
ing known relations for the forcing terms the equations
of motion can be rewritten. The unknown forces will be
neglected. However, unknown forces can have a major
impact on the motion of particles. This is for example
observed for a particle settling along the axis of rotation
in a solid body rotation, e.g. see Candelier et al. (2005);
Candelier (2008).
For the derivation of the forces, it will be assumed that
oil droplets can be approximated by solid spherical par-
ticles. As the particle Reynolds number will generally
be low for settling conditions, deformations due to parti-
cle induced flow will be negligibly small. The particles
used in settling experiments are solid spheres, therefore
no assumptions have to be made with respect to geome-
try or substance.
The volumetric forces such as gravity are written as:

volume = f (4)

The forcing due to the divergence of the stress of the
undisturbed flow field at the location of the particle V o


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


has been rewritten using the Navier-Stokes equations for
one-way coupling in non-conservation formulation:

Vao 1Du
fVstress- f, (5)
PP Dt
with the material derivative defined as:
D a
D +u.V. (6)
Dt Ot
The effect of the drag is written as a function of the
drag factor f, the particle relaxation time TV and the rel-
ative velocity w.

fdrag -- (7)
TV
A larger particle relaxation time implies more slip and
therefore smaller influence of the drag on the particle
acceleration. The drag factor is related to the drag co-
efficient by Stokes' drag. It shows the influence of non-
potential flow effects on the drag:

Re C'
f = D 2- 6- fStokes 1, (8)
24 CD, Stokes

with the Reynolds number Re:

Dlwl
Re (9)
Vc
In the present research the relation for the drag factor
found by Clift and Gauvin (1970) is used, although the
relation for the drag factor found by Schiller and Nau-
mann (1933), a good approximation for Re < 800 ac-
cording to Crowe et al. (1998), can be used for most
cases considered. The drag factor, from the drag coeffi-
cient found by Clift and Gauvin (1970), is:


f 1 + 0.15Re0 68


0.0175Re
1 + 4.25 x 104Re 116 (10)


Figure 1 shows four common expressions for the drag
coefficient. It can be noted that Stokes drag underesti-
mates the actual drag for all Reynolds numbers, there-
fore f > 1. However, for Re < 1 Stokes' drag is a very
good approximation. The drag factor can be corrected
for the case of a non-solid particle, i.e. a droplet or a
bubble. For the present research such a correction is not
used, as the viability of such a correction for the non-
Stokes regime does not appear to be proven. Moreover,
in case of contaminated bubbles such a correction would
underestimate the drag force on a bubble significantly.
Also, in the case of oil droplets in water, the viscosity of
the droplet is generally larger than the viscosity of water
and therefore the correction is quite small.
The Saffman and Magnus lift forces capture phenom-
ena related to rotation of particle and flow field. Figure







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


10


10-1


101 0'5
100 10
Re

Figure 1: The drag coefficient CD for several expres-
sions found in literature.


FL,S


U


Figure 2: Saffman force effect(left) and Magnus force
effect(right).


2 shows the two flow phenomena taken into account by
these forces. The Saffman lift captures the effect of a
higher fluid velocity on one side of the particle with re-
spect to the other side of the particle due to a flow field
velocity gradient. This effect causes a pressure differ-
ence over the particle resulting in a lift force that de-
pends on the direction of the rotation and the direction
of the relative velocity. The Magnus lift force represents
the effect of the flow field being curved by a rotating
particle. Rotation of the particle causes a change in an-
gle between the flow ahead of the particle and the flow
in the wake of the particle. This is caused by a force
of the particle on the flow and vice versa. Crowe et al.
(1998) use the relative rotation for the Magnus lift, while
Loth and Dorgan (2009) use the absolute rotation of the
particle. For the cases studied, both methods have been
investigated, and differences turn out to be small in cur-
rent research. An approach using the absolute rotation
will therefore be used.
The Saffman lift force per unit mass is:

9.66 n^
fSaffman LS x w6, (11)
sann Di) V~ I W| 0


with the rotation of the flow field:


W = VX U.


The Saffman lift coefficient CL can be written as, see
Mei (1992); Crowe et al. (1998):


L,S o0.33143 (1
CL, 0524 Re
0.0524 7Re


Ro+ \
Re
+e 10


Re < 40,

Re > 40,


where
D I
= 0.005 < 3 < 0.4, (14)
Iw|
This is a correction for the lift force which in the limit
of low Reynolds numbers, for which the influence of
B becomes negligible, reduces to the relation found by
Saffman (1965), and for high Reynolds numbers be-
comes the relation from numerical results by Dandy and
Dwyer (1990).
The Magnus lift force per unit mass is:
31
fMagnus 4 3 CL,MWp X W. (15)

The Magnus lift coefficient CL,M is a Reynolds number
and spin correction on the Magnus lift as derived by Ru-
binow and Keller (1961). Loth and Dorgan (2009) use
a relation that gives the lift coefficient in the same order
of magnitude as measurements. This relation is:


CL = 1 +0.15) tanh (0.18Re)

-0.15tanh(0.28(a- 2))tanh (0.18vRe), (16)

with
c =D |p (17)
w|
When a particle is moving with respect to the flow,
a similar volume of fluid has to be moved in opposite
direction. In case of an acceleration of either the contin-
uous phase, or the dispersed phase, a resulting force will
remain. This effect is written as:


C 1Du
added mass Cam
p Dt


dv
dt)


The expression for the added mass force coefficient pro-
posed by Auton et al. (1988) is Cam = 1, and is ac-
curate for a fairly large range of Reynolds numbers, see
Loth and Dorgan (2009). There has been some debate
whether or not the use of the material derivative is cor-
rect, but general consensus appears to be that the ma-
terial derivative should be used. However, in case the


3


+ L,M







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


particle velocity does not differ too much from the fluid
velocity, the derivative of the flow velocity along the par-
ticle trajectory is a good approximation of the material
derivative.
The history force will be presented in more detail after
introducing the equations of motion. It can be written as
an integral over time, having the dimensionless kernel
K (t, T), with T the integration time variable:


fh T K dT.
-C


Using the relations for the forcing terms, neglecting
the unknown forces, and rewriting parts for convenience,
the equation of motion for a particle used in present re-
search is written as:


dv
(p+ Cam)
dt


SPf
(P 1) f ---w
TV


Du 3
+ (1 + Ca,) Du+ 3-CL,Mp X W
9.66 L
DCLS1 W X Ww
7nD V FW o


t

TV


din
K (t, T) dw (20)
dT


The history term is related to Stokes' second problem
for the boundary layer on a flat plate oscillating in the
direction tangential to the plate. The drag force in equa-
tion (7) is based on a developed flow field around the
particle. However, the boundary layer will only be fully
developed after a certain time span. To correct for this
phenomenon the history term will have to be included,
describing the non-instantaneous effects. The kernel that
was originally derived for this force, also called the Bas-
set force, will be called the Basset kernel and can be
written as:

d1 Td 1
KBasset (t, T) = (21)

with the diffusive time-scale:


The kernel is based on the solution for the case of
constant acceleration in the low Reynolds limit. This
solution involves the derivative of an error function. For
a variable acceleration infinitesimal steps can be com-
bined to the resulting kernel with a 1 behaviour.
It was observed that the decay of influence of the ac-
celeration over time was more rapid for the long time


100
(t-t)/d [-]


Figure 3: Reproduction of figure in Loth and Dorgan
(2009) of the history kernel with cl = 2.5 and
C2 0.2 for three Reynolds numbers. The
history time TH is as defined in equation (25).


behaviour than is predicted by the Basset kernel. There-
fore alternative kernels have been proposed. An alter-
native kernel, proposed by Mei et al. (1991); Mei and
Adrian (1992), has a square root decay for short time be-
haviour, and a more rapid quadratic decay for long time
behaviour:


K (t, T)


1 +


KBasset (t, T)

( Re
KBa ,e(t,T))


1] cl (23)

TW c


with the history factor fH as:

fH = (0.75 + caRe)3,


and cl and c2 dimensionless empirical coefficients that
have been taken cl = 2.5 and c2 0.2 by Loth and
Dorgan (2009), using several sources of settling experi-
ments. However, Mei and Adrian (1992) originally pro-
posed cl = 2 and c2 0.105 on semi-analytical basis.
In the line of Loth and Dorgan (2009) a history time is
now defined to indicate the distance in time from time t
at which the influence of the square root behaviour is
equal to the influence of the quadratic behaviour, i.e.
when the denominator in equation (23) is equal to 2.
Therefore:

H fH Re2 (25)
7 Re
The history time is only influenced by the choice of the
parameter c2, the Reynolds number Re and the particle
diffusive time-scale Td. In figure 3 the history kernel is
presented including the history time. For the limiting
case of Re = 0 the history time would be infinite, as











no transition to the quadratic decay regime will occur.
The parameter cl does not influence TH as defined in
equation (25). However, cl does influence the size of
the transition region of the curve between the square root
and quadratic decay regime, i.e. a larger cl results in a
larger transition region.
The kernel was subsequently refined by Lawrence and
Mei (1995), leading to a drag factor depending kernel,
with history factor fH:

3
fH (t, T) 3 (f (t) + Re (t) f' (t)) f' () (26)
4
For low Reynolds numbers Lawrence and Mei (1995)
used the drag factor

f =1 + Re. (27)
16
It can be derived that using this expression for the drag
coefficient in equation (26), the drag coefficient based
history factor and the empirically determined history
factor in equation (24) only slightly deviate in values up
to order 0 (Re).
where f' is the derivative of f to Re. Lawrence and
Mei (1995) also found that the quadratic time compo-
nent can be described by an integral of the velocity from
T to t. In present research, after comparison of numer-
ical results with experimental data, the original expres-
sion for the history kernel is here adapted with an added
drag factor correction leading to:


K (t, T) 4fH (t, T)
7T


x Re (t)j(t, + T (t,r) (28)

with the history factor fH as defined in equation (26)
and dimensionless time T as:


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


there will be some history prior to to, the time at which
initial information is known. Therefore the integral is
split into a part for which data can be obtained, either
exactly or by interpolation, and a part for which data is
not available:


t

dw
-00O


to t
L dw f dw
K dT dr+ K dT. (30)
dT J dT
-CO tO


Subsequently the first part of the integral can be
rewritten using integration by parts, as


to
K dT
j dr
-00


to
f dK
[Kw] o,- J w dT.


The integral in the right hand side is still unknown.
However, the derivative of the kernel has a larger rate
of decay than the kernel itself, while relative velocities
can be expected to be near zero. The kernel times the
velocity at minus infinity will be zero, and therefore the
integral over the entire time domain can be written as if
the particle obtains its velocity instantanuously at t = to:


t
f dw
K dr
-00


SKdw dT + Kowo.
S dr


The torque on the particle is also required for solu-
tions of the equations of motion of the particle. The
torque on a particle tends to force the rate of rotation
of the particle to be the same as half the rate of ro-
tation of the flow. In Crowe et al. (1998) a relation
for the torque vector is given, and subsequently a rela-
tion for the torque in a quiescent fluid depending on the
Reynolds number is provided. It will be assumed that
the Reynolds number correction will be approximately
valid for the relation with a fluid flow, leading to:


SRe (') ,
Tdf (T) d .


The expressions in equations (28),(29) differ from the
equation proposed by Lawrence and Mei (1995) due to
added drag coefficient corrections and a smooth switch
from the square root to the quadratic decay regime.
These corrections are based on a comparison of numeri-
cal results with experimental results of Mordant and Pin-
ton (2000), this will be shown in more detail in the next
section.
It must be noted that the lower limit of integration for
the history term is minus infinity. If the particle and fluid
are initially, at time to, at rest, this is not a problem as
the contribution prior to to is zero. If this is not the case,


T = 2.01/D3 cwe


where:
1 )


p (1 + 0.201 Re,


p D2
4)P 7we[ D2


Analytic solution of the equations of motion is only
possible in very specific situations and often the assump-
tion of Stokes flow. For a more generic case a numerical
integration scheme is required to obtain an approximate
particle trajectory. In the present research two numerical
methods have been used. An explicit algorithm using the
classical Runge-Kutta scheme, which is globally fourth
order accurate in time, is used for most calculations. It
is especially useful for the comparison of the numerical
solution with experimental data and fully resolved flow











solutions. The Runge-Kutta scheme has good capabil-
ities for adaptive time-stepping which is used for some
of the calculations. For particle tracking in a generic
swirling flow field, an implicit Crank-Nicholson scheme
is used, which is globally second-order accurate in time.
The major advantage is that in this scheme there is no
strict limitation to step size for a stable solution. How-
ever, if flow time scales are small, solutions generally do
not show these high frequency effects. A major disad-
vantage of an implicit scheme is the requirement to solve
a system of equations. In the case of the present equa-
tions, this will involve a non-linear system of equations.
However, iterative methods can generally approximate
solutions quite efficiently with only a few iterations.
The history term requires a numerical integration
scheme. A Simpson quadrature method is used to ob-
tain a proper order of accuracy. However, one of the
major problems for a standard integration scheme is the
singularity at T t in the kernel for both the Basset ker-
nel and the alternative kernel. One of the solutions for
this problem, which is only possible for the Basset force,
is to write the entire integral as a fractional derivative,
see e.g. Bombardelli et al. (2008); Coimbra and Rangel
(1998). The fractional derivative makes use of the divi-
sion by the square root, and therefore cannot be used for
the kernel involving quadratic decay. Another solution,
also only applicable to the Basset kernel, is to assume a
constant acceleration in the final time step and use the
analytic integral for the Basset kernel for this situation.
For the alternative kernel it can be assumed that in the
limit of very small time steps (t < t/ which
implies that the kernel reduces to the Basset kernel in the
final step. This makes the second method, in which an
analytic solution for the Basset kernel is used valid for
the calculation of the history term. Therefore the method
used in the present research for calculation of the history
term will have an altered integral:


dw
/K d dr
j dr
-CO


t -t
f dwi
lirn Kdw dT
StRo J dT
to


dw
dt
+ d- | .KBasst drT+ Kowo. (35)
2 t-6t

Settling sphere

In 2000 Mordant and Pinton (2000) published experi-
mental results of experiments of spheres settling in qui-
escent water. These experiments have been used to adapt
the kernel in the history term in such a way that par-
ticle behaviour is predicted fairly accurate for a large
range of Reynolds numbers. The results show that pre-


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


dictions without a history term overestimate the initial
settling velocity and underestimate the time required for
the sphere to reach its terminal velocity. On the other
hand, predictions using the Basset kernel significantly
underestimate the long term velocity and overestimate
the time required for the sphere to reach terminal veloc-
ity. Both predictions sometimes overestimate the veloc-
ity of the sphere in part of the time domain for which
the Reynolds number is high. Therefore a kernel should
not only show a more rapid decay than the Basset kernel,
but it should also result in larger values than the Basset
kernel in specific situations.
Mordant and Pinton (2000) calculated the history
force for a specific case. The calculated force is the re-
sult of the force determined by the acceleration of the
particle, with known forces such as the drag and buoy-
ancy substracted. Figure 4 shows the resulting history
force. It includes the history force caclulated by Mor-
dant and Pinton (2000) using the Basset kernel. It was
observed that while the peak force is of the proper or-
der of magnitude, however, the long term decay of the
Basset kernel is much too slow. This case has been
used to compare the numerically obtained history force
with experimental data. It has been assumed in the nu-
merical simulations that the experiments have been per-
formed in Lyon. Using the altitude and latitude of Lyon
the gravitational acceleration has been calculated to be
g z 9.807 [ms 2]. Mordant and Pinton (2000) de-
termined water conditions to be conditions at 25 [ C].
At that temperature, density and kinematic viscosity are
pc 997 [kgm 3] and v = 8.9 10 [m2 s- ], re-
spectively.
The numerically determined history force using the
kernel as defined in equation (28) is also shown in fig-
ure 4 for a visual comparison. The general trend of the
force is predicted correctly, although the peak value of
the force is a bit lower than that of the experimentally
determined force. In figure 5 the history force predicted
using several kernels are shown. While the Basset ker-
nel overestimates the force significantly, the kernels pro-
posed by Loth and Dorgan (2009) and the one proposed
by Mei and Adrian (1992) predict a lower history force
in the part of the domain where accellerations are rela-
tively large.
Another quantity that can be used to compare numeri-
cal results with experimental results is the velocity. Usu-
ally the velocity time history is provided in experimental
data found in literature. Mordant and Pinton (2000) pro-
vide velocity histories for several cases. The first case
considered here is the case for which the history force
was determined. Figure 6 shows the velocity history
for this case including the velocity predicted using the
kernel of equation (28). While the numerical prediction
slightly overestimates the velocity reaching terminal ve-







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


z ........4 .. ..... . ......" --............t-......-







0 0.1 0.2 03 04 0.5 0.6

time (s)

Figure 4: Experimentally determined history force, ir-
regular solid black line, and predicted history
force using the Basset kernel, dashed line, by
Mordant and Pinton (2000). Prediction using
the kernel of equation (28) is added for com-
parison, solid blue line. pp 7850 [kg m3],
D 1 10-3 [m].








x10-6
7
-Basset
Mel et al
SLoth and Dorgan
6 . Equation 28


4

S3

2





00 01 02 03 04 05 06
t[s]

Figure 5: Numerical results of the history force for sev-
eral kernels. pp 7850 [kg m-3], D
1 10-3 [m].


600


Figure 6: Experimental obtained velocity profile, solid
black line, numerical predictions without his-
tory force, dashed line, and with history force
using the Basset kernel by Mordant and Pin-
ton (2000). Prediction using the kernel of
equation (28) is added for comparison, solid
blue line. pp 7850 [kg m-3], D 1 1
10-3 [m].



locity, and slightly underpredicts the initial velocity, re-
sults are quite satisfactory.
Another case from the experiments has been selected
to present a case for which a sphere attains a high
Reynolds number during the settling process. The high-
est Reynolds number that was observed in the experi-
ments was Re = 7700. This Reynolds number implies
that the drag coefficient is in the Newton regime, the
regime with a near constant drag coefficient. Figure 7
shows the result of the experiment. It can be noted that
the experimentally obtained velocity is lower for part of
the domain than the calculation both without and with
the history force, using the Basset kernel. Numerical
prediction using the kernel of equation (28) does not
show the initial overprediction of the velocity. Due to
the drag correction in the kernel, the history force can
obtain larger values than when using the Basset kernel.
However, the current kernel fails to correctly predict the
particle velocity near its terminal velocity and further in-
vestigation is required.
One of the causes for different particle behaviour
could be the kinematic viscosity used in the numeri-
cal calculations. The Reynolds number listed by Mor-
dant and Pinton (2000) can be obtained using either a
kinematic viscosity of v, 9.02 10 7 [in2 s1] or
a diameter of D 5.6 10-3 [m], while the kine-
matic viscosity used in numerical calculations was set
to v, 8.92 10 7 [m2s-1] and the diameter to
D 6 10-3 [m]. Although such a difference in kine-







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


0 200 400 600 800
time ( ms )

Figure 7: Experimental obtained velocity profile, solid
black line, numerical predictions without his-
tory force, dashed line, and with history force
using the Basset kernel by Mordant and Pin-
ton (2000). Prediction using the kernel of
equation (28) is added for comparison, solid
blue line. pp 7750 [kg m3], D 6 6
10-3[m].


matic viscosity or diameter seems minor, it has a sig-
nificant impact on the obtained Reynolds number and
subsequently on the obtained drag factor, equation (10),
and history factor, equation (26). As the flow is in the
Newton regime, a near constant drag coefficient implies
that the drag factor increases almost linearly with the
Reynolds number.
The above two cases are two important cases, one
for a relatively low Reynolds number, and one for a
very large Reynolds number. The kernel introduced in
the present research is able to predict particle behaviour
more accurately than the Basset kernel and agrees with
the experimentally determined history force. In the re-
mainder of this paper the obtained kernel, equation (28),
will be used in the prediction of particle behaviour in
rotating and swirling flow fields.

Swirling flow field

One of the properties of swirling flow, is the rotating
flow field. Therefore, a flow field consisting of purely
solid body rotation has been studied to further improve
the prediction of particle trajectories. In Candelier et al.
(2005) and Candelier (2008) such a case of solid body
rotation has been presented, and predictions of particle
behaviour using the present model will be compared to
results obtained by Candelier (2008).
A solid body rotation is described with a single pa-
rameter w for the rotation. Using this parameter the flow


can be described by:


S(


0WY
0o


Important properties of such a flow field for particle
tracking include the material derivative, the rotation of
the flow field and the Laplacian of the flow field:

DuW 2 Y ;wc 0 ;V2U =0,
0 2w
(37)
respectively. The material derivative is used in both the
added mass force and the description of the pressure gra-
dient and stress divergence forces. The rotation is re-
quired for the lift forces. When the Laplacian of the
flow field is not zero, corrections have to be added to
the forces, known as the Faxen force corrections, see e.g.
Crowe et al. (1998). As the Laplacian of a solid body ro-
tation is zero, such corrections are not required and will
therefore simplify the equations of motion to the form
presented in equation (20). It will also be important to
notice that u x w = -2.D
Candelier (2008) studied two situations for the case of
solid body rotation. The first case is a particle released
in the flow field with gravity working in the direction of
the axis of rotation. The second case, only numerically
studied, is a particle released in the flow field with grav-
ity working in the direction perpendicular to the axis of
rotation. The first case shows that the rate of rotation has
an influence on the particle settling velocity. In equa-
tion (20) all terms involving rotation only act in a plane
perpendicular to the axis of rotation, and therefore the
currently used equation fails to describe this particle be-
haviour.
In the second case gravity is acting in the plane per-
pendicular to the axis of rotation. There will be no move-
ment in a direction parallel to the axis of rotation, and
therefore an influence of the rotation to the drag in the
direction parallel to the axis of rotation will not influence
results. If a particle is released at the axis of rotation in
such a flow field, the particle will be forced away from
the axis of rotation due to gravity. However, centrifu-
gal forces, caused by the flow field, will force the parti-
cle back to the axis of rotation. In the case studied, the
forces acting on the particle will cause the particle to at-
tain equilibrium at a specific position with respect to the
axis of rotation. The location of that position will first be
predicted using the assumption of steady state. The par-
ticle velocity and its derivatives will then be zero at this
equilibrium point. The steady state of the particle also
implies that rotation of the particle will be equal to half
the rotation rate of the flow field, because the torque on







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


the particle will than be zero. The reduced equation of
motion for the particle in terms of x, y, w, g, with gravity
working in negative y-direction will then become:


0(1


w2[A ( -Y
XA -


B(x


)13
(38)


with the dimensionless factor A (x):

A (x) W T
A=TV


and dimensionless factor B:

3 9.66 v
B(x)= 1 +Ca ,M -"2 9 L,S c. (40)
4 7D 2w

Equation (38) provides a relation for the equilibrium
location of the particle:


Depending on the sign of A, B and 1 p, the point
where steady state is attained is within a specific quad-
rant with respect to the axis of rotation. A can never
be zero, and therefore the equilibrium point can only be
on the line x = 0 if the density ratio = 1. This im-
plies that the equilibrium point can only be at the origin
if the equilibrium point is located at x 0. Candelier
(2008) used as parameters D = 3 10-3 [m], = 0.1,
v, = 2 104 [m2 s ] and w 15 [rads- ]. Gravity
is chosen to be the same as in the numerical simulations
for a settling sphere. Using an iterative approach, the
equilibrium location for this situation is then obtained
as:
X x ( 0.45
XD D 0.067 ) (42)

However, the location of the equilibrium position shown
by Candelier (2008) has a positive value for the y coor-
dinate. The sign of B determines the sign of the y coor-
dinate, and if lift forces are taken into account in such a
way that B > 0 the y-coordinate will become positive.
In figure 8 the particle trajectory predicted by the
equations of motion in equation (20) is presented. The
trajectory converges to the equilibrium point predicted in
equation (42). In figure 9 the result of a simulation tak-
ing the history force into account has been compared to
the result without taking the history force into account.
The history force decreases the maximum excursion of
the particle in x- and y-direction. Therefore, the particle
trajectory without history force will initially be a curve
with a larger radius than the trajectory with history force.
After some time however, the trajectory will converge to


-05 -04
XD[-


00


Figure 8: Prediction of a particle trajectory in a solid
body rotation flow field with gravity acting
in a plane perpendicular to the axis of rota-
tion. D = 3 10-3 [m], 0.1, vc = 2
104 [m2 s1] and w 15 [rads- ]. x and
y coordinates have been divided by the par-
ticle diameter D to obtain non-dimensional
variables.


-08 -07 -06 -05 -04
XD[-


00


Figure 9: Predictions of particle trajectories in the solid
body rotation flow field also used in figure
8 with and without taking the history force
into account. x and y coordinates have been
divided by the particle diameter D to obtain
non-dimensional variables.


_ ( A (x)) 1-
-B (x) W 2 (A2 (x) + B2 (X)) 9.







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


Table 1: Summary of fluid, flow field and particle prop-
erties used in computations for swirling flow.
Variable Value Units

U 2 ms 1
F 2 m2s 1
TV 500 s
Pc 997 kgm 3
vc 8.9 10-7 m2s 1
pp 800 kg m3
D 1 i 10-4 m
Ro 5 10-2 [m]



the equilibrium position faster than the trajectory pre-
dicted taking history force into account.
Forces in the direction parallel to the axis of rotation
caused by rotation still need to be investigated. How-
ever, present research first focuses on the effect of the
history force for particles in a swirling flow field, i.e. a
flow field with both rotation and an axial velocity. To in-
vestigate the influence of the history force for such flow
fields, a generic swirling flow field is described as:


U (


(1
(1
2rrr2 i


2 )
2 )


where F is a parameter determining the circulation of
the vortical flow with dimension [m2s-1], and T, is a
parameter determining the decay in time of the circula-
tion away from the axis of rotation with dimension [s].
The rotating part of the velocity is a Lamb-Oseen vortex,
however to provide a steady flow field the decay in time
has been frozen by taking T, constant. The parameters
for the performed simulations are summarized in table
1.
The parameters for the flow field and the particle have
been chosen to represent conditions similar to swirling
flow produced water treatment. During such a separation
process oil droplets up to 1 10-4 [m] are found. As
a reduction in diameter significantly increases the time
for an oil droplet to reach the axis of rotation, a droplet
size has been chosen close to the upper limit of often
occurring droplet diameters.
Figures 10 and 11 show the resulting trajectory of a
droplet, predicted taking the history force into account
as well as the trajectory of a droplet without taking the
history force into account. The droplet is released at
x Ro with the velocity of the flow field at that lo-
cation. Distances have been divided by the initial radius
with respect to the axis of rotation to obtain dimension-
less distances. The trajectories follow the rotating flow


25---- ---------- I
S-- With history
I---- Without history

20-





20 - --- ---





00 01 02 03 04 05 06 07 08 09 10
R-1

Figure 10: Dimensionless distance travelled in axial di-
rection versus dimensionless radial position
of the particle trajectory for a swirling flow
field. Radius and axial distance have been
divided by the initial radius Ro to obtain
non-dimensional variables. Table 1 summa-
rizes the flow conditions. Density ratio of
particle and fluid p 0.8


1 00 L
-1 00


NxlO

-080 -060 -0 40 -0 20 0 00 020 040 060 080 1 00
X-1


Figure 11: Particle trajectories projected on the plane
perpendicular to the axis of rotation. x and y
position have been divided by the initial ra-
dius Ro to obtain non-dimensional variables.
Table 1 summarizes the flow conditions











and gradually move to the axis of rotation. Initially the
effect of the history force is not very large. However,
the distance required for the particle to move within one
tenth of the initial radius shows a difference of about
two times the initial radius. This case shows that a sig-
nificant influence of the history force for swirling flow
separation is to be expected.


Conclusions

To provide an adequate description of the trajectories of
particles in separation processes a Lagrangian particle
tracking model is used. First the settling of particles in
quiescent fluids is investigated, numerical predictions
are compared with available experimental results. The
results showed the importance of a history force acting
on the particle. A square root decay, resulting from the
Basset kernel, is not able to describe the quadratic decay
found for the long time behaviour of the history force.
Therefore, a kernel initially proposed by Lawrence and
Mei (1995) is adapted to provide a kernel for the history
term with a smooth transition from square root decay to
a quadratic decay. The kernel used in the present re-
search depends among others on the drag factor, and is
therefore able to adequately predict particle behaviour
for a large range of Reynolds numbers However, in
some cases predictions of the particle reaching terminal
velocity do not agree completely with experimental data.
The influence of the history term is subsequently in-
vestigated in a solid body rotation flow field. A numeri-
cal study performed by Candelier (2008) is used to pro-
vide for a case of a solid rotation flow field. For the case
of solid body rotation flow with gravity acting parallel
to the axis of rotation a significant influence of rotation
rate on the particle trajectory can be found. This cannot
be predicted by equations of motion, used in present re-
search, that include only Saffman and Magnus lift forces
as rotation dependant forces. However, the observed ef-
fect will be relatively small for swirling flow, as axial
velocity differences will generally be small.
The second case studied is the case of solid body ro-
tation flow with gravity acting perpendicular to the axis
of rotation. Simplification of the equations of motion for
a steady state situation provides a prediction of the point
of equilibrium for a steady state situation. Both the par-
ticle trajectory taking history force into account and the
particle trajectory without taking history force into ac-
count were found to converge to the predicted equilib-
rium point. However, history force significantly influ-
ences the maximum excursion of the particle as well as
the time in which the particle converges.
Finally the case of swirling flow is investigated. The
history force shows a significant influence under condi-
tions occurring during produced water treatment. Dif-


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


ferences in the order of pipe diameter can be observed
for the distance it takes to reach a radius will outside the
axis of rotation.


Acknowledgements

The research presented is made possible by FMC Tech-
nologies/CDS Separation Systems.


References

T. R. Auton, J. C. R. Hunt, and M. Prud'Homme.
The force exerted on a body in inviscid unsteady non-
uniform rotational flow. Journal of Fluid Mechanics,
197:241-257, December 1988.

Fabian A. Bombardelli, Andrea E. Gonzalez, and
Yarko I. Nino. Computation of the particle basset force
with a fractional-derivative approach. Journal of Hy-
draulic Engineering, 134(10):1513-1520, 2008.

F. Candelier. Time-dependent force acting on a particle
moving arbitrarily in a rotating flow, at small Reynolds
and Taylor numbers. Journal of Fluid Mechanics, 608:
319-336, August 2008.

F Candelier, J.R. Angilella, and M. Souhar. On the ef-
fect of inertia and history forces on the slow motion of a
spherical solid or gaseous inclusion in a solid-body ro-
tation flow. Journal of Fluid Mechanics, 545:113-139,
December 2005.

R. Clift and W.H. Gauvin. The motion of particles in
turbulent gas streams. Proceedings of Chemeca, 1:14-
28, August 1970.

C. F M. Coimbra and R. H. Rangel. General solution
of the particle momentum equation in unsteady Stokes
flows. Journal ofFluidMechanics, 370:53-72, Septem-
ber 1998.

Clayton T Crowe, Martin Sommerfeld, and Yutaka
Tsuji. Multiphase flows with droplets and particles.
CRC Press, 6000 Broken Sound Parkway, NW, (Suite
300) Boca Raton, Florida 33487, USA, 1998. URL
www.crcpress.com.

David S. Dandy and Harry A. Dwyer. A sphere in shear
flow at finite reynolds number: effect of shear on particle
lift, drag, and heat transfer. Journal ofFluidMechanics,
216:381-410, July 1990.

Christopher J. Lawrence and Renwei Mei. Long-time
behaviour of the drag on a body in impulsive mo-
tion. Journal ofFluidMechanics, 283:307-327, January
1995.







7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 -June 4, 2010


E. Loth and A. J. Dorgan. An equation of motion for
particles of finite Reynolds number and size. Environ-
mental FluidMechanics, 9(2):187-206, April 2009.

Renwei Mei. An approximate expression for the shear
lift force on a spherical particle at finite Reynolds num-
ber. International Journal of Multiphase Flow, 18(1):
145-147, January 1992.

Renwei Mei and Ronald J. Adrian. Flow past a sphere
with an oscillation in the free-stream velocity and un-
stready drag at finite Reynolds number. Journal ofFluid
Mechanics, 237:323-341, April 1992.

Renwei Mei, Christopher J. Lawrence, and Ronald J.
Adrian. Unsteady drag on a sphere at finite Reynolds-
number with small fluctuations in the free-stream veloc-
ity. Journal ofFluid Mechanics, 233:613-631, Decem-
ber 1991.

N. Mordant and J.-F. Pinton. Velocity measurement of
a settling sphere. European Physical Journal B, 18(2):
343-352, November 2000.

S. I. Rubinow and Joseph B. Keller. The transverse force
on a spinning sphere moving in a viscous fluid. Journal
ofFluid Mechanics, 11(03):447-459, November 1961.

P. G. Saffman. The lift on a small sphere in a slow shear
flow. Journal ofFluid Mechanics, 22(2):385-400, June
1965.

L. Schiller and A. Naumann. Uber die grundlegenden
berechungen bei der schwerkraftaufbereitung. Vereines
DeutscherIngenieure, 77:318-320, 1933.




University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - Version 2.9.7 - mvs