Langmuir 2000, 16, 2363-2368
2363
Experimental Evidence of Several Time Scales in Drop Spreading Michel J. de Ruijter, Magali Charlot, Michel Voue´, and Joe¨l De Coninck* Universite´ de Mons-Hainaut, Centre de Recherche en Mode´ lisation Mole´ culaire, Place du Parc 20, B-7000 Mons, Belgium Received June 15, 1999. In Final Form: October 21, 1999 We studied the dynamics of spontaneously spreading drops of viscous liquids over several decades of time. We show that there exist at least two different time scale behaviors which are compatible with the recently proposed combined theory of spreading dynamics.
Introduction The dynamics of wetting, or the displacement of one fluid in contact with a solid by another fluid, is central in many industrial applications, such as coating deposition,1 printing,2 and enhanced oil recovery.3 It has frequently been characterized by the spontaneous spreading of droplets.4-17 Within this example of spontaneous wetting, the variation in time of the contact angle θ, base radius R, and/or contact line velocities ν is in general measured by means of optical inspection. The advantages of this method are that only small samples of flat solid substrates are required, that a large range of contact angles is covered by a single experiment, and that the geometry of the drops is well-known. One disadvantage exists in the fact that with spontaneous spreading, the dynamics are governed by the physics of the system, not by the user. Several models have been proposed in the literature to describe the dynamics of the associated contact angle or base radius. These theories differ essentially by the way they describe the dissipation of energy during the drop shape transformation. Among them, we have the friction dissipation theory due to Blake et al.18,19 also known as the molecular kinetic model and the hydrodynamic dissipation theory due to Voinov7 and many others. For vanishing contact angles, it has been shown20 that the hydrodynamic regime may be characterized (1) Blake, T. D.; Clarke, A.; Ruschak, K. J. AIChE J. 1994, 40, 229. (2) Karnik, A. R. J. Photogr. Sci. 1977, 25, 197. (3) Paterson, A.; Robin, M.; Fermigier, M.; Jenffer, P.; Hulin, J. P. JPSE 1998, 20, 133. (4) Strella, S. J. Appl. Phys. 1970, 41, 4242. (5) Dussan V., E. B.; Davis, S. H. J. Fluid Mech. 1974, 65, 71. (6) Ogarev, V. A.; Timonina, T. N.; Arslanov, V. V.; Trapeznikov, A. A. J. Adhes. 1974, 6, 337. (7) Voinov, O. V. Fluid Dyn. 1976, 11, 714. (8) Tanner, L. H. J. Phys. D 1979, 12, 1473. (9) Chen, Q.; Rame´, E.; Garoff, S. J. Fluid Mech. 1997, 337, 49. (10) Dodge, F. T. J. Colloid Interface Sci. 1988, 121, 154. (11) Extrand, C. W. J. Colloid Interface Sci. 1993, 157, 72. (12) Zosel, A. Colloid Polym. Sci. 1993, 271, 680. (13) Brochard-Wyart, F.; de Gennes, P. G. Langmuir 1994, 10, 2440. (14) Seaver, A. S.; Berg, J. C. J. Appl. Polym. Sci. 1994, 52, 431. (15) Cazabat, A. M.; Gerdes, S.; Valignat, M. P.; Vilette, S. Interface Sci. 1997, 5, 129. (16) de Ruijter, M. J.; De Coninck, J.; Blake, T. D.; Clarke, A.; Rankin, A. Langmuir 1997, 13, 7293. (17) de Ruijter, M. J.; De Coninck, J.;Oshanin, G. Langmuir 1999, 15, 2209-2216. (18) Blake, T. D. The Contact Angle and Two-phase Flow. Ph.D. Thesis, University of Bristol, 1968. (19) Blake, T. D.; Haynes, J. M. J. Colloid Interface Sci. 1969, 30, 421. (20) de Gennes, P. G. Rev. Mod. Phys. 1985, 57, 827.
asymptotically by:
R ∼ t1/10
(1)
θ ∼ t-3/10
(2)
and
However, within the same conditions, the friction regime leads to different time scales, i.e.
R ∼ t1/7
(3)
θ ∼ t-3/7
(4)
and
In a recent article,17 a combined model was presented to describe the dynamics of spreading drops, allowing dissipation in both channels: friction against the substrate and hydrodynamic effects. It is predicted that, in the case of droplet spreading, a hydrodynamic regime is preceded by a regime where physicochemical aspects of the liquid and the solid are dominating. The identification of these two regimes is very important to assess the influence of the different parameters (viscosity, friction, ...) on the dynamics of wetting. It is precisely the aim of this article to verify experimentally the new model and to measure the influence of the parameters. Theoretical Background Let us first introduce the combined model of wetting.17 Consider a nonvolatile, sessile liquid drop of volume V, which is placed on an ideal, horizontal solid substrate, exposed to a neutral gas phase, for example, air. At time t ) 0 the initial shape of the drop is characterized by the base radius R and the contact angle θ (Figure 1). In accordance with the experiments presented later, we will suppose in what follows that θ is less than π/2. We are only concerned with the case of partial wetting. Consequently, the liquid drop will wet the solid substrate until it reaches its equilibrium shape with base radius R0 and contact angle θ0 (θ0 > 0) which obeys Young’s equation. It is also supposed that the drop, which is placed initially in a nonequilibrium configuration, retains the form of an ideal spherical cap at any moment in time. This implies that the instantaneous configuration of the drop can be totally described by a single parameter: the timedependent base radius R(t) or the dynamic contact angle
10.1021/la990769t CCC: $19.00 © 2000 American Chemical Society Published on Web 02/08/2000
2364
Langmuir, Vol. 16, No. 5, 2000
de Ruijter et al.
Figure 1. The schematic profile of a sessile drop, spreading on a solid substrate.
θ(t). For instance, when the dynamic contact angle is measured, the corresponding values of the instantaneous base radius can be calculated by virtue of the conservation of volume condition, which gives
R3(t) )
(3Vπ )Φ[θ(t)]
(5)
where
Φ[θ(t)] )
sin3 θ(t) 2 - 3 cos θ(t) + cos3 θ(t)
(6)
The drop is considered as a purely mechanical system, and the wetting is related to the total dissipation occurring during the change in the drop radius from the value R(t) to the value R(t) + δR(t). Within the spherical cap approximation, the change of the associated free energy δF is given by
δF[R(t)] ) -2pR(t)γ/(cos θ0 - cos θ(t))δR(t)
(7)
with γ the liquid-gas interfacial tension. Now, the relation between the driving force of spreading, expressed by the free energy change, and the energy dissipation, is provided by the standard mechanical description of dissipative systems dynamics.21 Since F(R) is not explicitly dependent on time, the corresponding dynamic equation will contain only two terms
∂F[R(t)] ∂T[R˙ (t);R(t)] )∂R˙ (t) ∂R(t)
(8)
where T[R˙ (t);R(t)] denotes the dissipation function, which describes the total dissipation in a circular liquid drop. The term on the right-hand-side determines the driving force of spreadingsthe derivative of the liquid/solid potential energy with respect to the base radius:
∂F[R(t)] ) -2πR(t)γ(cos θ0 - cos θ(t)) ∂R(t)
(9)
The dissipation function is represented as a sum of two different components
F[R(t)] ) Fl[R(t)] + Fw[R(t)]
(10)
where the first term on the right-hand side describes the dissipation occurring in the immediate vicinity of the contact line and the second one the losses due to viscous flow in the rest of the drop. We do not take into account in this approach the existence of a precursor film, which (21) Landau, L. D.; Lifschitz, E. M. Mecanique, 3rd ed.; MIR: Moscow, 1969.
Figure 2. Seaver-Berg approximation of the flow pattern in a sessile drop.
may appear in the complete wetting case. We have also checked in the time scale considered for our experiments, that there was no evidence of a precursor film during spreading. Dissipation in the vicinity of the contact line results from various physicochemical processes, which lead to the attachment of liquid molecules to the solid. Such a dissipation channel has been considered already in refs 18, 19, and 22 and subsequently scrutinized by several authors (see, e.g., refs 21 and 23). It is in general referenced to as the molecular-kinetic theory of wetting. At low velocities ν, this theory predicts that
R˙ (t) ) ν(t) )
γ (cos θ0 - cos θ(t)) ζ0
(11)
where ζ0 is considered to be a friction coefficient per unit length of the contact line.17 The leading contribution to Fl is then
Fl[R(t)] ) 2πR(t)
ν(t) γ(cos θ0 - cos θ(t)) ) 2 πR(t)ζ0Rν2(t) (12)
To calculate the contribution of the viscous flows to the total dissipation function, we followed the scheme introduced by Seaver and Berg.14 They assumed that the fluid dynamics of the spreading spherical cap can be approximated by that of a spreading cylindrical disk with height h* and volume equal to the volume of the drop (Figure 2). The height of the disk h* is not an independent parameter but is fixed by the volume and the radius
h* )
R(t) 3Φ[θ(t)]
(13)
where φ has been defined in eq 6. Under the assumption that the radial component vr of the fluid velocity is much greater than the vertical component, it was found that the dissipation due to radial flows is given by
Fw[R(t)] ) η
∫V( ∂hr) ∂ν
2
dV
(14)
with η the liquid viscosity. Integrating over the drop volume resulted in (22) Cerry, B. W.; Holmes, C. M. J. Colloid Interface Sci. 1969, 29, 17. (23) Ruckenstein, E.; Dunn, M. J. Colloid Interface Sci. 1976, 59, 137.
Evidence of Several Time Scales in Drop Spreading
( )
Fw[R(t)] ) 6πR(t)ηΦν(t)2 ln
R(t) a
Langmuir, Vol. 16, No. 5, 2000 2365
(15)
where the parameter a denotes the lower cutoff value of the radial coordinate. Without this cutoff the viscous dissipation would diverge, because the velocity at the symmetry axis of the drop has a singularity. This singularity is quite artificial, however, and arises only because of the simplified problem treated by Seaver and Berg. Clearly, instead of the singular behavior predicted by the equations, one may expect that the radial velocity is zero at the symmetry axis and remains small within some cylindrical region of radius a centered at this axis. This is precisely the meaning of the parameter asit is the radius of the core region, which is stagnant with respect to the radial dilation. Or, in other words, it is the region where the dissipation is negligible. Combining the two types of dissipation in eq 8, the desired dynamic equation was found17
ν(t) )
γ(cos θ0 - cos θ(t)) ζ0 + 6ηΦ[θ(t)] ln(R(t)/a)
(16)
or, in dimensionless units
Ca )
cos θ0 - cos θ(t)
(17)
δ + 6Φ[θ(t)] ln(R(t)/a)
where Ca ()ηv/γ) is the capillary number and the parameter δ is the ratio of the friction coefficient, for motion of the contact line on the solid substrate, to the bulk viscosity, δ ) ζ0/η. It is easily seen that if δ is large compared to 6Φ ln(R/a), this final equation reduces to the molecular-kinetic theory of wetting (eq 11). This will be the case if the frictional forces which act in the vicinity of the contact line are more important than the viscous forces. If, on the other hand, δ is small, the viscous forces are dominating. By expansion under the condition of small contact angles, eq 17 can be approximated by
θ3 - θθ20 ∝ Ca ln(R/a)
(18)
Figure 3. Schematic of the analysis assembly of spreading drop.
considered are 20 cP (referred to PDMS), 45-85 cP (referred to PDMSOH1), and 90-130 cP (referred to PDMSOH2). The liquids were chosen because they are nonvolatile at 3D. Throughout, we used monocrystalline (100) silicon wafers, provided by ACM. To prepare reproducible substrates, we followed a rigorous cleaning procedure.24-26 It consists of five steps. First, the wafers were sonicated during 3 min in two successive ultrasound chloroform baths. Then, they were placed under an ultraviolet lamp in a reactor saturated with ozone, to remove all organic contamination from the surface. Third, the inorganic contaminants are removed by immersion in a 7:3 mixture of sulfuric acid and hydrogen peroxide, at 150 °C. Afterward the wafers are rinsed with MilliQ water and dried under a flux of nitrogen. Finally, the wafers are placed a second time in the UV ozone reactor under a flux of oxygen saturated by water. The wafers were always cleaned just prior to the experiments. In Voue´ et al.,27 it is shown by complete spreading experiments that the wetting properties of wafers prepared using this procedure do not change in the time scale of our experiments. Methods
This result is typical for a hydrodynamic regime (see, for example, ref 8 or 25). Analyzing the asymptotic solutions of the dynamical equation (eq 16), we get that R(t) will behave as t1/7 (and accordingly θ(t) will behave as t-3/7) for times less than a certain t*, which is a function of several variables such as the volume of the drop the surface tension ... Above t*, we obtain that R(t) will behave as t1/10 (and accordingly θ(t) will behave as t-3/10). That means that for the spreading of a drop, a distinct molecular-kinetic regime is followed by a hydrodynamic regime. The characteristic time of the changeover between these two regimes is given by
( ( ))
V 1 t* ≈ γ-1V1/3(ζ0)10/3η-7/3 ln 3 2 a
-7/3
(19)
In the following sections, we present our experimental methods and results. Afterward we discuss the results in terms of the combined model. Materials The experiments were performed by spreading poly(dimethylsiloxane)s, with trimethyl or hydroxyl ends on top of bare silicon wafers. The liquids were provided by ABCR GmbH Germany and the viscosities of the liquids
To measure the contact angle relaxation of spreading drops over a wide range of times, we combined two methods (Figure 3). For the short times, the spreading of the drops was captured in profile by a video recorder connected to a CCD camera and a stereo binocular microscope (×8, Zeiss Stemi SV 11). The drops were back-lighted. With this method, we were able to capture 50 images of the relaxing drop per second. From each consecutive image, the edge of the drop and its reflection was located by an efficient edge detector.16,17 The drop profile was then interpolated with a Laplacian fit. The contact angle was determined at the point where the fit to the drop and its reflection intersected. This procedure was done automatically for each image. With this arrangement it was possible to measure contact angle relaxation directly from about 90° down to 5°, with a precision of about 1°. Simultaneously, the advancing contact line was captured on video from a top view with an interference (24) Dodge, F. T. J. Colloid Interface Sci. 1988, 121, 155. (25) Cox, R.G. J. Fluid Mech. 1986, 168, 169. (26) Broszka, J. B.; Shahidzadeh, N.; Rodelez, F. Nature 1992, 360, 719. (27) Voue´, M.; Valignat, M. P.; Oshanin, G.; Cazabat, A. M.; De Coninck, J. Langmuir 1998, 14, 5951.
2366
Langmuir, Vol. 16, No. 5, 2000
de Ruijter et al.
Figure 4. log θ versus log t for the spreading of PDMS.
Figure 6. log θ versus log t for the spreading of PDMSOH2.
between the contact angle θ and R (eqs 5 and 6)
R ) (3V/π)1/3Φ(θ) avec Φ(θ) )
(22)
sin θ (2 - 3 cos θ + cos3 θ)1/3
Here it is assumed that the shape of the drop is spherical and that the volume V is constant (the liquids have low volatility). These assumptions are valid for the experiments analyzed here.29 Differentiating this equation we find
dR 3V )dt π
(1 - cos θ)2
1/3
( )
3
(2 - 3 cos θ + cos θ)
Figure 5. log θ versus log t for the spreading of PDMSOH1.
4/3
dθ dt
(23)
Combining eqs 21 and 23, we obtain microscope (CF Plan NIKON 20× objective and standard 10× oculars) in combination with a CCD camera. Since the silicon surface is highly reflecting and the PDMS and PDMSOH are transparent, equal thickness interference fringes were observed near the edge. The distance between the fringes is a measure for the contact angle
tan θ ) pλr/2nx
Results In Figures 4-6, we give the contact angle relaxation for PDMS, PDMSOH1, and PDMSOH2 versus time. Discussion All the equations describing the different models can be expressed as
(21)
where f is independent or at most a weak function of the base radius R. In addition, we have the following relation (28) Greivenkamp, J. E., Jr. Physical Optics. In Handbook of Optics; Bass, M., Ed.; McGraw-Hill: NewYork, 1995; pp 2.1-2.43.
1/3
( )
Ω(θ)f(θ,θ0,γLV,η,parameters) (24)
where
Ω(θ) )
(20)
where p is the number of fringes, λr the wavelength of the reflected light, n the refractive index of the liquid, and x the distance between the fringes under consideration.28 The images were analyzed using our software. With this method contact angles were measured below 6°, with a precision of about 0.15°. We take great care that in the small overlap around 6°-5°, the two methods provided the same result.
ν ) dR/dt ) f(θ,θ0,γLV,η,parameters)
π dθ )dt 3V
(2 - 3 cos θ + cos3 θ)4/3 (1 - cos θ)2
To fit the experimental data to the combined model, the general equation used is
dθ π )dt 3V
1/3
( )
with
{
(cos A - cos θ) B + C Φ(θ)
Ω(θ)
(25)
A ) θ0 ζ0 B) γLV 6η R C) ln γLV a
()
where A, B, and C and the initial condition are the parameters to fit. In this equation, two parameters are a priori unknown, namely ζ0 and a. This differential equation is not analytically solvable. We therefore used the following calculation scheme: Initially, we guess the values of ζ0 and a. For these values, the numerical solution of eq 25 is calculated. The resulting curve is then compared to the experimental (29) Shikhmurzaev, Y. D. J. Fluid Mech. 1997, 334, 2-11.
Evidence of Several Time Scales in Drop Spreading
Langmuir, Vol. 16, No. 5, 2000 2367
Table 1. Parameters Obtained by Applying the Combined Model to the Data for PDMS, PDMSOH1, and PDMSOH2 η(cst)a PDMS PDMSOH1 PDMSOH2
20 87 111
ln(R/a)
ξ0
θ0 (deg)
t* (s)
0.38 ( 0.04 3.24 ( 0.15 0.1 ( 0.01 0.47 1.94( 0.06 4.69 ( 0.1 0.1 ( 0.01 0.24 4.6 ( 0.1 1.49 ( 0.1 0.1 ( 0.01 17
Table 2. Parameters Obtained by Applying the Molecular Kinetic Model, the Hydrodynamic Model, and the Combined Model to the Data for PDMS, PDMSOH1, and PDMSOH2
PDMS PDMSOH1 PDMSOH2
molecular kinetic ξ0
hydrodynamic ln(R/a)
1.16 6.81 8.29
1.68 2.81 1.03
combined model ξ0 ln(R/a) 0.38 1.94 4.6
3.24 4.69 1.49
contact angle relaxation under consideration. The values of ζ0 and a are modified in such a way that the distance between the numerical solution and the experiment is minimal. We wrote a FORTRAN program based on an adaptive step-size Runge-Kutta algorithm combined with a downhill Simplex method30 to perform this fitting procedure automatically. The errors on the best fitting parameters are calculated with the bootstrap method,30 using a standard deviation of the measured dynamic contact angle of 1° using G-contact and 0.15° using the interferometric method. With this powerful method, the original experimental data are used as the basis for a Monte Carlo simulation. From the original data set, 37% (≈1/e) randomly chosen points are replaced by duplicates. The duplicate points are chosen according to a normal distribution function, with the original data as average and the expected error on each datum point as the standard deviation. In this way, we replace the original set of data with a new one, for which the corresponding parameters are then calculated with the method described above. Successive cycles result in a simulated set of values for each parameter. If enough cycles are used (>100), these sets turn out to be normally distributed, providing us with the mean value and standard deviation for each parameter. All the values of the parameters and their predicted errors given below were calculated with this method (Table 1). The resulting best fits are shown in Figures 4-6. It is seen that the friction between the liquid and the solid substrate increases with viscosity. The cutoff ln(R/ a) has a typical value of 1 to 4, which makes sense with many other experiments. The validity of this model is confirmed if we compared our data either to the molecular kinetic theory or to the hydrodynamic model as given in Table 2. The friction ξ0 calculated with the molecular kinetic model follows the same trend as the values calculated by the combined model, but the values are a factor 2 or 3 higher. This can be explained as follows: The molecular kinetic theory accounts only for the dissipation in the threephase zone; by omitting the dissipation due to the flow, the dissipation in that regime is thus overestimated. We observe that the values of the cutoff associated to the hydrodynamic model decrease from 10 -4 to 10 -6 m by taking into account the friction contribution. The quality of the fit is also sensitive to the considered theory. Let us introduce n
1 χ2 ) ni)1
∑
(
)
θiexp - θith σθ
2
Figure 7. cos θs - cos θd versus velocity for different PDMS liquids and the fit of the data with the combined theory (cos θs ) 1 for all data), from ref 30. (θs is the static contact angle and θd is the dynamic contact angle.) Table 3. Calculated χ2 for the Molecular Kinetic Model, the Hydrodynamic Model, and the Combined Model PDMS PDMSOH1 PDMSOH2
molecular kinetic
hydrodynamic
combined
12.64 31.26 15.47
1.67 1.89 7.75
1.48 1.61 2.46
Table 4. Parameters Resulting from Applying the Combined Model to the Data of Sauer et al.30 PDMS 3.9k PDMS 17k PDMS 71k
η (cP)
ξ0
ln(R/a)
50 500 1300
0.97 ( 0.007 10.6 ( 0.2 98 ( 4
0.40 ( 0.1 0.85 ( 0.1 0.85 ( 0.1
where θiexp (θith) refers to the experimental measurement (theorical prediction) and where σθ represents the estimated error of the measurement: 1° with G-contact and 0.15° with the interferometric method. We then get Table 3. The lower values of χ2 obtained using the combined theory support its validity. This combined model provides also an explanation for the different power laws, which have been reported in the literature for advancing contact angles. When concentrating on systems far from equilibrium, as is mostly the case in forced wetting experiments, the molecular-kinetic theory should fit the data very well. On the other hand, when systems are studied close to equilibrium, we are likely to find the power law predicted by the hydrodynamic models, or at very long times, to find an exponential behavior. For receding contact angles, the effect will be opposite. For systems far from equilibrium, the contact angle will be small, and the dissipation due to viscosity correspondingly more important. Closer to equilibrium, the molecular-kinetic model should fit the data better. This behavior has been reported in a recent article by Sauer et al.31 In a series of experiments with PDMS of different viscosity, the authors recognized an exponential decaying regime and a hydrodynamic regime. In their experiments, they also observed a regime which they supposed was due to “viscoelastic processes near the solid”. In this regime, the contact line velocity is more or less linear with the driving force. This behavior is characteristic of the molecular-kinetic model. Thus, our combined model should be able to predict all three different regimes using (30) Press: W. H.; Teukolsky, S. A.; Vetterling, W. T.; Flannery, B. P. Numerical Reciepes in Fortran, 2nd ed.; Cambridge University Press: Cambridge, 1992. (31) Sauer, B. B.; Kampert, W. G. J. Colloid Interface Sci. 1998, 199, 28.
2368
Langmuir, Vol. 16, No. 5, 2000
just one equation. We have thus applied our equation (eq 25) to their data. The results are given in Figure 7. The associated fitted parameters are given in Table 4. We again observe the relation between viscosity and friction. Let us here mention that the small values of ln(R/a) may be due to the cylindrical geometry which refers to the considered fibber but the quality of the fit is indeed quite remarkable (Figure 7). Conclusions We have thus shown using several independent experiments the validity of the combined model of spreading allowing dissipation of energy by friction against the solid
de Ruijter et al.
surface or by the hydrodynamic effects. The fact that we obtain good agreement between the proposed theory and the data over 6 decades of time is a very strong argument in favor of the theory. We have also shown how to estimate the time window in which each particular regime should hold preferentially. Acknowledgment. This research has been partially supported by the Ministe`re de la Region Wallonne. aThe viscosities were determined with a Brookfield cone/ plate viscometer and corresponded with the data given by the manufacturer. LA990769T