7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
A new Volume of Fluid model for handling wetting effects.
Application to the impact of emulsion O/W droplets on a moving plate
R. Guillaument*t, S. Vincent*, J.P. Caltagirone*
J. Duclost, M. Laugiert, P. Gardint
*TREFLE Laboratory (UMR CNRS 8508), ENSCPB, F33607 Pessac, France
t ArcelorMittal Research, Maizieres les Metz, 57283, France
guillaument@enscpb.fr
Keywords: Wetting effects, Threedimensional numerical simulation, Impact of droplets, Emulsion O/W, VOF methods
Abstract
The knowledge of oil thickness, difficult to measure, on the production line is an important industrial stake to control
of lubrication in cold rolling system. Lubrication is accomplished by an oilinwater emulsion spray, so it is necessary
to model the behavior of emulsion droplet impact on a moving plate.
The hydrodynamic phenomena are simulated in solving two phase Navier Stokes equations in threedimensional
Cartesian coordinates with a Finite Volume method.
To describe capillary effects, a model of contact angle and a three phase contact line model are developed in the
computational fluid dynamic (CFD) code.
The contact line model is validated according to the theory of the floating lens on a bath of fluid. The contact angle
model allows wetting effects to be simulated with a new VOF method (3), that was validated by the impact of water
droplets on normal (14) and on oblique substrates (8).
Thanks to the coupling of models developed in this paper, the impact of an oilinwater emulsion droplet on a wetting
plate is studied as a function of two parameters: the oil concentration and the impact velocity. No experiments are
found in the literature to validate the evolution of oil thickness.
Nomenclature
Roman symbols
g gravitational constant (m.s 2)
p pressure (N.m 2)
t time (s)
u local velocity (m.s 1)
C liquid volume fraction ()
do droplet initial diameter (m)
Fv surface tension force (N.m 3)
ni the normal to the interface( )
Greek symbols
p dynamic viscosity (Pa.s)
K local curvature of the interface (m 1)
6i a Dirac function indicating the interface (m 1)
a surface tension coefficient (N. m 1)
p density (kg.m3)
Superscipts
g Gas
d Droplet
1 Liquid
s Solid
Introduction
Impact and spreading of multiphase droplet on solid
wall is of fundamental importance in many engineer
ing and natural processes including jet printing, spray
coating and spray cooling. The three phase moving con
tact line and wettability properties are a notoriously diffi
cult problem involving highly complicated physical pro
cesses and offers a challenge for computational mod
els. Furthermore, velocity measurement by Particle Im
age Velocimetry (PIV) is complex, due to the interfaces
which induce optical difficulties. Although numerous
models and solutions to this problem have been pro
posed. A successful numerical modeling of the impact
and spreading of an emulsion droplet on partially wet
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Ki gpj
Figure 1: A floating lens, definition of total thickness
(e) and submerged thickness (e')
ting substrate should address : (i) resolve the interface
between the several fluid, (ii) accounting for the effects
of interfacial tension, (iii) resolving the contact line, (iv)
incorporating the effect of the solid wettability. The pa
per is organized as follows: the foundation is described
in the next section, the modeling is explained in section
2, the validation and results are presented and discussed
in section 3 and some conclusions are drawn in section
4.
Foundation
The contact angle plays a part capillary effects on wall.
The deposition of droplet on a wetting surface or non
wetting, when it impacts, are controlled by the contact
angle ec defined by the Young law (5):
coseO,
Ysg Ysl
Where y,, is a surface coefficient between solid and gas,
sy, is a surface coefficient between solid and liquid,;yl
is a surface coefficient between liquid and gas. This
law leads the surface energy of droplet where the sur
face area between every phase and the total area on the
wall are summed, for small droplets, gravity effects are
neglected so the droplets are spherical. The contact an
gle is a local property that do not depend on droplets
diameter but it depends on state of surface. The contact
angle guides the velocities of contact line where a three
phase flow exists. This contact line is studied for the
floating lens case on liquid both in figure 1. The three
fluids are non miscibles where the fluid A and B are liq
uids and the fluid C is a gas. The thickness (e) is the
total thickness of the lens and the thickness (e') is the
submerged thickness of the lens. In the De Gennes the
ory (5), PA, PB and pc are the density of fluid A, B
and C respectively where pc is neglected. In function
of three surface tensions (YA, YB, TAB), the spreading
parameter S of A on B is determined by the following
equation: S = s (YAB +YA). The spreading is con
sidered total if S>0 when the fluid A spreads completely
on fluid B, else the spreading is partial for S<0 so the
fluid A floats as a lens on fluid B. The curvature radius
is defined by the following equation : where nK is the
curvature of phase i:
At equilibrium, the total force existed on the lens
called F(e) is defined by
1
F(c) [ pAg(c
2
1 2
S')2 + 1 C/g(pB
2
PA) S] (3)
Where V is total volume of lens and v is a surface of
the lens. Hence, by Archimedes law, pAge page'
with assumption that the density of fluid A is compared
to the density of fluid B. Consequently, e' is negligible
in front of e. So the equilibrium force is simplified and
defined by
1
F(c) [ Pm9 2
2
with the average density pm, (PB PA). The
assumption of conservation of fluid volume is developed
by De Gennes (5), so the asymptotic thickness (e) is
defined by
2S8
ef P=(m
(Omf
Physical Model
Conservation equations
The flow is assumed to obey the NavierStokes
equations for incompressible viscous fluids in three
dimensional Cartesian coordinate system maintained
at isothermal condition with constant thermofluid
properties. The roughness of the steel strip is not
explicitly taken into account. The model is based on
the convolution of the motion equations in each phase.
In the present model, the effects of gravity, viscosity,
surface tension at free surface, and wettability between
liquids and solid are taken into account. The mass and
momentum conservation equations, are given by
at
at u+uVu)
D( u
VP + pg + Au + Fv(6)
V.u = 0 (7)
(8)
Where t, u = (u,v,w) and p are time, velocity and pres
sure. Moreover p ,p and Fv denote apparent fluid
density, apparent viscosity and surface tension volume
B
force, respectively. The simulation and the modeling of
multiphase flow, where the characteristic scale of inter
face length are bigger than the space step of grid, are car
ried out by one fluid model (10). A color function (C)
that defines the interface evolution on time is updated by
an advection equation. A the interface, it is assumed that
is no sliding between the fluids (6). The velocity field is
continuous through the interface. For every fluid k, the
color function Ck is used as a heaviside function :
Ck t) 1, for r e phase,
C'(, ) 0, else .
The interface between phase i and other phases is de
fined by isosurface C(=0.5. The local quantities such
as density or viscosity are built with the color function
according to the following numerical mixing laws.
P = ko1CkPk
,t E=k0 Ckk
Hence, the three phase flow is modeled by an equiv
alent fluid and the NavierStokes equations in the one
fluid model (10), (12) are completed by material inter
face tracking equation.
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
the interface C = 0.5. The coefficient Ah is equal to
the local characteristic size of the grid cell (3). On a
numerical point of view, equation (11) is implicitly dis
cretized with finite volume and a centered scheme, to be
consistent with the NavierStokes equation approxima
tion. The linear system is solved with a direct MUMPS
L;iAh2
solver (7) If we define T* as  the numerical
algorithm for obtaining is (3):
algorithm for obtaining CS is (3):
CS'1 C
For k = 1..N 1, solve
V T*VC'Sk + CSk+1 CS
The Continuum Surface Force CSF developed by Brack
bill (2), consisting in summing surface tension forces at
an interface on the outline of the control volume, is used
to model Fv The capillary force reads by
Fv(M) = y(M)ni(M)6i
Where is a surface tension coefficient.
(n) and the Dirac function (68) indicating
calculated from the gradient of C,, as nS,
the curvature K of interface is obtained by
The normal
interface are
 VC, and
Vr Vc
K v I
S+ u VCk
at
In the evolution of paper, the different models of surface
tension effects, devoted in the closing of the one fluid
model, are presented. The methods used to resolve con
servation equations are based on finite volume approach,
structured grid and PLIC Volume Of Fluid (VOF) ap
proach (15) which are presented and are validated in the
paper Vincent (et.al) (13) and Youngs (15). In general,
VOF methods used are discontinu on the one grid space.
A the consequence, the curvature calcul requested to sur
face tension simulation is difficult numerically.
So, a new function is build, corresponding to a smooth
function of VOF called Cs (Smooth VOF, SVOF). This
function is used only for the curvature calculation (3)
and (11). The color function C is not replaced by Cs in
order to keep properties of mass conservation by algo
rithm PLIC. By analogy with the previous concentration
developments, the SVOF method consists in building a
smooth VOF function Cs by solving an Helmhotz equa
tion such as (??) with a source term initially equal to the
sharp VOF function C as follows:
V DVCS + CSn"+1 s'" (11)
where the diffusion coefficient D is equal to Li Ah2.
This parameter is fixed in order to ensure that the VOF
function Cs spreads on a distance Li on each side of
Wetted model
With the SVOF method, a color diffusion is applied to
the surface in the first cell. A correlation is established
between the contact angle and the value of the smooth
VOF function coefficient called wetted coefficient (Cp).
Thanks to control of diffusion, different properties are
modeled by the equation following :
V DCS + CS + Bik(Cs Cp) = C (17)
With Bik is penalization coefficient that strives toward
infinity on the wall and strives toward zero in the fluid
domain. In figure 2, the color function C is presented by
line and the Smooth VOF Cs to presented by a flood
color. It is observed that penalizing Cs on the wall
(flood color) allows to impose a contact angle and that
Cs is smoother than C, Cs =0.5 being almost the same
as C 0.5, the laws between contact angle and the wet
ting coefficient Cp are defined in figure 3 As previ
ously mentioned, the surface properties are not modeled
by the color function so the contact angle is not defined
whereas with a new model SVOF, the penalization on
the wall allows to keep the continuity between surface
and droplet with known contact angle.
Figure 2: The comparison of color function C defined
by line and wetted coefficient Cp defined by
the flood color
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
v (v cC 1 c.
Line triple
PhaseS
Figure 4: definition of normal for contact line model
Figure 5: The comparison of color function C defined
by line and Smooth VOF Cs defined by the
flood color in the case of the floating lens
Validation and Results
1,(Cp)= 75 Cp 170
SConstantS00250 EXPO256.128
* EXPos00 20 *Constant256128
8 (Cp) 286Cp+286
wetting coefficient (Cp)
Figure 3: Evolution of the numerically measured con
tact angle 0c according to imposed values of
Cp on a flat surface.
Contact line model
The volume force Fv, appearing in Eq 6 thought
(15 16), is generalizing by the model of Cranga in fig
ure 4 (4) who postulates that at every point, Fv is the
sum of the contribution of the capillary forces for each
couple of fluids :
Fv
Si=
3,j=3
.. i. VCtj
1Ej 1
where
C
Cij 
(3 +t cj
This model is improved by using C, instead of Cj,
which illustrated in figure 5.
Floating lens
The case of validation is the experience performed
in the thesis of Bonometti (1), with the following
characteristics:
surface tension between fluids A and C : 0.148N.m
surface tension between fluids B and C : 0.148N.m
surface tension between fluids A and B : 0.296N.m
initial diameter of drop(d) : 0.04m
step time (dt): 0.5e 3s
gravity (g): 9.81m.s 2
We calculate the ratio between e and e ', according
to the law of Archimedes, the density for fluid A is
10(nil.,,,, 3 and that for fluid B is r5lllli. i ,, 3.
We obtain a ratio of 0.2 in the theory and 0.18 in the
simulation for a 90x60x90 mesh. In the figure 6, the
result of simulation is equal to asymptotic solution
with the SVOF method, in figure shows the good
resolving of curvature is necessary to obtain to the equi
librium between the pressure and the capillary pressure.
Droplet impact on the wall
The study of a water drop impacting a hydrophilic
surface (v, lliiii I or a hydrophobic surface (no v, lliI i"
is characterized by simulating the conditions of Wang
experiments (14). The droplet diameter is 2.15mm and
the impact velocity is 0.514m/s for a Reynolds number
of 1121 and a Weber number of 8.
The definition of diameter and height measurement
are defined in figure 7 for the horizontal and oblique
cases.
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
a) NORMAL COLLISION
[E]=e/d
o,9 *
o
0,6 *
05* 
,7 
o,. +
0,4
.s 4 simulation
0,  DeGenne [T]=t(d/g)^O.
S2 4 10 12 4 16 18 20
d
h
d
Figure 6: case of test of the floating lens simulation
(blue symbol), red line asymptotical theory of
De Gennes
S U Ih
d
b) OBLIQUE COLLISION
Threedimensional simulations are leaded in a calcu
lation domain whose dimensions are 6dp x 1.5dc x 6dp
where dc represents the initial diameter of the droplet.
The grid is exponential with 120 x 80 x 120 calculation
points in each directions. The smaller space step in each
direction is 20.10 ,,.. The time step is 4.106s.
The first case concerns the study of the droplet impact
on a wetting surface. The experimentally measured con
tact angle Os is in the range 20 40 According to the
SVOF numerical law, a value of Cp = 0.9 is imposed on
CS on the flat surface, such that Oc = 28.6 The simu
lated results are proposed in figure 8. The time evolution
of the dimensionless droplet height during time is well
predicted by the SVOF approach. Concerning the di
mensionless diameter, differences are observed between
experiments and simulations, certainly due to the defini
tion of the diameter and the experimental technique used
to measure it which are not described in the publication
of Wang (14). Numerically, the maximum of the droplet
diameter is undertaken.
A second test case involving a nonwetting surface is
simulated to verify the ability of the SVOF method to
correctly predict contact angles. Following Wang exper
iments (14), the contact angle Oc varies between 80
and 130 Numerically, Cp = 0.55 is imposed on the
impact surface such that a contact angle of 128 is ver
ified in the simulations. In figure 9, a numerical and
experimental droplet dimensionless diameter and height
are presented according to time. As for the wetting case,
a very good agreement is demonstrated between the sim
ulations and the measurements of Wang (14). In the hy
drophobic case, the droplet spreads on the flat surfaces
until t 4.103s while a receding motion is observed
for t < 4.103s corresponding to a droplet oscillation
during impact.
The third test case, impact of a droplet on an oblique
surface is now undertaken. The orientation of the
surface is accounted for by changing the orientation
of the gravity. A water droplet of initial diameter
Figure 7: Definition of droplet diameter, droplet height
and contact angle for a droplet impact on a flat
surface
dp 0.56.10 3m is considered. The gravity effects
are negligible in this case as dp is lower than the
capillary length equal to 2.7.103m. The size of the
simulation domain is 6dp x 2dp x 1.5dc while the
calculation grid is exponential with 200 x 100 x 200 in
each Cartesian direction. The size of the smaller space
step is 12.1, ,. The time step is chosen such that
the CourantFriedrichsLevy CFL conditions is less or
equal to 0.2. The experiments of Fujimoto (8) are used
to validate the simulations. The results are presented
in figures 10. Again, a very good agreement is found
between simulation and experiments.
Emulsion Droplet impact on the wall
A spherical emulsion droplet impinges a horizontal
solid substrate at a temperature of 20 C (figure 11). Wet
ting coefficients are applied at the same time for water
and oil, the contact angle being respectively 150 and
20 The respective densities of air, water and oil oil are
1.2 Kg.m3, 1000 Kg.m3 and 900 Kg.mM3. The di
ameters of water and droplets are 100 pm and 10 pm.
The thickness of oil as a function of concentration, ve
locity of spray and surface tension are studied. The nu
merical domain is 400pm x 400pm x 400pm with a
space scale of 2 pm. The boundary conditions are slid
ing on the steel strip, periodicity on left and right slides
and symmetry for other directions. The definition of a
percentage of set oil, which is the oil volume occupying
the ten first micrometers over the steel strip, is consid
ered in comparison with the total oil volume in the cal
culation domain in figure 12 et 13.
The results of figure 12 show the percentage of set oil
when the surface is hydrophilic without precursor film
on its surface initially. The first impact of droplet is re
alized at 25 pus and the second impact at 75 Ips. The
behavior of emulsion is similar for concentration from
0.2 % at 2 %, in the first time, the droplet impact al
lows augmenting concentration in oil. Next time, the
crown around the zone of impact is created causing a hy
drodynamic flow against the direction of impact, draw
ing away the oil droplet. To finish the second and third
droplet impacts have little effect upon deposition of oil.
The study investigates the impact velocities of emul
sion drop to estimate the behaviour of oil droplet on the
steel strip. The results are showed in figure 13. For
all velocities, almost 10% of the oil of emulsion impact
on the plate. After impact, the water droplet is sheared
transporting the oil. The larger velocity is low and the
larger shearing is weak so the migration of oil is slower.
The Reynolds number and the Weber number are vari
able from 200 at 1000 and from 117 at 2935, respec
tively.
Conclusions
Development of a smooth VOF function allows to sim
ulate curvature correctly so capillary effects are deter
mined precisely. The coupling of the model of contact
angle and contact line allows to model impacts of emul
sion droplets on a wall with different properties of wet
tability. The validation of model are conclusive on dy
namic and static contact angle cases.
Acknowledgements
The authors thank the Aquitaine Regional Council for
the financial support dedicated to a 256processor cluster
investment, located in the TREFLE laboratory. We are
grateful for access to the computational facilities of the
French CINES (National computing center for higher
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
education) and CCRT (Center of Research Calcul and
Technology) under project number x2009026115 and
gen6115 rescpectly.
References
[1] T. Bonometti Developpement d'une m6thode de
simulation d'6coulements a bulles et a gouttes.
PhD thesis, Institut National Polytechnique de
Toulouse, 2005.
[2] J.U. Brackbill, B.D. Kothe, and C. Zemach. A con
tinuum method for modelling surface tension. J
Comput. Phys., 100(335):335354, 1992.
[3] J.P. Caltagirone, S. Vincent, R. Guillaument. A
new numerical method for handling wetting effects
in VOF approaches. submitted to Computer and
Fluids, 2010.
[4] J. Cranga Simulation num6rique direct des 6coule
ments ditri phasiques engendr6s par l'injection de
bulles dans un bain de m6tal liquid. PhD thesis,
Institute National Polytechnique de Toulouse, 2002.
[5] P.G. DeGennes F B. W. and Qu6r6 D. Gouttes,
bulles, perles, et ondes. Belin, 2005.
[6] J.M. Delhaye. Jump conditions and entropy
sources in twophase systems. local instant formu
lation. Int. J Multiphase Flow, 1:395409, 1974.
[7] I.S. Duff, P.R. Amestoy, and J.Y. L'Excellent.
Multifrontal parallel distributed symmetric and
unsymmetric solvers. Comput. Methods in Appl.
Mech. Eng., 184:501520, 2000.
[8] H. Fujimoto, Y. \ lII.II, and A.Y. Tong. Three
dimensional numerical analysis of the deformation
behavior of droplets impinging onto a solid sub
strate. International Journal of Multiphase Flow,
33:317332, 2007.
[9] I. Gustafsson. On first and second order symmet
ricfactorization methods for the solution of elliptic
difference equations. Phd dissertation, Chalmers
University of Technology, 1978.
[10] I. Kataoka. Local instant formulation of twophase
flow. Int. J Multiphase Flows, 12(5):745758,
1986.
[11] G. Pianet, S. Vincent, J. Leboi, J.P. Caltagirone
and M. Anderhuber. Simulating compressible gas
bubbles with a smooth volume tracking 1Fluid
method. Int. J Multiphase Flows, 36:273283,
2010.
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
[12] R. Scardovelli and S. Zaleski. Direct numeri
cal simulation of freesurface and interfacial flow.
Ann. Review Fluid Mech., 31:567603, 1999.
[13] S. Vincent and J.P. Caltagirone. A one cell lo
cal multigrid method for solving unsteady incom 2 .
pressible multiphase flows. J. Comput. Phys.,  *
163:172215, 2000. .
[14] M.J. Wang, F.H. Lin, J.Y. Ong, and S.Y; Lin. Dy  1, s
namic behaviors of droplet impact and spreading ; 0
water on glass and paraffin. Colloids and sur
face A : Physicochemical and Engineering aspect,
339:224231, 2009. _
[15] D.L. Youngs, K.W. Morton, and M.J. Baines. tlme.v./d,
Timedependant multimaterial flow with large fluid
distorsion. In Numerical Method for Fluid Dynam
ics, Academic Press, New York, 1982.
between 0.4 ms and 3.2 ms
between 3.6 ms and 6.4 ms
between 6.8 ms and 9.6 ms
between 10 ms and 11.6 ms
Figure 8: Relaxation of droplet wetting diameter(d),
droplet height(h), and contact angle (c0) for
water droplet impinging onto glass (wetting
surface
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Sdexpermental dmaxCFD
hmaxCFD hexpenmen. l
n, 1 t5 2
tme. tv/d,
between 0.4 ms and 3.2 ms

between 3.6 ms and 6.4 ms
between 6.8 ms and 9.6 ms
al a a a
ae ae a a
af
between 10 ms and 11.6 ms
Figure 9: Relaxation of droplet wetting diameter(d),
droplet height(h), and contact angle (c0) for
water droplet impinging onto paraffin (no
wetting surface)
Shcf dd h[fulmoto071
4, .....T ........
2 't .. t ..............
S02,5
T o,o
Time,tVo/dp
Figure 10: Time evolution of spreading
length and height of liquid for
(d vo,6 0)=i, .1.,.,. 2.8m/s, 45)
Figure 11: Impact of emulsion droplet on oil wetting
surface, the upper part shows the tree phases
where the oil is in yellow color, the water
is in blue color and air is transparent. The
lower part shows only the phase oil.
5% oil
1% oil
S 2% oil
0.4% oil .: ./ o.2% il

Figure 12: The percentage of set oil volume on a dry
hydrophilic surface for several concentration
of oil in water
2,5 1 3,5
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Figure 13: The percentage of set oil volume on a dry
hydrophilic surface for several velocities
