Paper No 7t" International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Lift and drag on agglomerates attached to walls
Robyn Larsen1, Dmitry Eskin2, Jos Derksen1
'Chemical & Materials Engineering, University of Alberta, Edmonton, AB, Canada T6G 2V4, rlarsen@ualberta.ca,
jos@ualberta.ca
2Schlumberger DBR Technology Center, Edmonton, AB, Canada T6N 1M9, deskin@slb.com
Keywords: agglomerates, simulation, latticeBoltzmann, immersed boundaries, shear flow, asphaltenes
Abstract
With a view to modeling asphaltene deposition on pipeline walls we perform direct numerical flow simulations revealing drag
and lift forces on spheres and agglomerates in shear flow attached to walls. The direct simulations are based on the
latticeBoltzmann method for solving the flow equations, combined with an immersed boundary method for representing the
noslip conditions at the solidliquid interfaces. For single spheres good agreement with analytical, low Reynolds number
results of lift and drag due to Leighton & Acrivos (1985) is observed. If scaled with the agglomerate size, nondimensional
drag on agglomerates is roughly 30% lower than for a solid sphere, and drag varies over a range of 30% without clear trends in
terms of the number of primary spheres in the agglomerate or its voidage. At low Reynolds numbers, dimensionless lift is
much weaker than drag and is relatively more sensitive to the (random) relative placement if the primary spheres forming the
agglomerate.
Introduction
of the computational domain on the hydrodynamic forces.
Asphaltene particle deposition in turbulent pipe flow is a
complex, multiscale process involving various
transportrelated phenomena: the formation of primary
particles, their growth as a result of agglomeration (Kusters,
1991), the migration towards and subsequent deposition of
agglomerates on the wall as well as the potential detachment
of agglomerates from the wall as a result of hydrodynamic
forces. This work focuses on the latter aspect as it may be
crucial in the prevention of agglomerate deposition.
We perform direct numerical simulations (DNS) of shear
flow around agglomerates attached to walls. The
agglomerates are composed of equally sized spheres. For
assessing the potential of the shear flow to detach the
agglomerates, in the simulations we measure the drag and
lift force on the agglomerates. We assume that the
agglomerates have a nearspherical shape, i.e. high fractal
dimension.
The latticeBoltzmann flow simulator (LBM) (Chen &
Doolen, 1998; Eggels & Somers, 1995; Somers, 1993)
combined with an immersed boundary method (IBM)
(Derksen & Van den Akker, 1999) to represent the
spherically shaped solidliquid interfaces was employed.
Prior to agglomerate simulations, the computational method
was calibrated and verified by single sphere computations.
In the first place, a calibration of the hydrodynamic radius
for the primary spheres was done according to the approach
of Ladd (1994). In the second place we used the
lowReynolds number analytical solution of Leighton &
Acrivos (1985) for the lift and drag force acting on a single
sphere attached to a flat in shear flow wall to validate our
numerical approach. This exercise also provided
quantitative information regarding the influence of the size
Nomenclature
a,aa
F,F*
L,L*
Lx,L,,L
n
ppP
Re, Rea
primary sphere radius, agglomerate radius
deformation tensor
drag force, nondimensional drag force
lift force, nondimensional lift force
flow domain dimensions
number of primary spheres in agglomerate
pressure, dimensionless pressure
Reynolds number, agglomerate Reynolds
number
moving wall velocity
u, u fluid velocity component, fluid velocity
magnitude
x,y,z Cartesian coordinates
Greek letters
y, yf overall shear rate, deformation rate in fluid
E voidage
v kinematic viscosity
p fluid density
Background on asphaltenes
Crude oil is a colloidal system consisting of saturates,
asphaltenes, resins and aromatics. Asphaltenes are colloidal
particles in the oil phase, a structure of which has not been
completely understood yet. The asphaltene particles
precipitate from oil due to depressurization below the
Paper No
socalled asphaltene onset precipitation pressure. In practice
the pressure gradually drops along a transport pipeline;
asphaltene particles precipitate and grow forming
agglomerates, which partially deposit on the pipe walls.
Asphaltene deposition is generally considered a negative
phenomenon that may cause significant technical problems
and economical losses (Hammami & Ratulowski, 2007). In
extreme cases, the asphaltene deposition may lead to
complete disruption of the pipeline. It is therefore important
to model the asphaltene deposition process to anticipate and
prevent (negative) critical process conditions. Attempts have
been made to provide methods and guidelines for predicting
and evaluating asphaltene flocculation and deposition in oil
reservoirs (Alkafeef et al 2005).
In this work, the focus is primarily on the hydrodynamic
aspects of asphaltene agglomeration and deposition. In
addition to surface effects, hydrodynamics is considered
crucial for understanding agglomerate deposition. Direct
numerical simulations of the simple shear flow around
spheres and agglomerates assembled of equally sized
primary spheres attached to flat walls have been performed.
Lift and drag on spheres and agglomerates have been
determined.
The dimensionless coordinates of the parameter space we
will be partly exploring in this paper are (1) a Reynolds
number, (2) the size of the asphaltene agglomerates
expressed in the number of primary spheres it is made of,
and (3) the void fraction of the agglomerate. The way we
define Reynolds number and void fraction will be discussed
below. In addition to physical parameters there are
numerical parameters (such as spatial resolution of the
simulations and computational domain size) the effects of
which will be explored as well.
Mesoscopic modeling of multiphase flow
Mesoscopic modeling can be defined in terms of the length
scales it represents. The Greek word mesos means middle.
The meso scale is considered to be in between the
macroscopic scale of the world we live in (for our purpose
this typically is the diameter of the oil pipeline), and the
microscopic scale of individual molecules. For multiphase
flow processes, meso scales relate to the scales of bubbles,
particles, and droplets and their mutual interactions
(Derksen & Sundaresan, 2007).
Mesoscopic modeling is usually based on some sort of
kinetic representation of the physics involved. Instead of
considering the continuum differential equations,
mesoscopic modeling considers a representation based on
(sometimes real, sometimes fictitious) particles that interact
according to simple rules ("particlebased modeling"). The
latticeBoltzmann method (LBM) is an example of a
mesoscopic method. It is (supplemented with an immersed
boundary method, IBM) the primary tool used in the
research presented here. The simple update rules and
locality of operations makes the LBM highly
computationally efficient and inherently suitable to be run
on (massive) parallel computer platforms, allowing for
simulations with high spatial and temporal resolution (Succi,
2001).
The LBM is a relatively novel way to do Computational
Fluid Dynamics (CFD). The fundamental idea of the LBM
71" International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
is to construct a simplified kinetic model that incorporate
the essential physics of mesoscopic processes in the fluid.
The kinetic equation provides many of the advantages of
molecular dynamics, including clear physical pictures, easy
implementation of boundary conditions, and fully parallel
algorithms (Chen & Doolen, 1998).
The kinetic model comprises fictitious fluid particles
residing on a uniform, cubic lattice. The fictitious particles
move in a specific set of directions, and collide to mimic the
behavior of a viscous fluid. The kinetic nature of the LBM
introduces three important features as identified by Chen &
Doolen (1998). These are: (1) the convection operation of
the LBM in phase (or velocity) space is linear; (2) the
incompressible NavierStokes equations can be obtained in
the nearly incompressible limit; (3) the LBM utilizes a
minimal set of velocities in phase space.
From an initial state the timestepping is achieved in two
sequential phases (a) streaming, where each fluid particle
moves to the nearest node in the direction of its velocity;
and (b) collision, in which particles arriving at a node
interact and exchange momentum according to scattering
rules (Chen & Doolen, 1998).
As discussed above, the LBM employs a uniform, cubic
grid. In order to impose noslip conditions at the surfaces of
solid, spherical particles we use an immersed boundary
method (IBM). The spheres are defined by a collection of
points on their surface with mutual spacing a little less than
the latticespacing. The noslip condition is achieved by
calculating forces acting on the liquid such that the linearly
interpolated velocity at the points defining the sphere's
surface gets zero (Ten Cate et al., 2002).
In our simulations, the spherical particles forming the
agglomerates have fixed positions. The forces determined in
the LBM are forces on the fluid required to maintain noslip.
The opposite of the sum of these forces can thus be
interpreted as the force exerted by fluid on the particles to
keep the latter stationary.
U0 _
T 
 >
2a
 a
L
F
a
z
t#x
Figure 1: Agglomerate consisting of equally sized primary
spheres of radius a attached to a wall in a simple shear flow
generated by a moving upper wall. Lift force (L) and drag
force (F).
A calibration procedure (as introduced by Ladd, 1994) is
needed to determine the effective radius of the spheres
(commonly referred to as the hydrodynamic radius). The
hydrodynamic radius is recognized with the symbol a and is
given in lattice units. In our study a radius a=6 has been
used. Grid refinement performed in previous studies
(Derksen, 2008; Derksen & Sundaresan, 2007) clearly
indicates that for low to modest spheresizebased Reynolds
numbers (up to order 10) an a=6 resolution is sufficient to
4 .
Paper No
capture the essential physics of spherefluid interactions.
Problem definition
We perform direct numerical simulations of the simple shear
flow around agglomerates consisting of equally sized,
spherical primary particles attached to walls (see Figure 1
for a typical flow system) and determine the lift and drag
force on the agglomerates, as well as the distribution of
hydrodynamic forces on the primary spheres throughout the
agglomerate. Shear is generated by placing a wall parallel to
the wall the agglomerate sits on at distance Lz and moving
that wall with constant velocity uo so that in steady state
no
y= The agglomerates considered have a nearly
Lz
spherical shape. Given the (small) size of asphaltene
agglomerates (typically, not larger than a few tens of
micrometers) and the slow flow at such close distances to
the pipeline wall, we only consider low Reynolds number
flows. The (nonzero) value of the Reynolds number does
matter given the inertial nature of the lift force and its
proportionality to Re for spheres attached to walls (Leighton
& Acrivos, 1985) and to Re away from walls (Saffman,
1965 & 1968; McLaughlin, 1993).
71" International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
center of the domain. They collide inelastically and
eventually form a spherically shaped agglomerate. For n>16
typically such agglomerates have a voidage of cE0.62,
where the voidage is defined as void space volume fraction
in the smallest possible sphere (with radius a,) fully
enclosing the agglomerate: =n Less dense
aa
agglomerates are made by contracting spheres with a
slightly larger radius (a + 8) and subsequently changing the
radius back to a.
In the agglomerate simulations the Reynolds number based
2a
on the agglomerate radius (Rea = aa ) was kept to a
v
constant (low) value of Re, = 0.09. Also dimensions of the
flow
domain
Lx Ly, L,
(streamwisespanwisewallnormal) relative to a, (i.e. the
L L, L
aspect ratios x ,y, ) were kept constant based on
aa aa aa
results of singlesphere simulations (see the next section).
The two variables in the agglomerate simulations were the
number of primary spheres in the agglomerate (in the range
2 to 64), and the voidage of the agglomerate (E in the range
0.62 to 0.82).
*
*.~ '^A*
o ..*. .
t o 
,, g o
", ',, 6 ,,
Lx
Figure 2: A sphere attached to a wall in a simple shear flow.
Direct numerical simulations were first preformed on single
spherical particles attached to walls (see Figure 2 for the
flow configuration). For this situation at low Reynolds
ka2
numbers (defined as Re = ) there is an analytical
solution of Leighton & Acrivos (1985). The lift and drag
solution of Leighton & Acrivos (1985). The lift and drag
force stemming from this solution are
L
= 9.22Re
vpya2
F
and V =32.1 respectively. The analytical solution
vpya2
considers an unbounded flow domain. We use the analytical
solution to validate our numerical procedure, and to
establish aspect ratios (computational domain size in
streamwise, spanwise, and wallnormal direction over
particle radius) for which the numerical solution approaches
independency of domain size. Once appropriate aspect ratios
(i.e. aspect ratios beyond which the hydrodynamic forces are
largely insensitive of the domain size) were found, we were
able to move from single spherical particles to agglomerates
attached to walls.
The agglomerates are made by randomly releasing a preset
number (n) of equally sized spheres in a threedimensional
domain (as shown in Figure 3), and attracting them to the
4 I:
",, .
Figure 3: Randomly distributed primary particles are
contracted to the center of a cubic domain to generate an
agglomerate.
Results & discussion
Forces on a single spherical particle on aflat wall
To validate our numerical approach we started with
simulations of single, spherical particles attached to a flat
wall in a simple shear flow at a low Reynolds number
(Re=0.012). The results of these simulations can be
compared to the analytical results due to Leighton &
Acrivos (1985).
Illustrations of the flow around a single sphere in terms of
velocity vectors and the pressure field are given in Figures 4
and 5 respectively. Dimensionless pressure p* is defined
asp 
p('L )2
As expected, the numerical results in terms of the lift and
the drag significantly depend on the size of the
computational domain, see Figure 6. In this figure we
introduce dimensionless lift and drag according to
L F
L F so that L =9.22 and
vpya2 Re
vpya2
F = 32.1 represents the Leighton & Acrivos solution.
Paper No
To get to a situation for which the lift and the drag are
virtually insensitive of domain size we first fixed the
domain size in streamwise direction (L,=64), and
wallnormal direction (L7=31) and systematically increased
the spanwise size (L,). For Ly>64 drag and lift get
insensitive of L, (see Figure 5, top panels). Subsequently we
fixed L, to 71 and checked the sensitivity with respect to L,.
At L, around 80 lift and drag tend to constant values. Then
the wallnormal distance was increased and at L=100 lift
and drag converged to virtually constant values of 7.8 and
34.6 respectively. These numbers are within 15% of the
analytical solution.

       
      
Figure 4: Velocity vector field of the flow around a single
sphere attached to a wall in simple shear flow.
 x
x p* m
00
Figure 5: Nondimensional pressure p* around a single
sphere attached to a wall at Re=0.012. Top panel: xzcross
section through sphere center; bottom panel: xy cross
section through sphere center.
Reasons for the deviation from the analytical solution relate
to the (still) finite Reynolds number and the relatively
modest spatial resolution: the sphere radius a corresponds to
6 lattice spacings. This resolution was chosen because of the
large numbers of spheres per agglomerate (and related large
71" International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
domain sizes) we will be dealing with in the agglomerate
simulations.
10
: 5
0 5 10 15 20
4y1e
39
36 *
33
30
0 5 10 15 20
Ly/a
S 30 *
* s o
*3 =,20 *
10
0 4 0
0 4 8 12 16 20 0 5 10 15 20 25
Lx/a L.',,
36.0
35.5
7.5 *
7.0
6.5 34.5
6 10 14 18 6 10 14 18
Lz/. Le/a
Figure 6: Normalized lift force L (left) and drag force
F* as a function of the size of the flow domain as defined
in Figure 2. The top panels have L,=64, L3=31, while Ly
varies. The middle panels have Ly=71, L,=31, while L,
varies. The bottom panels have L,=81, Ly=71, L, varies.
Re=0.012; a=6.
Agglomerates attached to aflat wall
Agglomerates are characterized by the number of primary
particles (n) they consist of, and by their porosity. In order
to determine the porosity the radius of the agglomerate aa
(the radius of the smallest sphere that fully encloses the
agglomerate) needs to be determined. Given an agglomerate
with n primary particles and voidage E, the flow conditions
are fully defined by the Reynolds number based on aa:
Re = ya (with a =a ); Re=0.09. The
v 1e
agglomerate simulations have constant domain size relative
to a,: ,, =(10.8,10.2,13.7). The singlesphere
aa aa aa
results indicated that for these aspect ratios the drag and the
lift do not depend much on the size of the flow domain.
As an example we now analyze the flow at Re=0.09 around
an agglomerate that consists of n=64 primary spheres and
that has a porosity of E0.75, see Figures 7 to 10. In Figure
7 we show the velocity field in the xzcross section. The
upper panel (showing the entire cross section) illustrates the
laminar nature of the flow. The large domain size in
zdirection (L,) is required to avoid too significant impact of
the upper wall on the lift and the drag force. The lower
panel of Figure 7 zooms in on the region around the
agglomerate and has a logarithmic velocity magnitude scale
so that it becomes clear that the velocities inside the
agglomerate are extremely small. The latter is reiterated in
the velocity vector plot in the xyplane in Figure 8.
3 5 0
7.0
Paper No
S 0 0 0 0 0
r
Figure 8: Velocity vectors in the xyplane at z=4a above the
bottom plate. Agglomerate with n=64, e=0.75.
x I 
0 o0
0 0
Figure 7: Velocity magnitude contours in the xz cross
section around an agglomerate with n=64 primary spheres,
and e=0.75. Top panel: the entire xz cross section & linear
color scale; bottom panel: lower portion of cross section &
logarithmic color scale.
The pressure distribution around (i.e. outside) the
agglomerate (Figure 9) is very similar to the distribution
around a solid sphere (Figure 5). Inside the agglomerate
small irregularities in the pressure field can be observed that
are the result of the specific relative placement of the
primary spheres in the agglomerate. Finally, liquid
deformation (and thus viscous stress) inside the agglomerate
is very weak, see Figure 10. In this figure y' denotes the
generalized deformation rate in the fluid 3f = 2 ,
with d +I a, ) the deformation rate tensor.
2 ax ax
These flow visualizations indicate that the flow around this
agglomerate can be largely considered to be the flow around
a solid sphere (with a rough surface) of radius aa.
The results of all our simulations in terms of dimensionless
lift L* and drag F* are given in Figure 11. Dimensionless
drag is in the range 17 to 26, i.e. some 60% of the drag on a
solid sphere attached to a flat wall in shear flow. There is
not a clear trend with respect to the number of primary
spheres in the agglomerate or its voidage. The essential
parameter is the size of the agglomerate while its internal
.& x p* m
.oX P
On
ch C)
x p*
Figure 9: Nondimensional pressure p* around an
agglomerate with n=64, e0.75. Top panel: lower part of the
xzplane through the center of the agglomerate; bottom
panel: entire xyplane at z=4a.
71" International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
structure is a secondary effect. Due to the limited flow
through the agglomerates (Figures 7, 8) this is not
surprising.
At low Reynolds numbers, lift is a much weaker force than
the drag force, and lift will only marginally impact
detachment of agglomerates attached to walls. Compared to
solid spheres, agglomerates feel a lift force that is roughly
50% of the lift force a solid sphere with radius aa would
experience.
000.
..w;
S*
...................................
............. ...... :.:: =;
............
 
..= ... .......... .. ... ...:. .
   
.. .. . .. . .
Paper No
x C
Figure 10: Fluid deformation rate
the agglomerate with n=64, E0.75.
SX
x
C o
7f around and inside
0 0
a3
o I  0
0 13 26 39 52 65
Nam,.rw fPhniclt
+ 6064% 0 6567% 7173% 7475% a 7881%
30
e o a
o X
o x
15 '
0 13 26 39 52 65
NraheraffP.'ddi
+ 60M4% 656 7% 7173% X 7475% A 7s81%
Figure 11: Dimensionless lift (top) and drag on
agglomerates as a function of the number of primary
particles; voidages are as indicated.
Conclusions
Motivated by the potential role of hydrodynamic forces on
the detachment of agglomerates off solid walls we
performed simulations of the lift and the drag force on
agglomerates in shear flow over a flat surface. The
simulation procedure was first calibrated by modeling a
single, solid sphere and comparing the simulation results
with the analytical solution of Leighton & Acrivos (1985).
These calibration simulations revealed the significant role of
the size of the flow domain on the lift and the drag. For
sufficiently large domains our simulations showed fairly
good accuracy (errors in the 15% range).
Agglomerates, even the ones with loose structure (voidage
of 75%) showed limited flow through them and therefore
the drag and the lift force could be well characterized by the
agglomerate size.
Acknowledgements
Support of the Schlumberger DBR Technology Center is
gratefully acknowledged.
71" International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
References
Chen, S., Doolen, G.D. Lattice Boltzmann method for fluid
flows. Annual Rev. Fluid Mech., Vol. 30, 329 (1998)
Derksen, J. & Van den Akker H.E.A. Largeeddy
simulations on the flow driven by a Rushton turbine. AIChE
J., Vol. 45, 209 (1999)
Derksen, J.J. & Sundaresan, S. Direct numerical simulations
of dense suspensions: wave instabilities in liquidfluidized
beds. J. Fluid Mech., Vol. 587, 303 (2007)
Derksen, J.J. Flow induced forces in sphere doublets. J.
Fluid Mech., Vol. 608, 337 (2008)
Eggels, J.G.M. & Somers, J.A. Numerical simulation of free
convective flow using the lattice Boltzmann scheme. Int. J.
Heat Fluid Flow, Vol. 16, 357 (1995)
Hammami, A. & Ratulowski, J. Precipitation and deposition
of asphaltenes in production systems: A flow assurance
overview. In Asphaltenes, Heavy Oils, and Petroleomics;
Mullins, O. C., Sheu, E.Y., Hammami, A. & Marshall, A. G.,
Eds., Springer, New York (2007)
Kusters, K.A. The Influence of Turbulence on Aggregation
of Small Particles in Agitated Vessels. PhD Thesis,
Eindhoven University of Technology, Netherlands (1991)
Ladd, A.J.C. Numerical simulations of particle suspensions
via a discretized Boltzmann equation. Part 2. Numerical
results. J. Fluid Mech., Vol. 271, 311 (1994)
Leighton, D. & Acrivos, A. The lift on a small sphere
touching a plane in the presence of simple shear flow.
ZAMP, Vol. 36, 174 (1985)
McLaughlin, J.B. The lift on a small sphere in wallbounded
linear shear flows. J. Fluid Mech., Vol. 246, 249 (1993)
Alkafeef, S.F., AlMedhadi F., AlShammari, A.D. A
Symplified Method to Predict and Prevent Asphaltene
Deposition in Oilwell Tubings: Field Case, SPE Production
and Facilities, May, 126132 (2005).
Saffman, P.G. The lift on a small sphere in a slow shear flow.
J. Fluid Mech., Vol. 22, 385 (1965) & Correction. J. Fluid
Mech., Vol. 31, 624 (1968)
Somers, J. A. Direct simulation of fluid flow with cellular
automata and the latticeBoltzmann equation. App. Sci. Res.,
Vol. 51, 127 (1993)
Succi, S. The lattice Boltzmann equation for fluid dynamics
and beyond, Clarendon Press, Oxford (2001)
Ten Cate, A., Nieuwstad, C.H., Derksen, J.J. & Van den
Akker, H.E.A. PIV experiments and latticeBoltzmann
simulations on a single sphere settling under gravity. Phys.
Fluids, Vol. 14, 4012 (2002)
