4748
Langmuir 2004, 20, 4748-4755
Low Inertia Impact Dynamics for Nanodrops F. Gentner, R. Rioboo, J. P. Baland, and J. De Coninck* Center for Research in Molecular Modeling, Materia Nova/Universite´ de Mons-Hainaut, Avenue Copernic, 1, 7000 Mons, Belgium Received October 24, 2003. In Final Form: March 10, 2004 Simulations of a droplet impacting a flat solid surface with a small initial speed have been studied using molecular dynamics. Approximating the shape of the drop by a spheroid, spreading radii, and dynamic contact angles are measured. The data reproduce well experimental results from literature. We show that the difference between the equilibrium and the dynamic contact angle cosines, that is, the spontaneous driving force, versus the spreading velocity of the three-phase line varies with impact speed and consists of two distinct regimes which can be described by existing models of moving contact lines.
1. Introduction The dynamics of spontaneous drop spreading is a very active research subject. Within this topic, the problem is to find explicitly an efficient way to describe or even predict the base radius or the drop contact angle versus time. It is known that we can use a few good models for that such as the hydrodynamic model, the molecular kinetic model, the combined model, and so forth (see refs 1-3 and references therein) on the basis of some parameters. When a drop impacts the solid surface, it can easily be understood that the impact will enhance the spreading. It has indeed been possible to study such systems by experimental methods,4-7 numerical computations based on Navier-Stokes equations,8-10 or theoretical modeling.10,11 Recently,5,8,11,12 it has been shown experimentally that the wettability of the substrate can play an important role in the phenomenon. A few years ago, we successfully used large scale molecular dynamics (MD) to study the details of the spontaneous spreading dynamics. It is now our aim to extend these techniques to consider also the impact case. This technique has the great advantage of allowing us to change a single parameter individually, which is not always possible in real experiments. In this way, we hope to contribute to the understanding at the microscopic scale of the mechanisms controlling the impact dynamics. The cases studied here have relatively low inertia in the sense that the associated Reynolds and Weber numbers (1) Brochard-Wyart, F.; de Gennes, P. G. Adv. Colloid Interface Sci. 1992, 39, 1. (2) Gentner, F.; Ogonowski, G.; De Coninck, J. Langmuir 2003, 19, 3396. (3) Oshanin, G.; de Ruijter, M. J.; De Coninck, J. Langmuir 1999, 15, 2209. (4) Chandra, S.; Avedisian, C. Proc. R. Soc. London, Ser. A 1991, 432, 13. (5) Mao, T.; Kuhn, D. C. S.; Tran, H. AIChE J. 1997, 43, 2169. (6) Rioboo, R.; Marengo, M.; Tropea, C. Atomization Sprays 2001, 11, 155. (7) Rioboo, R.; Marengo, M.; Tropea, C. Exp. Fluids 2002, 33, 112. (8) Fukai, J.; Miyatake, O.; Shiiba, Y.; Yamamoto, T.; Poulikakos, T.; Megaridis, C.; Zhao, Z. Phys. Fluids 1995, 7, 236. (9) Geldorp, W. I. Numerical study of drop impact on solid surfaces. Msc. Thesis, Delft University of Technology, The Netherlands, and Darmstadt University of Technology, Germany, 1999. (10) Pasandideh-Fard, M.; Qiao, Y. M.; Chandra, S.; Mostaghimi, J. Phys. Fluids 1996, 8, 650. (11) Rioboo, R. Impact de gouttes sur surfaces solides et se`ches. Ph.D. Thesis, University Pierre et Marie Curie, Paris 6, France, February 2001. (12) Prunet-Foch, B.; Legay, F.; Vignes-Adler, M.; Delmotte, C. J. Colloid Interface Sci. 1998, 199, 151.
(Re ) FDV/µ and We ) FDV2/σ where D is the drop diameter, V is its impact velocity, and µ, σ, F are respectively the liquid viscosity, surface tension, and density) are very low for nanodrops. Thus, we will consider here the limiting case of an impact where inertia is supposed to play “little” role on the spreading. More general studies will be published elsewhere. The paper is organized as follows. Section 2 is devoted to the MD model used for these simulations. The results showing the spreading radius and dynamic contact angles follow in section 3 with the spontaneous case first and then with the influence of impact speed on the phenomenon. Interpretation and concluding remarks end the paper. 2. Simulation Model In this section, we present briefly for completeness the basis of the MD simulations and the definition of the model used to simulate the behavior of a droplet impacting a solid substrate. In a MD simulation, all atoms are considered as point masses defined by a position and its time derivatives. Each particle is interacting with the other particles through interaction forces derived from the interaction potentials detailed in this section, and time evolution is governed by Newtonian mechanics. At each time step, particle accelerations are calculated using Newton’s second law and their positions are then updated. This method is completely deterministic and has been widely developed in the literature.13-15 To detail the model, the main step is to define the force field containing all potentials of interaction between atoms. The basic interaction applied between all atoms, solid as well as liquid, is the pairwise Lennard-Jones (LJ) 12-6 interactions:
[( ) ( ) ]
Vij(r) ) 4ij
σij r
12
-
σij r
6
(1)
where r is the distance between any pair of atoms i and j. The parameters ij and σij are related in the usual manner to the depth of the potential well and the effective molecular diameter, respectively. Translated into reduced (13) Haile, J. M. Molecular Dynamics Simulation: Elementary methods; Wiley-Interscience: New York, 1997. (14) Allen, M. P.; Tildesley, D. J. Computer simulation of liquids; Oxford University Press: New York, 1987. (15) Israelachvili, J. Intermolecular and surface forces, 2nd ed.; Academic Press: New York, 1992.
10.1021/la030393q CCC: $27.50 © 2004 American Chemical Society Published on Web 04/27/2004
Impact Dynamics for Nanodrops
Langmuir, Vol. 20, No. 11, 2004 4749
(dimensionless) units,14 eq 1 becomes
Vij*(r*) ) 4
(
Cij
r*
12
-
Dij r*6
)
(2)
where the asterisk stands for reduced units. For simplicity, Cij and Dij are chosen the same for each type of atom. Following previous publications for such liquids and solids,16-19 we choose here Cff ) Dff )1.0, Css ) Dss ) 1.0, and Csf ) Dsf ) 0.6, where the subscripts stand for fluid/ fluid (ff), solid/solid (ss), and solid/fluid (sf) interactions. The intrafluid coefficients are standard and the solidsolid coefficients are chosen to produce a stable lattice structure at the temperature of interest. The chosen solid/ fluid interactions are characteristic of partial wetting.16 For computational performance, we apply a spherical cutoff at 2.5σij for the LJ interactions; that is, we only consider short-range interactions in these simulations. In addition to the basic interaction, specific potentials are used to describe the liquid or the solid. We simulate fluid chains into the drop by including a strong elastic bond between adjacent atoms within a molecule, of the form Vconf ) Dconf(r/σ)6 with r being the distance separating two consecutive atoms in the chain, σ being the LJ diameter, and Dconf ) 1.0. This intramolecular interaction forces the atoms within one molecule to stay together and reduces evaporation considerably. We apply a harmonic potential on the solid atoms so that they are strongly pinned on their reference face-centered cubid (fcc) lattice configuration, to give a realistic atomic representation of the solid surface. Our system is made of a liquid droplet composed by 1600 16-atom chains on top of a solid substrate composed by one layer of fcc unit cells (57 600 atoms). As the initial configuration, we consider a liquid drop far above the solid surface that we allow to equilibrate in a vacuum. This takes around 500 000 time steps. Before impact, the system is indeed equilibrated and looks like a sphere in a vacuum, far away from the solid surface. Once the two systems, liquid as well as solid, are equilibrated, the droplet is translated close to the solid at a distance of d ) 2.5σ small enough to initiate the contact, before impact. Then, we add an instantaneous downward initial velocity to the liquid atoms to mimic the impact of a drop issued from a spray. The cylindrical symmetry of the problem is used and the calculations are done in two spatial dimensions (r, z). During the following dynamics, we observe the positions of the atoms and molecules as a function of time. Once the droplet begins to spread, we only keep the temperature of the solid constant by rescaling the speed of the atoms, so as to allow the dissipation of energy between the solid and the liquid. This technique is indeed reminiscent of real experiments and has already led to interesting results within the spontaneous spreading case on a flat surface.17 We always apply a time step of 5 fs during our simulations with σ ) 3.5 Å. The values chosen for these parameters do not affect the spreading behavior, but they are usually used to compare the measured contact angles with true experimental results. (16) de Ruijter, M. J. A Microscopic Approach to Partial Wetting: Statics and Dynamics. Ph.D. Thesis, University of Mons-Hainaut, Belgium, September 1998. (17) de Ruijter, M. J.; De Coninck, J.; Blake, T. D.; Clarke, A.; Rankin, A. Langmuir 1997, 13, 7293. (18) de Ruijter, M. J.; Blake, T. D.; De Coninck, J. Langmuir 1999, 15, 7836. (19) De Coninck, J.; Dunlop, F.; Menu, F. Phys. Rev. E 1993, 47, 3.
Figure 1. Snapshots of the impact for Vimp ) 18 m/s. Times from left to right: t ) 0 ps, 0.05 ns, 0.15 ns, 0.3 ns, and 0.8 ns.
Figure 2. Side view of the droplet profile (open circles are positions where the density function is 0.5) during the spreading dynamics fitted by a circle (solid line). The inclination of the tangent to the circle defines the contact angle θ at the level of the solid layer, the atoms of which are represented by full circles at the bottom of the figure. The cylindrical symmetry is used.
Simulations have been launched for a range of impact velocities from 0 to 18 m/s. As a function of this impact velocity, we expect different behaviors of the liquid that we will detail in the following section. Figure 1 shows a typical series of snapshots at different times during the dynamics of the impact spreading case. The impact velocities studied in this article are qualified as low velocities (lower than 18 m/s). As we will see in the next section, these velocities are sufficient to induce an effect due to the impact, but we mean that inertia is low enough to allow the contact line to reach but not overshoot its equilibrium value. Considering that there is here no hysteresis here and that the equilibrium contact angle is of 105°, the equilibrium diameter is 1.16 times the initial drop diameter. 3. Analysis and Results Determination of the Shape of the Drop. To study the dynamical parameters and compute the associated contact angles, we must first establish the shape of the liquid/vapor interface in the vicinity of the solid. First, we subdivide the liquid droplet into several horizontal layers of an arbitrary thickness. The modification of this thickness does not influence the position or shape of the profile. The constraint on the number of layers is provided by the need to maximize the number of layers while ensuring that each layer contains enough molecules to give a uniform density. Then, each layer is analyzed and the density function is calculated in each small volume unit. We then locate the extremity of the layer as the distance where the density falls below a cutoff value of 0.5 times the liquid density. This method enables us to construct the complete profile of the drop (see Figure 2) and to determine how it evolves with time. Because surface tension effects dominate here, the resulting shape, whenever the impact speed is zero, has to be approximated by an arc of a circle. Only very close to the solid surface do we expect the profile to be perturbed by the solid for energetic and entropic reasons.19
4750
Langmuir, Vol. 20, No. 11, 2004
Gentner et al.
a2 volellips_cap ) πh2 2(3b - h) 3b
Figure 3. Scheme of the ellipsoid that approximates the shape of the drop during impact.
To avoid this problem, we investigated the profile as a function of the number of removed layers near the solid. To compare results of spreading radii or dynamic contact angles with those of macroscopic experiments, it is necessary to consider enough layers at the top of the drop but also to stay sufficiently far from the substrate. The circular fit using all the experimental points except the last four close to the substrate leads to stable results for drops with more than 20 000 atoms. The associated contact angle is then given by the tangent to the circle extrapolated to the solid surface as shown in Figure 2. Whenever the impact speed is different from zero, we may expect some hydrodynamic enhancement of the spreading dynamics. There is then no more reason to expect a spherical shape for the drop. The next easiest picture we may imagine that would correspond to real experiments is the spheroid (cylindrical ellipsoid). The flattened shape of the droplet shown in Figure 1 corroborates this hypothesis. This shape is indeed seen experimentally for low inertia impacts.7,20 We have, thus, studied in detail how our drop may be approximated by a spheroid. Using the cylindrical symmetry of our system, we may consider its two-dimensional projection into the X-Z plane. The associated profile should obey the following equation of an ellipse:
(x/a)2 + (z/b)2 ) 1
(3)
where a and b are the semi-axes respectively along the X and the Z axes. Here, x and z are the coordinates of the profile starting from the center of the ellipse. In Figure 3 is presented a scheme of the spheroid. Because the volume of the drop is approximately conserved during the impact process, we may use the formula of an ellipsoidal cap:
(4)
with h being the distance from the apex of the drop to the surface (see Figure 3). We can easily calculate the volume of the drop before impact when it has a spherical shape. The radius of the sphere is then 6.52 nm, and, thus, the associated volume 1161 nm3. We then fit the shape of the drop with an ellipse of the same volume V0 by fitting the parameters a, b, and h defined on Figure 3. Studying the regression coefficient R2 for the elliptic fit, we can see that the fit is good because R2 is always above 0.99 with about 60 points. From these fits, we measure the dynamic contact angle as the inclination of the line tangent to the ellipse at the point in contact with the solid surface. These procedures allow us to measure simultaneously the contact radius and the contact angle versus time for all the considered snapshots. Figure 4 shows the spreading radius and the dynamic contact angle versus time for the null impact case. In the same way, Figure 5 shows the dynamic contact angles and spreading radii for different low impact speeds. These results are reminiscent of real macroscopic experimental observations.5,7,8,11 The impact deforms the drop during the process, but it can be noticed that all curves in Figure 5 converge to the same equilibrium position. This implies a null or very low hysteresis. It can also be seen that increasing the impact speed results in an enhancement of the spreading in terms of decreasing the time to reach equilibrium. Flow Fields in the Drop During Spreading. Figure 6 shows the flux in the drop during its spreading for the null impact speed case and a low impact speed case (Vimp ) 18 m/s). The flow fields at times t ) 0.075 and t ) 0.125 ns show that the flow is accelerated from the contact line in both cases. Thus, for low impact velocity cases, the spreading is not due to impact pressure20 like in other studies6,21 on drop impact but to capillary effects. The case of spreading with an impact speed shows a velocity field with larger velocities and gradients in most of the drop volume. Comparing the two dynamics in terms of flow fields, we observe that the spreading with impact speed is faster than the one with null impact speed. In particular, at time t ) 0.275 ns for the null impact speed, the velocity field shows that the spreading is far from being completed. Moreover, strong gradients are still present in the drop. On the contrary, for the impact speed case, most of the spreading is already accomplished and the velocity field presents low velocities and low gradients compared with those of previous times.
Figure 4. Spreading radius (left) and dynamic contact angle (right) versus time.
Impact Dynamics for Nanodrops
Langmuir, Vol. 20, No. 11, 2004 4751
to slip freely across a solid surface. Adsorption of liquid molecules, immobilizing them to some extent, is unavoidable. The mobility of these liquid molecules in the vicinity of a solid has been described by Blake and Haynes23 using the MKT of liquids, developed by Eyring and coworkers.24-28 One introduces an equilibrium frequency of molecular displacement K, related to the molar activation free energy ∆G*:
K)
(
)
kBT -∆G* exp h NkBT
(5)
where kB and h are Boltzmann’s and Planck’s constants, N is Avogadro’s number, and T is the absolute temperature. In bulk liquid and at equilibrium, this frequency term should be isotropic. Blake and Haynes assumed that near the solid the frequencies are anisotropic during spreading, with K+ and K- being the frequencies of displacements parallel to the solid. The velocity of the contact line is then
v ) λ(K+ - K-)
(6)
with λ being the characteristic length of displacement. At equilibrium, the velocity of the contact line is zero, and on average, it can be written that
K+ - K- ) K0
Figure 5. Dynamic contact angles (top) and spreading radii (bottom) for four impact velocities [top to bottom (top) and bottom to top (bottom): 0, 8, 12, and 18 m/s].
We also observe on Figure 6 a bump at the times t ) 0.225 ns and t ) 0.275 ns. It reveals a small deformation of the drop due to the competition between impact and spreading. It corresponds to the moment when the vertical arrow’s length falls down and the impact effects become smaller, whereas the spreading process is dominant in the direction of the contact line. The region of this bump is even almost static at t ) 0.275 ns: this part of the drop is stopped and relaxes very slowly while a ribbon of flow shows that the spreading process is still acting in the direction of the contact line. The bump is retracting slowly and has completely disappeared at the time t ) 0.525 ns. 4. Interpretation Several models exist to predict and sketch the behavior of the drop in the case of spontaneous spreading.16,22 For a static contact angle higher than 90°, two models can be used to describe the spreading: the hydrodynamics theory if the speed is low1,25-27 and the Blake-Haynes model (molecular kinetic theory, MKT) at a higher velocity.2 When we consider an impact speed, the speed of the contact line is increased. Moreover, as we consider the first instants when the spreading speed is high, the MKT is more appropriate. Beside the large range of validity of this model in terms of angles and spreading velocities, it should be recalled here that this model is based on molecular displacements and, thus, is rather well-adapted to a nanoscopic study. In reality, a liquid is never allowed
(7)
with K0 being an equilibrium constant. If the contact line is in relative motion to the solid, the frequencies parallel to the solid can no longer be in balance (eq 6). If we assume for example that K+ is greater than K-, then this imbalance, whatever its source, facilitates the motion in that direction and makes it more difficult for molecules to move in the opposite direction:
K+ )
( (
) )
kBT ∆GW* w exp + h NkBT 2nkBT
kBT ∆GW* w exp K ) h NkBT 2nkBT -
(8)
with ∆GW* being the molar activation free energy of wetting, w the work per unit of area done by the driving force, and n the number of sites per unit area of solid at which this work is dissipated. Combining this with eq 6 results in
v ) 2K0λ sinh
(
)
w 2nkBT
(9)
In a final step, it was assumed that the driving force is the out-of-balance surface tension force acting on the (20) Schiaffino, S.; Sonin, A. A. Phys. Fluids 1997, 9, 3172. (21) Field, J. E.; Lesser, M. B.; Dear, J. P. Proc. R. Soc. London, Ser. A 1985, 225, 401. (22) Kistler, S. F. In Wettability; Berg, J. C., Ed.; Marcel Dekker: New York, 1993. (23) Blake, T. D.; Haynes, J. J. Colloid Interface Sci. 1969, 30, 421. (24) Blake, T. D. In Wettability; Berg, J. C., Ed.; Marcel Dekker: New York, 1993. (25) Hoffman, R. L. J. Colloid Interface Sci. 1975, 50, 228. (26) Tanner, L. J. Phys. D 1979, 12, 1473. (27) Voinov, O. V. J. Colloid Interface Sci. 1976, 11, 714. (28) Lesser, M. B.; Field, J. E. Annu. Rev. Fluid Mech. 1983, 97, 15.
4752
Langmuir, Vol. 20, No. 11, 2004
Gentner et al.
Figure 6. Velocity fields calculated over a 0.05-ns duration at the beginning, during, and at the end of the spreading process for a null impact speed (a, c, e, g, i) and a low impact speed Vimp ) 18 m/s (b, d, f, h, j) for representative steps [0.075 ns (a, b); 0.125 ns (c, d); 0.225 ns (e, f); 0.275 ns (g, h); and 0.525 ns (i, j)].
Impact Dynamics for Nanodrops
Langmuir, Vol. 20, No. 11, 2004 4753
Figure 7. Fitting curves using MKT applied on the spontaneous spreading case compared with simulation data. The solid line for the radius curve is calculated using the fitted parameters obtained from the contact angle data.
contact line. Per unit area swept out by the wetting line, this work can be written as
w ) γLV(cos θ0 - cos θd)
(10)
with γLV being the liquid-vapor surface tension and θ0 being the equilibrium contact angle. The dynamic contact angle used here is defined by thermodynamic considerations and is necessarily macroscopic. This results in the final equation
Vspr ) 2K0λ sinh
[
]
γLV(cos θ0 - cos θd) 2nkBT
(11)
Assuming the simplest case for flat surfaces n ) 1/λ2, we get the following fit for the zero impact speed case in Figure 7. Recall here that the static contact angle is high (above 90°) and, thus, the use of some linearized theory16 cannot be applied. From this fitting procedure, we can deduce the parameters K0 ) 50.0 × 107 ( 7.08 × 107 Hz and λ ) 3.5 ( 0.70 Å. We will see in the next section how the data with impact deviate from this model. Driving Force versus Spreading Velocity. The MKT cannot be applied directly to the impact cases because it does not take into account the role of inertia. Thus, an added force should be implemented in the model. Our starting point being the spontaneous spreading model based on the MKT,16 we, thus, need to express the problem in terms of the relevant variables: the driving force [cos θ0 - cos θd] versus the triple line speed Vspr(t). The first step is to compute the instantaneous velocity for the null impact speed directly with the analytical expression of the velocity versus the contact angle from eq 11. Figure 8 shows the case of the null impact speed in terms of the difference of cosines versus spreading velocity (we just recall here that the shape of the drop is approximated by a sphere). The first simulation data of the dynamics provide high spreading velocities (above 30 m/s) which confirm the necessity to choose the MKT. Let us now explain the way this instantaneous velocity is calculated for the impact cases where no explicit expression exists. Because the raw data results show some scatter due to the individual movement of the atoms (Brownian motion), it is necessary to “smooth” out the curves to focus and study the general behavior of the drop. The instantaneous spreading velocity is determined by fitting the radius represented in a log-log plot versus time, with a function, that we differentiate to get the triple line instantaneous speed. For all the impact velocities, the same shape of the curves indicates that the logarithm of the radius follows two distinct linear regimes with different slopes, except
Figure 8. Difference of cosines of dynamic and static contact angles versus spreading velocity for Vimp ) 0 m/s. The instantaneous velocity is calculated directly as an analytical function of the contact angle. The data follow the MKT model with K0 ) 50.0 × 107 ( 7.08 × 107 Hz and λ ) 3.5 ( 0.70 Å.
in the region where they cross each other and at the extreme end of the simulation, when we reach equilibrium. Hence, we fit both regimes of the logarithm of the radius linearly and do not take into account the estimated instantaneous velocities in the regions where the linear fits are of poor quality. As an example, Figure 9 shows for the impact velocity of 18 m/s the linear fits of the radius in the log-log plot. The quality of these linear fits with simulation raw data is very high (R2 > 0.98). For all the other impact velocities, the quality of the fit was always good in the sense that the smoothing function passes always in the near vicinity of the experimental points and the R2 were ranging between 0.98 and 1. Figure 10 shows the driving force versus the triple line speed for low impact speeds and again the case of the null impact speed as a comparison. When the impact speed increases, the overall shape of these curves [cos θ0 - cos θd(t)] versus Vspr(t) changes. The main difference between the impact cases appears for high contact line speeds, that is, in the beginning of spreading. At t ) 0, simple geometric arguments indicate that the “spreading” or contact line velocity is infinite as a result of the curvature of the drop.7 The advancing of the contact line is not due to the out-of-balance but to the movement of the drop toward the solid by its own impact velocity and subsequent inertia.7,21 Considering a very simple model of a sphere that is cut by a plane moving toward the sphere at Vimp, we can consider in two
4754
Langmuir, Vol. 20, No. 11, 2004
Gentner et al. Table 1. Values for the Fitted Ordinate and Slope from Equation 15 Compared with the Theoretical Values from the Geometric Model
Vimp
fitted ordinate
5 8 10 12 15 18
0.991 ( 0.016 1.080 ( 0.028 1.182 ( 0.024 1.098 ( 0.029 1.203 ( 0.019 1.171 ( 0.017
1/
theoretical ordinate 2 log(2RVimp) 0.907 1.009 1.058 1.097 1.146 1.185
fitted slope 0.503 ( 0.014 0.558 ( 0.021 0.607 ( 0.019 0.516 ( 0.021 0.571 ( 0.015 0.504 ( 0.013
theoreti- correlacal tion slope R2 0.5 0.5 0.5 0.5 0.5 0.5
0.989 0.986 0.991 0.982 0.994 0.993
drop radius before impact, x is the abscissa of the contact line, and y is its ordinate. We define ∆cos θ ≡ cos θ0 - cos θ d. At y ) 0, we get
x ) xVimpt(2R - Vimpt) ) R sinθd Figure 9. Droplet spreading radius versus time in a log-log plot for Vimp ) 18 m/s. The squares represent the data, and the solid lines respresent the linear fits.
(13)
that we derive to get the expression of the contact line speed Vspr:
Vspr(t) ≡ Vimp
dx ) dt
∆cos θ - cos θ0
x(1 - ∆cos θ + cos θ0) - (1 - ∆cos θ + cos θ0)2 R cos θd
Figure 10. Difference [cos θ0 - cos θd(t)] versus the triple line speed Vspr(t) for Vimp ) 0, 5, 8, 10, 12, 15, and 18 m/s respectively from left to right.
Figure 11. Description scheme of the geometrical model at the beginning of the impact process.
dimensions the intersection of the line with a circle (which is the contact line) and calculate its evolution. An explicit scheme of this geometrical model28 is shown in Figure 11. The contact line position is described by the equation of a circle:
x2 + (y - y0)2 ) R2
(12)
where y0 is defined as y0 ) R - Vimpt ) -R cos θd, R is the
)
dθd (14) dt
Let us remark that Vspr diverges at contact time Vspr ≈tf0+ (RVimp/2t)1/2. This behavior is characteristic of a drop that is only driven by inertia and falls vertically as if there is no solid surface. This corroborates the analysis of the flux in Figure 6, where the flow is mainly vertical in the two first figures at Vimp ) 18 m/s and the effect of the driving force, that is, large inclined vectors near the substrate, are visible because t ) 0.125 ns. Moreover, as the instantaneous velocity is derived from the radius and according to eq 13, we expect the following behavior for the radius at the beginning of the simulation: x ≈tf0+ RVimpt. This gives
1 1 log t + log(2RVimp) log x ≈ + 2 2 tf0
(15)
Hence, the slope in a log-log plot should be 0.5 and the ordinate should be 1/2 log(2RVimp). The log-log plot in Figure 9 as well as the linear smoothing of the radius used to estimate the instantaneous velocity is then justified. The parameters resulting from the fits of the first linear regime radius presented in the log-log plot in Figure 9 are presented in Table 1. In Table 1, the results corroborate the hypothesis of a linear regime for the radius when the time is small. The ordinate and slope data are remarkably close to those predicted by the geometric model. Inverting eq 14 gives an expression for ∆cos which is represented in Figure 12 for several impact speeds:
∆cos θ ≡ cos θ0 - cos θd ) cos θ0 +
Vspr2 Vspr2 + Vimp2
(16)
Of course this simple model does not take into account mass conservation, but the smaller the time after impact considered is, the smaller is the neglected mass and the more accurate is the analysis.
Impact Dynamics for Nanodrops
Figure 12. Difference of cosines versus contact line speed compared with the geometric model (valid for high contact line speeds) and the MKT (valid for medium and low contact line speeds).
Figure 12 shows the comparison at high contact line speeds between the simulation data for the difference of cosines at different impact velocities and the simple geometrical model. Let us note that the agreement improves with the contact line speed. However, as time increases, or in other words as the contact line speed decreases, the simulation data deviate from the geometrical model and tend to get closer to the spreading law corresponding to MKT. At low spreading velocities, the scatter of the data in Figure 13 represents a difference of the contact angle of about 10°, which is low considering the raw data. Another representation of the geometric model is shown in Figure 13 in terms of radius versus time. The agreement of this asymptotic model with the data is better at a higher impact speed. This model is valid at the beginning of the simulation (solid line) when the radius is low and the contact line speed is high, up to a threshold at about 5.5 nm and 0.3 ns. This threshold seems to depend only on the impact velocity, at least for these simulations. Conclusions We have, thus, developed MD simulations to study the drop impact dynamics with a small initial speed. During all our analyses, we have compared the impact case with the spontaneous spreading dynamics. We have shown that the shape of the drop can be fitted by a sphere when the impact speed is zero, that is, the spontaneous case. Whenever the initial speed is not zero, we have established that the drop shape can be well approximated by a spheroid, at least in our small impact speed case.
Langmuir, Vol. 20, No. 11, 2004 4755
Figure 13. Droplet spreading radius versus time compared with the geometric model prediction. The squares, circles, and crosses represent the data respectively for the velocities of 5, 10, and 18 m/s. The solid lines represent the geometric model in the range of small radii corresponding to high velocities, and the dashed lines are the continuity of these curves where this model is not valid anymore.
For all the considered impact speeds, we have, thus, measured the associated base radii and contact angles versus time. The corresponding flux lines inside the drop have also been derived. All our data show that the spreading velocity is a unique function of the driving capillary force γ(cos θ0 - cos θd), which differs only at high speed (just after impact) in agreement with the geometric model. Our simulations support, thus, the idea that the geometric model is able to describe the impact of the drop at a very high speed (just after impact) and the MKT model works well when the spreading velocity becomes small enough. Just after impact, the advancing of the contact line is, thus, due to the movement of the drop toward the solid by its own impact velocity and subsequent inertia. After some time, the friction between the liquid and the solid surface takes place and the driving force becomes the out-of-balance of the capillary force. The crossover time between the two mechanisms should of course be a function of the affinity between the solid surface and the liquid or, in other words, its wettability. This missing piece toward a complete theory is now under investigation. Of course, all these simulations have been considered in a vacuum but we expect the results to hold, at least qualitatively, when a fluid is present in place of the vacuum. Indeed in the presence of air, the drop would decelerate and slightly deform as it approached the wall leading to a small change in the crossover time between the two regimes. Acknowledgment. The authors acknowledge Daniel Weiss for stimulating discussions on the subject. This research has been partially supported by the Structural European Funds and by the Re´gion Wallonne. LA030393Q
