Paper No 7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Numerical Simulation of a Compound Droplet by
ThreeFluid FrontTracking Method
Shunji Homma, Muneyuki Yokotsuka, Toshihiro Tanaka,
Jiro Koga and Gretar Tryggvasont
Division of Materials Science, Graduate School of Science and Engineering, Saitama University
255 ShimoOkubo, Saitama, 3388570, Japan
honma@al, s..ii.,in.,u .ik ip
Department of Mechanical Engineering, Worcester Polytechnic Institute
100 Institute Road, Worcester, MA 016092280, USA
gretar@wpi.edu
Keywords: three immiscible fluids, interface, fronttracking method, compound droplets
Abstract
We develop a threefluid fronttracking method in order to simulate the motion of a compound droplet that consists of three
immiscible fluids separated by two different interfaces. The two interfaces of the compound droplet are represented by two
different sets of the fronttracking elements immersed on the Eulerian grid mesh, where the velocities and the pressure are
calculated. The density and viscosity profiles with jumps at the interfaces are successfully determined from the location and
the connection information of the fronttracking elements. The motion of a compound droplet is simulated on axisymmetric
cylindrical coordinates. The results show that the threefluid front tracking method works appropriately for the motion of
compound droplets. Stability of a spherical compound droplet is examined numerically and the stable position is found
from the simulations.
Introduction
Although the flow involving three distinct fluids is rather
complex, it is often encountered in industrial processes.
Double emulsion is one example (Gart, 1997). This is
highly structured fluids of emulsion drops that contain
smaller droplets inside. The emulsion drop of this type is
usually referred as a compound droplet, which consists of
two immiscible phases in yet another immiscible
continuous phase. The fluid dynamics of the compound
droplets has been studied (Johnson and Sadhal, 1985), and
it has still attracted considerable attention, partly because
of its potential to manufacture highly structural materials
such as micro capsules for drug delivery vehicle (Prankerd
and Stella, 1990).
In this study, we develop a threefluid fronttracking
method based on the original fronttracking method
(Unverdi and Tryggvason, 1992) in order to simulate the
motion of a compound droplet that consists of three
immiscible fluids separated by two different interfaces.
The two interfaces of the compound droplet are
represented by two different sets of the fronttracking
elements immersed on the Eulerian grid mesh, where the
velocities and the pressure are calculated.
This paper is organized as follows: After describing how
the fronttracking method is modified, the density and
viscosity profiles are shown to verify if the modified
fronttracking method successfully determines the
densities and viscosities with jumps at the interfaces. The
motion of a compound droplet is discussed by comparing
with the theoretical and numerical results for a spherical
compound droplet in a uniform flow (Sadhal and Oguz,
1985).
Nomenclature
area
distance of the center between core and shell
droplets (m)
gravitational acceleration (ms2)
indicator function
pressure (Nm2)
radius (m)
radial coordinate
time (s)
translation velocity of shell droplet (ms1)
velocity vector
translation velocity of core droplet (ms1)
position vector
Paper No
z axial coordinate
Greek letters
S delta function
e eccentricity
K curvature (mn)
p viscosity (Pas)
p density (kg m3)
a interfacial tension coefficient (N m1)
Supersripts
T transpose
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
equations for onefluid formulation (Unverdi &Tryggvason,
1992) are the continuity equation, the NavierStokes
equation, and the equations of state for density and
viscosity:
Vu= 0,
pu + V puu
t D
D D
Dt Dr
VP + V u(Vu+ VuT)
+pg+ f orKnf(x xf)dAf,
f
Subsripts
1 external fluid
2 shell fluid
3 core fluid
f front
z 10R1
2 P
z
I 
g T
0 r r=R,
Figure 1: Computational domain for simulation of a
compound droplet.
Governing Equations and Numerical Method
Figure 1 shows a schematic of the problem. An initially
spherical compound droplet of radius R2 is placed in the
quiescent external fluid of viscosity ul and density pi.
The compound droplet consists of a shell (u2 and pz) and a
core (/3 and p3) whose radius is R3. Under constant
gravity, the entire droplet and the core move depending on
the density difference between the external fluid and the
droplet, or on that between the shell and core fluids. The
droplet eventually reaches the terminal velocity balancing
the buoyancy force with the viscous drag force. The
entire droplet and the core droplet can be deformed
depending on the relative contribution of the interfacial
tension to the other forces.
If the external fluid and the compound droplet are assumed
to be incompressible and Newtonian, the governing
Here, u is the velocity vector, P the pressure, p the density,
p the viscosity, a the interfacial tension coefficient, and t
the time. The interfacial tension is found using twice the
mean curvature, K, the unit normal, nf, and a delta function,
5(xxf), which is zero everywhere except at the interface,
Xf.
Equations (1) and (2) are discretized on the axisymmetric
cylindrical coordinate by finite difference approximations
with second order. The boundary conditions of the
leftside of the domain (zaxis) are symmetry, and the other
boundaries are noslipwall. The initial velocities are
zero (u 0) as the fluid is at rest. The motion of the
interface is traced by a fronttracking method (Unverdi and
Tryggvason, 1992). The details of the FrontTracking
method were described by Tryggvason et al. (2001), and
the axisymmetric version of the FrontTracking method
was examined in jet breakup problems by Homma et al.
(2006).
For the fronttracking method by Tryggvason et al. (2001),
the density and viscosity fields are determined by an
indicator function, defined, for example, by I = 1 for fluid
1 and I = 0 for fluid 2. Since the indicator function
depends on the location of interfaces, it is necessary to
determine every time step from the locations of the front
elements. The gradient of I has a value only at the
interface and can be expressed as
VI = AInf(x x)dAf .
f
Taking divergence of Eq. (4) results in the Poisson
equation
V2 =V.f n (x x)dA .
f
The right hand side can be calculated from the information
of the front elements and Al is unity for the above example.
The lefthand side is approximated by standard centered
differences, and solving the Poisson equation with the
appropriate boundary conditions yields the field of the
indicator function everywhere.
For threefluid front tracking method, we use an indicator
function defined by
Paper No
1 a fluid with maximum density (MXD)
I(x,t) = 0 a fluid with intermediate density (IMD) (6)
1 a fluid with minimum density (MND)
In this case, Al is 1 for the interface between MXD and
IMD and is 1 for that between IMD and MND. Once the
indicator function is determined, we can calculate the
density from
p(x,t) = max(pp,p,3) max(I(x,t),0)
+(p, + p2 + 3) max(p,p2,p,3)
min(p1,p2,p3) (1 sgn(I(x,t)) I(x,t))
min(p, p2,p3) min(I(x,t),0).
The viscosity can be determined by a similar manner to the
density.
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
lowReynoldsnumber flow. The solution is limited to
small capillary number ranges, where the interfacial
tension overcomes the viscous forces, thereby keeping the
compound droplet spherical during the translation. Since
the aim of this section is to compare our numerical results
with the theoretical ones, the condition similar to their
study was used in our simulations.
One important result of their study is stability of the
position of the core droplet inside the shell droplet. The
position depends on the density ratios between core and
shell droplets. Figure 3 shows the relative translation
velocity of the core droplet to the shell droplet (V) as a
function of eccentricity, which is the position relative to
the compound droplet itself (shell droplet). The
eccentricity is defined by e = d / (R2 R3), see Fig. 4.
The intersections between the curves and the line V/U = 0,
where U is the translation velocity of the compound
droplet itself or the velocity of the external flow, show the
equilibrium points. There are three kinds of equilibrium
points: stable, unstable, and metastable (e = 0). In this
study, whether the stable equilibrium is attainable, is tested
by our numerical simulations.
V/U
0 0.2 0.4 0.6 0.8
r/R1
0 0.2 0.4 0.6
r/R,
3 4
a2
0 0
0 0.2 0.4 0.6 0.8 0 0.2 0.4 0.6 0.8
r/R1 r/R1
Figure 2: Density (left) and viscosity (right) distributions
determined from the position of the fronttracking
elements.
O
s=1
0
Figure 3: Schematic of the core droplet velocity as a
function of its relative position (Sadhal and Oguz, 1985):
*stable equilibrium, X unstable equilibrium.
Density and viscosity field
In order to check whether the density and viscosity are
obtained correctly, those distributions are plotted for a
spherical compound droplet (R2 = 0.3 and R3 = 0.2) as
shown in Fig. 2. The ratios of density and viscosity in
this case are as follows: p2lp1 = 2, p3lp1 = 0.5, p/2//1 = 2
and p3//l1 = 3. The crosssectional distributions along the
center (a broken line in the top figures) of the compound
droplet show that the density and viscosity with jumps at
the interfaces are successfully determined (bottom figures).
Thus the algorithm based on Eqs. (5)(7) works
appropriately in the threefluid front tracking method.
Motion of a spherical compound droplet
Sadhal and Oguz (1985) studied theoretically the
translator motion of a compound droplet for
Figure 4: Schematic of the definition of eccentricity.
N 50
Results and Discussion
'oor 'o00'00000'0'0'1
Y
~
E~C = 0
configurations for /2/1 1=/u3//U1=l and
Figure 5 shows three eccentric configurations for RJR3 =
0.5 and A2/li = u3/ul = 1. The density ratio between
external fluid and shell droplet (p2/pi) is fixed at 1.2. The
configuration depends on the density ratio between
external fluid and core droplet (P3/P1). In the case of
P3/P1 = 1.01, the core droplet falls with shell droplet
because the average density of the compound droplet is
larger than the external fluid. In this case, e = 1 and it is
unstable equilibrium (the most left point in Fig. 3). For
P3/P1 = 0.89, the core droplet does not touch the top of the
shell droplet and the entire droplet falls holding this
configuration. This case shows the stable equilibrium,
where buoyancy forces balance with viscous forces. The
flow pattern and the dynamics of eccentricity in this case
will be shown. In the case of P3/Pl = 0.8, a similar
eccentric configuration to p3/lp = 0.89 is observed, but the
core droplet touches the top interface of the shell droplet.
This is also stable equilibrium configuration.
Figure 6 shows the dynamics of eccentricity and the
velocities of core and shell droplets for p3/lp = 0.89. The
velocity of the shell droplet increases slightly faster than
that of the core droplet, and those velocities eventually
equalize. Thus, V/U becomes 0 and the compound
droplet attains an equilibrium configuration. The
eccentricity increases and reaches around 0.8.
Figure 7 shows the comparison of the flow stream lines
around the compound droplet between exact solution
(Sadhal and Oguz, 1985) and fronttracking simulation.
The eccentricities of the exact solution and the
fronttracking simulation are 0.9 and 0.89, respectively.
There are double vortices inside the compound droplet.
Although the stream lines in the external fluid are slightly
different, the flow patterns of both the exact solution and
the fronttracking simulation are almost the same. The
threefluid fronttracking method developed in this study,
therefore, can reasonably simulate the translation of a
spherical compound droplet.
Paper No
p, / P =0.89
p, / =0.80
p,/p, =1.01
Figure 5: Eccentric
P2/pi = 1.2.
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
1.0 0.16
0.8
0.12
S 0.08 >
u, 0.4 
Eccentricity
0.2 *Shell Drop 0.04
ACore Drop
0.0 0.00
0 10 20 30 40 50
time
Figure 6: Eccentricity and velocities for stable equilibrium
configuration (p3/p = 0.89).
Figure 7: Comparison of the flow stream lines between
exact solution (left: e = 0.9, Sadhal and Oguz, 1985) and
fronttracking simulation (right: = 0.89).
Conclusions
Threefluid fronttracking method is developed in order to
simulate the motion of a compound droplet that consists of
three immiscible fluids separated by two different
interfaces. The two interfaces of the compound droplet
are represented by two different sets of the fronttracking
elements immersed on the Eulerian grid mesh, where the
velocities and the pressure are calculated. The density
and viscosity profiles with jumps at the interfaces are
successfully determined from the location and the
connection information of the fronttracking elements.
The motion of a compound droplet is simulated on
Paper No 7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
axisymmetric cylindrical coordinates. The results show
that the threefluid front tracking method works
appropriately for the translation of a compound droplet.
Stability of a spherical compound droplet is examined
numerically and the stable position is found from the
simulations.
References
Gart, N., Double emulsions scope, limitations and new
achievements. Colloids and Surfaces A: Physicochemical
and Engineering Aspects 123124, pp. 233246 (1997)
Homma, S., J. Koga, S. Matsumoto, M. Song and G.
Tryggvason, Breakup mode of an axisymmetric liquid jet
injected into another immiscible liquid, Chemical
Engineering Science, Vol. 61, pp. 39863996 (2006)
Johnson, R. E. and Sadhal, S. S., Fluid mechanics of
compound multiphase drops and bubbles. Annual Review
of Fluid Mechanics, Vol. 17, pp. 289320 (1985)
Prankerd, R. J. and Stella, V. J., The use of oilinwater
emulsions as a vehicle for parenteral drug administration.
Journal of Parenteral Science and Technology, Vol. 44, No.
3, pp. 139149 (1990)
Sadhal, S. and Oguz, H. N., Stokes flow past compound
multiphase drops: the case of completely engulfed
drops/bubbkes, Journal of Fluid Mechanics, pp. 511529
(1985)
Tryggvason, G., B. Bunner, A. Esmaeeli, D. Juric, N.
AlRawahi, W. Tauber, J. Han, S. Nas and Y. Jan, A
fronttracking method for the computations of multiphase
flow, Journal of Computational Physics, Vol. 169, pp.
708759 (2001)
Unverdi, S. 0. and Tryggvason, G., A fronttracking
method for viscous, incompressible, multifluid flows.
Journal of Computational Physics, Vol. 100, pp. 2537
(1992)
