Paper No 7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Modelling drop interactions with topological changes in liquidliquid dispersions
Fabian Krause and Udo Fritsching
Particles and Process Engineering Department, University of Bremen,
Badgasteiner Str. 3, 28357 Bremen, Germany
ufri @iwt.unibremen.de
Keywords: multiphase flow, Computational Fluid Dynamics CFD, Volume of Fluid VOF, droplet coalescence,
emulsification, topological changes
Abstract
Emulsification is an important dispersion step in many liquid processes for production of dispersed liquid systems. During a
membrane emulsification process the droplets are formed at the individual pores on a porous layer. New interface area is
produced thereby while also the new formed drops may get in contact. This drop surface interaction may lead to relative
movement and deformation (sliding), repulsion or coalescence of the droplets. Main influence of droplet coalescence processes
on the molecular scale is due to interface properties determined e.g. by surface active substances while also hydrodynamic
influences determine the interface interaction and contact. This paper presents a modification of the standard VolumeofFluid
model VOF that allows describing droplet interactions including contact, coalescence and rebound, depending on the local flow
conditions, droplet contact time and contact area.
Nomenclature
Emulsification of a liquid in another liquid is used in
processes as e.g. in the production of components in
pharmaceuticals, cosmetics or foods (Schramm 2006).
Conventionally, emulsions are produced in e.g. rotorstator
systems or high pressure homogenizers (Schubert 1989).
Here the liquid disruption occurs in the high energy shear
flow of the continuous phase to form new liquid entities
(droplets). An alternative liquid dispersion technique is
membrane emulsification, which works with much lower
energy input. During the membrane emulsification process
the droplets get formed at the pores of a porous solid structure
as well as during the flow through the pore structure inside
the membrane. The size distribution of the generated droplets
is determined mainly by the pore size and the pore geometry
as well as by the micro fluidic conditions in the structure,
caused by the integral outer process parameters (pressure
difference and velocity of continuous respectively disperse
phase flow) (Schroder 1999). The newly formed entities /
drops at nearby pores may get in contact, as to be seen in Fig
1. The droplet surface interaction may lead to relative
movement and deformation (sliding) of the interfaces, to
repulsion or coalescence of the droplets, respectively. Main
influence of this surface interaction process is on the
molecular level where the liquid interface properties are
determined e.g. by surface active substances (Bibette 1999).
Also hydrodynamic conditions during the liquid interface
contact determine the dynamical contact process and its
outcome (Danner2001).
A area (m2)
F force field (N)
g gravity (m s2)
G Gibbs enthalpy (J kg1 Kl)
p pressure (Nm 2)
r radius (m)
T temperature (K)
t time (s)
v velocity (m s')
x length (m)
Ca Capillary number
Co Courant number
Greek letters
a phase fraction
p viscosity (Pa s)
p density (kg m'3)
C surface tension (N m1)
Subscripts
A constant Area
m mixture
p constant pressure
T constant Temperature
Introduction
Paper No
piese 1
1
Ir~lertece
.phase 2
I h;L. ri
oriliinuous phae D
phnaB n+1
Figure 1: Phases in the computational model
Surface Tension
In a liquid/liquid emulsion the area of the interface between
the phases is relevant for the inner energy of the system
(Wedler 2007). For this the equation for the free Gibbs
enthalpy:
dG=I dT+ l dp(1)
T T)p,A 'p3 )T,A
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
following subsequent steps: The first step is the approach of
the droplets. Between the droplets the continuous phase gets
displaced. After a further approach only a thin film between
the droplets exist. This film ruptures, if a critical thickness is
reached. Estimations of the critical film thickness vary in
literature from 10 nm (Duineveld 1998) up to more than 100
nm (Witelski, Bernoff, 2000).
The behaviour of the approaching droplets is influenced by
the hydrodynamic forces due to the viscous effects of the
continuous phase, as well as by the intermolecular forces
(van der Waals, electrostatic, steric, etc) inside the surface
film. The balance of these forces determines the next
coalescence step: either rebounce of the droplets,
coalescence or remaining of the droplets in this state. In the
case of coalescence the liquid film between the droplets
ruptures and the two droplets fuse. The surface area tries to
minimize to reduce the Gibbs free energy (Eq. 2) and the new
droplet tries to reach a spherical form.
Shearinduced coalescence in emulsions has been
experimentally investigated by means of microscopic
observations by Burkhart et al. (2001). They found that the
coalescence efficiency is depending on the capillary number.
The capillary number Ca represents the relative effect of
viscous forces versus surface tension acting across an
interface between a liquid and a gas, or between two
immiscible liquids. It is defined as
is expanded to:
dp+ i T ,
.aAwT,p
aGp
A T,A
The additional term, the influence of the interface area to the
Gibbs free enthalpy is the surface tension c. A liquid system
always tends to a minimal inner energy, so it aims to
minimize the interface area. For emulsions this results in the
following effects:
* To prepare an emulsion, it is necessary to bring in at
least that amount of energy that is equal to the extended
surface multiplied by the surface tension.
A disperse fluid volume tends to a spherical form, due to
its optimum relation of surface area to volume.
When droplets get into contact, they try to minimize the
surface area by forming a bigger spherical drop (when
they are not influenced or hindered by other forces like
van der Waals, electrostatic, steric, etc).
Coalescence in emulsions
Emulsions are thermodynamically unstable, because of there
extend surface. The surfaces area tries to minimize to reduce
the Gibbs free energy. Coalescence is the effect of two
contact lines merging, leading to coarser drop size
distribution or total breaking of emulsions (complete phase
separation). This is a main destabilization effect in emulsions
(Bibette 1999). Coalescence of two droplets occurs in the
where g is the viscosity of the liquid, u is a characteristic
velocity and C is the surface or interfacial tension between
the two fluid phases.
The coalescence process of two droplets with a diameter of
several gm pending on the tips of capillaries has been
investigated by Danner (2001). He found for surfactant
stabilized droplets contact times in the range of several
seconds.
The destabilization process of an emulsion in a flow field has
been investigated in a micro fluidic device by Bremond et al.
(2008). They derived that the coalescence occurs during the
separation phase and not during the impact.
Integral approaches to describe the effect of coalescence
processes onto the droplet size distribution in a liquid/liquid
dispersion have been derived within stochastic methods like
the Bubble Population Balance Equation (BPBE) described
by Chen et al. (2005). With the Volume of Fluid (VOF)
Model developed by Hirt and Nichols (1981) it is possible to
resolve phase boundaries of droplets. Delnoij et al. (1998)
and Koebe (l"'rl described the coalescence of bubbles
based on the VOF model. For bubbles that coalesce, these
authors got good agreement with experimental data, but the
model is not able to predict possible rebouncing of the
droplets. A model for the viscous coalescence of amorphous
particles is described by Garabedian and Helble (2000). They
used a geometric function that approximates the surface of
the coalescing particles. The surface evolution and the fluid
dG(GI dT+
aT)p,A
Ca= U
Paper No
mechanics behaviour were decoupled, to reduce the
computational effort.
Liquid surface interaction scales
The relevant length scales within liquid surface interaction
are (as shown in Fig 2):
physicochemical interaction of surfactants
molecules at the liquidliquid interface e.g. steric
hindrance or electro static interaction in the range of
few nm
surface activation by demicellization, diffusion and
adsorption of surfactants (ntm)
drop interaction, coalescence and surrounding fluid
dynamic in the range up to several mm
This broad range of length scales (106) makes modelling
and simulation within one numerical scheme impossible.
Instead, here the macro scale behaviour is described where
the relevant micro and meso scale processes (and the sub
grid scale) are exemplary incorporated.
micro meso macro
phisico surfactants, liquid entities,
chemical demicellization flow field,
interaction diffusion/ adoption coalescense
Figure 2: Phases in the computational model
Numerical modelling of the macro scale
A finite volume method representing the membrane
emulsification process determines the temporal and spatial
pressure, velocity and phase fraction during the process of
drop generation at the pores. The numerical model describes
the multiphase flow by means of NavierStokes (NS)
equations (momentum and mass conservation) plus an
interface capturing model that is able to control the interface
dynamics and the coalescence of liquid interfaces.
Based on a standard VOF implementation the CFDcode
Fluent 12 is employed [Ansys 2009]. The VOFmodel uses
one continuity and momentum equation with locally phase
fraction properties as:
Continuity:
( )+ V m )= 0 (4)
m m m
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Momentum:
(Pm m)+ V(Pm m vm)=
8t
S n \T
Vp V[tm (Vm VTm)] (5)
n
Pmg + F + V (akPkVdr,kVdr,k)
k=l
The continuity and momentum equations use phase averaged
values for density, velocity and viscosity. The values are for
phase fraction i and mixed at the interface region and
weighted by density p as:
Pm=Z i1Pi
For every phase i, a volume fraction coefficient a, exists. The
volume fraction for the disperse phases a, is defined as:
continuous phase:
disperse phase:
interface: 0 < a, < 1
a, =0
a, = 1
The volume fraction of the continuous phase is defined by:
1 a,1,csp (8)
so that:
Co1 = 1 (9)
With the standard VOFmodel formulation it is possible to
resolve the location and movement of the phase boundary of
the growing drop at the exit of the capillary. The surface
tension is concerned in a source term in the momentum
equation using the CSF (Continuum Surface Force) model
(Brackbill et al., 1992). The pressure drop across the surface
depends as given by Young and Laplace (1805) on the
surface tension coefficient, C, and the main surface curvature
radii in orthogonal directions, r1 and r2:
(10)
Ap=o+
Irl r*2
The curvature of the surface is determined from the
derivation of the phase fraction gradient, like described by
Brackbill et al (1992). The conservation equations are solved
using a pressure based implicit PISO pressurevelocity
coupling.
If a cell filled with one single phase (liquid 1 or liquid 2), the
fluid flow over the cell faces is calculated with the standard
method. Is an interface inside the cell, the georeconstruction
method is used. This method used a piecewise linear
assumption for the interface (Youngs 1982).
Numerical test case
Membrane emulsification in a channel with a diameter of 14
gm and pores with a diameter of 4 gm is modelled. The
Paper No
simulation is performed in a three dimensional domain
reflecting three distinct pores. The resulting domain is shown
in Fig. 3.
10 pm
dispersed phase
5 PJm

10 pmr
d = 4 p
d = 4 pm
d =14
TVl
continuous phase
Figure 3: Computational mesh with two pores
The following boundary conditions are used:
Channel wall: fixed contact angle, zero pressure
gradient, fixed velocity = 0
Channel inlet: constant velocity, fixed volume
fraction = 0
Channel outlet: constant pressure, zero velocity
gradient, fixed volume fraction = 0
A variable time step is used for the unsteady simulations. The
time step is set for CourantNumber Co < 0.5
with Co=VAt
Ax
Where At is the time step, v is the magnitude of the velocity
through that cell and Ax is the cell size in the direction of the
velocity. The solver chose for a fixed grid the time step
depending on the resulting mean velocity.
By grid independency tests it has been confirmed that for a
correct prediction of the surface bending and the contact
angle at the wall, a minimum of 10 cells across the pore
diameter is needed.
For the testcase the following constant material conditions
are used:
Phase 1: water: Viscosity gt: 0.001 kg /(m s), Density
p 998; kg/m3,
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Phase 2: oil : Viscosity g: 0.003 kg /(m s), Density
p 773; kg/m3,
Surface tension of water / oil: C 0.03 kg/s2,
contact angle of oil: 1800 (no wetting with oil)
Definition of algorithms for surface interaction
Coalescence of droplets in emulsions or at nearby pores is an
effect occurring at broad length and time scales. The
coalescence effect depends on droplet interactions in the gtm
range down to the molecular layer of surfactants. The
timescale includes the stability of emulsion over weeks and
month down to the spontaneous film rupture (without
stabilization).
The analysis of the droplet behaviour in emulsions, dispersed
phases or in sprays is of broad interest. With the VOFmodel
it is possible to resolve directly the droplets with a diameter
of a few gtm but not to resolve the contact film thickness of
10 to 100 nm during droplet impact/contact. Consequently
the physical aspects in this contact region are not suitable
captured during a conventional computation. Therefore, a
droplet/droplet contact (or in a bubbly flow a bubble/bubble
contact) leads to immediate coalescence or no coalescence
depending on the model. It is not possible to model for a
droplet/droplet impact the previous phase boundary
deformation and retarded coalescence.
Coalescence retardation model:
In this contribution a modification of the VOFmodel is
introduced, that allows modelling multiple droplet
interaction including contact, coalescence and/or rebound,
depending on the local flow conditions, the droplet contact
time and the contact area. For fundamental investigation first
the case of two single colliding droplets is studied (binary
collision). The potential of the model to describe multi
droplet processes is outlined where the model is extended to
multiple droplets and droplet processes. In this case separate
numerical phases qn+1 (not necessarily with identical
physical properties) for n droplet phases (separated entities)
and the bulk liquid phase are used. The separated numerical
phases are not allowed to merge as prescribed by separate
marker functions (e.g. 0 for surrounding phase, 1 for drop 1,
2 for drop 2 and n for drop n) used for both droplets and the
bulk fluid (Fig. 1).
For a binary contact, the gradients of the phase fractions aql
and aq2 are calculated. Contact of the phases in a cell is
assumed, when the product of the phase fraction gradient
magnitude is greater than zero:
Igrad (aq) *grad (aq2) > 0
(12)
When the contact is spread over a sufficient number of grid
cells, the contact time approach is initialized.
The contact film dimension is in the sub grid scale, it is not
resolved explicitly, but the drainage is described by
appropriate boundary conditions for normal and tangential
stresses (with the possibility to include locally varying
surface tension). In case the contact gets lost before
Paper No
coalescence occurs (bouncing or sliding) the contact time
marker is reset again. Otherwise, after the defined contact
time the phase fraction of phase 2 (drop 2) is switched to
disperse phase 1 (drop 1), resulting in an immediate
disruption of the droplets contact area and the droplet
volumes and interfaces will merge together. A single
interface is provided of the combined entities.
Gql (tcoalesce) = ql + aGq2; aq2, new => 0
(13)
Due to surface tension the new droplet try to energetically
optimize the surface energy and change the form to spherical
or spheroid when it is connected at the ends.
The functionality of the coalescence retardation model is
illustrated by numerical reproduction of the binary drop
collision experiment of Danner (2001).Here two droplets
pinned to the tips of two adjacent needles and brought into a
defined contact by horizontal movement of the droplets. In
the experiments the different contact scenarios have been
investigated. The numerical simulation of this procedure is
shown in Figure 4 and 5. Figure 4 shows the initial phase
fractions for the continuous phase, and the two disperses
phases. Here red means a Phase fraction of 1 and blue a phase
fraction of 0. Other colours denote intermediate phase
fractions and indicate an interface in the numerical cell. The
age of the interface between the two disperse phases is show
in Figure 5. Blue means no interface between disperse phases
and the other colours indicate the "age" in seconds.
Drop 1
\Drop 2
^^^^^^^^^
Figure 4: Phase fraction of the phases i (red =1, blue =0)
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
A result of a binary droplet interaction simulation during
emulsification using the new phase interaction model is
shown in Fig. 6 as a sequence of surface contours. The
retarded coalescence of two contacting droplets that are
moved towards each other in a surrounding liquid is shown.
This case has been studied experimentally in Danner (2001)
as illustrated in Fig. 7. The contact plane between the
droplets is shown in light (yellow) colour. After contact of
the droplets, the contact area grows without coalescence of
the droplets and the droplets deform. The droplets do not
coalesce immediately. The breakup of the contact film
initiating the retarded coalescence process occurs after t =
14ms. After breaking the contact film, the coalesced droplets
tend to form a common shape due to surface tension
(influenced by the holding pins of the droplets).
t=0s
t =15 ms
t = 7 ms
t =18 ms
t =14 ms
t =24 ms
Figure 6: Retarded coalescence of two merging oil droplets
in water
nach 5.7G a an$ 0 1A
Figure 7: Microscopic photos of coalescence sequence,
from Danner (2001)
Results and Discussion
Results for the behaviour of droplets at nearby pores in
membrane emulsification process are examined. In the first
illustrated sequence in Fig. 8 the standard VOFmodel is used
and the droplets coalesce immediately, when contact occurs.
A single larger droplet is formed, redispersed at the pores and
recoalesce at the pore exit.
Figure 5: "Age" of the interface between phase oil1 and
phase oil2 (marker function: Phase oil1 in contact with
phase oil2, accrues at the interface between phase oil1 and
phase oil2)
tins
S2.00e4
S1.8084
1.7 e14
1.6O i4
1.5 0c 4
1.40 _4
1.3 e84
1 .20 4
1_101S4
1.0 0._
9.02o16
7.00e55
0.01o09
S4.0 15
3.0 a*5
2D00o1S
I O.OO1
Paper No
@
V
*
t= 0 s t= 3 s t= 18 ts
Figure 8: Droplets modelled with standard VOFmodel,
instantaneous coalescence
The second sequence in Fig. 9 is modelled with the extended
VOFmodel for retarded coalescence. A fully stabilized
system is assumed (coalescence time oo ), so no
coalescence occurs in case of drop contact. In the third
sequence in Fig. 10, also modelled with the extended model,
no coalescence of the droplets occurs. The droplets get in a
contact at a small area (smaller than the critical contact area)
and get redispersed at the pores.
SR
es
t= 12 s t= 18 ts t= 21 s
Figure 9: Droplets modelled with extended VOFmodel for
retarded coalescence, stabilized droplets assumed, no
coalescence occurs
*
U
is
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
In the last sequence (Fig. 11) the same droplets get in contact
again. The contact area here is greater than the critical
contact area. Film rupture occurs after the defined contact
time and afterwards the droplets coalesce.
S.
t= 15 ts
e
t= 18 ts
t= 21 ts
Figure 11: Droplets at modelled with extended VOFmodel
for retarded coalescence, contact and coalescence
Conclusion
This paper presents a modification of the standard
VOFmodel that uses a separate phase for every single
droplet. A dropinteraction zone (boundary film) is
determined from phase fraction gradients. The fluid dynamic
conditions as well as geometric and temporal factors
determine if the film is stable or will rupture leading to
coalescence of the nearby droplets. This model allows
simulating droplet interaction including contact, coalescence
and rebound, depending on the local flow conditions, droplet
contact time and contact area. The presented results shows
behaviour of drop interaction as expected from experimental
analysis.
The contact time model need to be suitably extended from
film drainage models in relation to experimental
investigations and detailed microsimulations. The analysis
of complex dispersed multiphase problems with multiple
dispersed droplets (or bubbles) is possible in this approach.
Acknowledgements
This work was supported by the German Research
Foundation (DFG) within the Research Training Group 1375
"Nonmetallic Porous Structures for PhysicalChemical
Functions".
References
Bibette, J., Calderon, F., Poulin, P., Emulsions: basic
principles, Rep. Prog. Phys. 62 pp 9691033, 1999
Brackbill, J. U., Kothe, D. B. and Zemach, C., A Continuum
Method for Modeling Surface Tension. J. Comput. Phys.,
100, pp 335354, 1992.
t = 0 ts t= 6 = s t = 12 gts
Figure 10: Droplets modelled with extended VOFmodel for
retarded coalescence, short contact, sliding and redispersion
Paper No 7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Bremond, N., Thiam, A. and Bibette, J., Decompressing Press, 1982
Emulsion Droplets Favors Coalescence, Physical review
letters, Vol. 100, 024501, 2008
Burkhart, B. E., Gopalkrishnan, P V, Hudson, S. D.,
Jamieson, A. M., Rother, M. A., Davis, R. H.,Droplet Growth
by Coalescence in Binary Fluid Mixtures, Physical Review
Letters,Vol. 87, 098304, 2001
Chen, P.; Dudukovic, M.P; Sanyal,J, Threedimensional
simulation of bubble column flows with bubble coalescence
and breakup, AICHE J., Vol. 51, 3, pp 696712, 2005
Danner, T, Tropfen Kollision in Emulsionen, Dissertation,
Karlsruhe, 2001
Delnoij, E., Kuipers, J. A. M. & van Swaaij, W. P. M.
Numerical simulation of bubble coalescence using a volume
of Fluid (VOF) model, Third International Conference on
Multiphase Flow, ICMF 98, Lyon, France, June 8th 12th,
1998
Duineveld, P C.,: Bouncing and coalescence of bubble pairs
rising at high Reynolds number in pure water or aqueous
surfactant solutions; Appl. Sc. Res. 58, pp 409439, 1998
Ansys Inc: Fluent 12 user manual, www.ansys.com, 2009
Garabedian, R. S., Joseph Helblel, J., A Model for the
Viscous Coalescence of Amorphous Particles, Journal of
Colloid and Interface Science 234, pp 248260, (2001)
Hirt, C. W., Nichols, B. D.: Volume of Fluid (VOF) method
for the dynamics of free boundaries; J. Comp. Phys. 39, pp
201225, 1981
Koebe, M.: Numerische Simulation aufsteigender Blasen mit
und ohne Stoffaustausch mittels der Volume of Fluid (VOF)
Methode; Ph.D. thesis, Universitat Paderborn, 2004
Schramm, L. L.: Emulsions, Foams, and Suspensions:
Fundamentals and Applications WileyVCH, Weinheim,
2006
Schroder, V: Herstellen von 01 in Wasser Emulsionen mit
mikroporosen Membranen; Ph.D. thesis, Universitat
Fridericiana Karlsruhe, 1999
Schubert, H.; Armbruster, H.: Prinzipien der Herstellung und
Stabilitat von Emulsionen; Chem. Ing. Tech. 61, Vol 9, pp
701711 (1989)
Witelski, T P.; Bernoff, A. J.: Dynamics of
ThreeDimensional Thin Film rupture; Physica D 147, pp
155176, 2000
Young, T, An Essay on the Cohesion of Fluids, Phil. Trans.
Roy. Soc. 1, pp 6587 (1805).
Young, Thomas (doi: 10.1098/rstl. 1805.0005)
Youngs, D. L.: TimeDependent MultiMaterial Flow with
Large Fluid Distortion. In K. W. Morton and M. J. Baines
(editors), Numerical Methods for Fluid Dynamics. Academic
