Langmuir 2004, 20, 8385-8390
8385
Liquid Coating of Moving Fiber at the Nanoscale David Seveno,* Gregory Ogonowski, and Joe¨l De Coninck Centre for Research in Molecular Modeling, Materia Nova, University of Mons-Hainaut, Av. Copernic, 1, B-7000 Mons, Belgium Received February 18, 2004. In Final Form: June 29, 2004 Using large scale molecular dynamics, we study the contact line motion of a liquid meniscus crossed by a moving nanofiber. Varying the amplitude of the liquid/solid interactions, we analyze the shape of the meniscus versus time for a range of velocities. The associated contact angles are estimated by fitting the profiles using the James equation. The corresponding flux lines describing the displacement of the liquid molecules inside the meniscus have also been measured. The analysis of the dynamic contact angle is in agreement with the molecular-kinetic theory and confirms the existence of an optimal speed for wetting.
1. Introduction Micro- or Nanotechnology has become progressively a very important area of research. Unfortunately, experimental techniques for the study of the nanoscopic characteristics of coating are rare. Concerning the flat susbtrate case, important progress has been achieved in particular using ellipsometric techniques1,2 and molecular dynamics.2-4 For other geometries, such as fibers, the problem is still open at that scale. The thermodynamic equilibrium properties of liquid films and droplets on fiber have already been examined theoretically and experimentally,5-6 but here we are interested in nonequilibrium phenomena. The case of the spontaneous rise of a liquid around a fiber has been examined experimentally7,8 and by very large scale molecular dynamics simulation.9 It has been shown in particular that during the process dissipation should occur both in the liquid (hydrodynamic description10-16), due to the reorganization of molecules as described by viscosity, and near the solid as described by some kind of friction (molecular-kinetic description17,18). The latter is necessarily dependent on the microscopic characteristics of the system. The range over which the different regimes operate may change depending on the specific macroscopic and microscopic characteristics of the system. Up to now, the (1) Cazabat, A. M.; Gerdes, S.; Valignat, M. P.; Vilette, S. Interface Sci. 1997, 5, 129. (2) Cazabat, A. M.; Valignat, M. P.; Vilette, S.; De Coninck, J.; Louche, F. Langmuir 1997, 13, 4754. (3) Blake, T. D.; Clarke, A.; De Coninck, J.; de Ruijter, M. Langmuir 1997, 13, 2164. (4) de Ruijter, M.; Blake, T. D.; De Coninck, J. Langmuir 1999, 15, 7836. (5) Neimark, A. J. Adhesion Sci. Technol. 1999, 13, 1137. (6) Que´re´, D.; Di Meglio, J. M.; Brochard-Wyart, F. Rev. Phys. Appl. 1988, 23, 1023 (7) Vega, M. J.; Seveno, D.; Lemaur, G.; Ada˜o, M. H.; De Coninck, J. Submitted for publication. (8) Clanet, C.; Que´re´, D. J. Fluid Mech. 2002, 460, 131. (9) Seveno, D.; De Coninck, J. Langmuir 2004, 20, 737. (10) de Gennes, P. G. Rev. Mod. Phys. 1985, 57, 827. (11) Brochard-Wyart, F.; de Gennes, P. G. Adv. Colloid Interface Sci. 1992, 1, 39. (12) Dussan, E. V. J. Fluid Mech. 1976, 77, 665. (13) Voinov, O. V. Fluid Dyn. 1976, 11, 714. (14) Tanner, L. H. J. Phys. D 1979, 2, 1473. (15) Hoffmann, R. J. Colloid Interface Sci. 1975, 50, 228. (16) Cox, R. G. J. Fluid. Mech. 1986, 16, 169. (17) Blake, T. D.; Haynes, J. M. J. Colloid Interface Sci. 1969, 30, 421. (18) Blake, T. D. Wettability; Berg, J. C., Ed.; Marcel Dekker: New York, 1993.
case of moving fiber has been considered, among others, by Inverarity,19 Que´re´,20 and Schneemilch et al.21 Our goal is to use very large scale molecular dynamics simulations to study in detail the dynamics of wetting on nanometric moving fibers and thus to study the optimization of the associated coating process. This technique has, in particular, the great advantage of allowing us to change a single parameter such as the liquid/solid affinity, which is never possible in a real experiment, since changing the solid will also change many other factors. The paper is organized as follows. Section 2 is devoted to the molecular modeling of the capillary rise around a fiber. The results associated with the dynamics of wetting are given in section 3. 2. The Model System In our simulations, all potentials between atoms, solid as well as liquid, are described by standard pairwise Lennard-Jones 12-6 interactions
(( ) ( ) ) σij d
Vij(d) ) 4ij
12
-
σij d
6
(1)
where d is the distance between any pair of atoms i and j. The parameters ij and σij are in the usual manner related to the depth of the potential well and the effective molecular diameter, respectively. Translated into reduced (dimensionless) units, eq 1 becomes
(
Vij*(d*) ) 4
Cij 12
d*
-
Dij d*6
)
(2)
where the asterisk stands for reduced units. For simplicity, Cij and Dij are chosen constant for each type of atom. We choose Cff ) Dff )1.0, Css ) Dss ) 1.0, and Csf ) Dsf ) 0.8, 0.9, 1.0, and 1.05 where the subscripts stand for fluid/ fluid (ff), solid/solid (ss), and solid/fluid (sf) interactions. These coefficients were chosen according to the previous studies of the spreading of a droplet on a flat surface4 and of spontaneous rise around a fiber.9 The intrafluid coefficients are standard, and the solid/solid coefficients are chosen to produce a stable lattice structure at the (19) Inverarity, G. Br. Polym. J. 1969, 1, 247. (20) Que´re´, D. J. Fluid Mech. 1999, 31, 347. (21) Schneemilch, M.; Hayes, R.; Petrov, J.; Ralston, J. Langmuir 1998, 14, 7047.
10.1021/la049574y CCC: $27.50 © 2004 American Chemical Society Published on Web 08/18/2004
8386
Langmuir, Vol. 20, No. 19, 2004
temperature of interest. The choice of the solid/fluid interactions ensures that the liquid wets the solid.3-4 The theoretical range of the Lennard-Jones 12-6 interactions extends to infinity. Strictly, one should therefore evaluate the interactions between all possible pairs in the system. Fortunately, the interaction potentials decrease rapidly as the distance becomes large. We therefore apply a spherical cutoff at 2.5σij, and shift the potential so that the energy and force are continuous at d* ) 2.5. As a result, we consider only short-range interactions in these simulations. We simulate a molecular structure for the liquid by including a strong elastic bond between adjacent atoms within a molecule, of the form Vconf ) Dconfd6 with Dconf ) 1.0. The liquid molecules are always 16 atoms long. This extra 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 initial face-centered cubic (fcc) lattice configuration, to give a realistic atomic representation of the solid fiber. Thanks to this added potential, the fiber structure is stable. To summarize, we consider a very simple chain-like liquid system, made by 16 monomers with spherical symmetry, in contact with a solid fiber made from a fcc solid lattice. We always apply a time step of 5 fs22 during our simulations with ij ) 0.276 × 103 J‚mol-1 and σij ) 3.5 Å, typical of carbon atoms. We consider a fiber with radius 20 Å and a meniscus of liquid molecules that we equilibrate inside an annulus of frozen liquid atoms. The simulation consists of a first part where the solid (16 400 atoms) and the liquid (40 288 atoms) are placed far from each other, so that they can equilibrate independently, which takes around 1 000 000 time steps. During this interval the liquid and the solid are kept at the same temperature. Then, through a small hole into the center of the annulus, we bring the fiber and the meniscus into contact and let all the elements equilibrate together. From that moment on, the temperature of the liquid is free to evolve, whereas the temperature of the solid is held constant to avoid heating of the system. Although the system is rather simple, it contains all the basic ingredients to model the wetting of a fiber. Once we have reached equilibrium, we then move the fiber through the meniscus at a constant speed and measure the dynamic contact angle around the fiber. Thanks to the periodic boundary conditions applied in the z-direction, we simulate a fiber of infinite length. There is no periodic boundary condition in the other directions. A typical side view of the process is given in Figure 1. To compute the associated contact angle for the liquid around the fiber, we first subdivide the liquid annulus into several concentric cylindrical shells of arbitrary thickness. The constraint on the number of shells is provided by the need to maximize their number while ensuring that each shell contains enough molecules to give a uniform density. For each shell, we compute the density of particles as a function of the distance x to the fiber. We then locate the extremity of the shell at the distance where the density falls below a cutoff value of 0.5 times the liquid density. To check the consistency of the method, different shell thickness and cutoff values were considered and these gave almost identical results. The James equation23 is then used to fit the meniscus profiles. This equation describes the equilibrium shape of (22) Allen, M. P.; Tildesley, D. J. Computer Simulations of Liquids; Oxford University Press: Oxford, 1987. (23) James, D. F. J. Fluid Mech. 1974, 63, 657.
Seveno et al.
Figure 1. Side view of fiber wetting with Csf ) Dsf ) 0.8 and a fiber velocity of 5 m‚s-1.
such menisci fairly well, and by using it we assume local equilibrium at each time step
(
yt(xt) ) r0 cos(θt)* ln(4) - E ln
( (( ) xt + r0
xt r0
2
) ))
- cos2(θt)
1/2
(3)
where yt is the height of the meniscus, θt is the contact angle, E is Euler’s constant, r0 is the fiber radius, and xt is the abscissa of the point considered. This method enables us to construct the complete profile of the meniscus and to determine how it evolves with time. The best fits through the profiles described by eq 3 were always located within the region where the density dropped from 0.75 to 0.25, except in the first few molecular layers in contact with the fiber. This confirms that the simulated menisci always retain their equilibrium form during wetting, except very close to the solid fiber where the fluctuations can be large. To avoid this problem, we investigated the averaged profile over 5 frames every 10 000 time steps as a function of the number of shells. We find that the fit using all the experimental points leads to stable results. By fitting the profiles with eq 3, we are thus able to measure the contact angle θ as a function of the number of time steps. 3. The Results In Figure 2, we have plotted the averaged profiles as described before for Csf ) Dsf ) 0.8 and a fiber velocity of 5 m‚s-1. We observe that the liquid shape, initially with a contact angle of 80°, is modified by the displacement of the fiber and that the global shape becomes elongated in the direction of the fiber displacement. Simulations were launched for a range of velocities from 0.02 to 15 m‚s-1. Advancing and receding angles (see Figure 2) were obtained by averaging the upper and lower profiles. The time evolution of the contact angle for several fiber speeds is plotted in Figure 3. After sufficient time at each fiber velocity, we can determine the stationary advancing and receding dynamic contact angles, by averaging the contact angle values over the last few nanoseconds (typically 2 ns) during which the contact angles do not vary by more than 5%. These angles
Liquid Coating of Moving Fiber
Figure 2. Interface profiles during fiber displacement. The upper meniscus exhibits an advancing angle with the fiber, the lower meniscus a receding one.
Langmuir, Vol. 20, No. 19, 2004 8387
Figure 4. Stationary contact angle versus fiber velocity with Csf ) Dsf ) 0.8 (circles), 0.9 (squares), and 1.0 (triangles). The error bars represent the standard deviation of the data averaged over the last 2 ns. Table 1. Dynamic Contact Angle for Systems of Increasing Sizesa no. of liquid atoms
dynamic contact angle (deg)
40 288 101 648 140 096
113.53 ( 6.44 111.54 ( 6.04 117.70 ( 5.40
a The error represents the standard deviation of the data averaged over the last 2 ns.
Figure 3. Dynamic contact angle versus time for different fiber velocities for Csf ) Dsf ) 0.8. Increasing contact angles are advancing ones with a fiber velocity of 10 m‚s-1 (triangles), 5 m‚s-1 (circles), and 2 m‚s-1 (squares) whereas decreasing ones correspond to receding contact angles with fiber velocities of 2 m‚s-1 (diamonds) and 5 m‚s-1 (stars).
are plotted versus the fiber velocity in Figure 4 for Csf ) Dsf ) 0.8, 0.9, and 1. As can be seen, all the advancing contact angles become more-or-less identical at sufficiently high speeds. This is in contrast to the equilibrium case where the contact angle is a clear function of the solid/liquid interaction.9 By extrapolation, these data also indicate that the limiting speed for coating before reaching an angle of 180° is, here, around 11 m‚s-1. To check the consistency of the method, we have compared our contact angle measurements with those obtained for larger systems, with 101 648 and 140 096 liquid atoms instead of 40 288. For a velocity of 4 m‚s-1 and Csf ) Dsf ) 0.8, we get the advancing dynamic contact angles given in Table 1. These data indeed confirm that we are considering a system large enough to avoid finite size effects. Let us now study in detail the flux lines in the meniscus. This will allow us to analyze the convection mechanisms in the vicinity of the fiber. Since our system is rotationally invariant, we project the volume of the meniscus onto the x-y plane. We then subdivide this plane into small square units or bins, so as to make a grid. The size of the units in the grid is constrained by the need to maximize the number of bins,
Figure 5. Velocity field between 0.75 and 1 ns.
while ensuring that every bin contains enough atoms to give reproducible results. For every bin we compute the center of mass of the atoms which are part of that bin. A short time later (typically 50 000 computer time steps, i.e., 0.25 ns), we determine the net displacement of the center of mass of each bin. This allows us to measure the velocity field of the meniscus in some coarse grained sense. Assuming the density in the meniscus is uniform, the velocity field may also be considered as lines of flux (streamlines). In Figures 5-8, velocity fields are shown for the menisci with Csf ) Dsf ) 0.8 and a fiber velocity of 5 m‚s-1 at different time steps. In Figures 5 and 6, the velocity field of the meniscus during the initial stages of wetting is shown. Every arrow represents the net displacement of a bin over a time period of 50 000 computer time steps (0.25 ns). The length of the arrow is a relative measure of the local velocity. The
8388
Langmuir, Vol. 20, No. 19, 2004
Figure 6. Velocity field between 2 and 2.25 ns.
Figure 7. Velocity field between 3 and 3.25 ns.
Seveno et al.
Figure 9. Cosine of the contact angle versus the fiber velocity with Csf ) Dsf ) 0.8 (lower triangles), 0.9 (squares), 1.0 (upper triangles), and 1.05 (circles). The full lines are the best fits corresponding to eq 4.
Figures 7 and 8, so that the average velocity field versus time for the highest layer remains similar from the beginning of the simulation to the stabilization of the meniscus. In summary, the flow of the meniscus is basically circular with not much variation in the magnitude of the velocity. Closer to the fiber surface, the velocity field shifts its direction toward that of fiber displacement. This results in a ribbon of flow, which follows the liquid-vapor interface. There is a gradient of velocity, almost zero in the middle of the meniscus, and maximum close to the contact line. It seems that the liquid is free to flow toward the three-phase zone and that the main dissipation occurs in this region. For high-speed wetting with large dynamic contact angles, it is expected that the molecular-kinetic theory should hold.17,18 Molecular-Kinetic Model. In its general form, the molecular-kinetic model relates the dynamic contact angle to the fiber velocity to the dynamic contact angle as follows (The inverse of the hyperbolic function that appears in ref 18 has been expressed analytically.)
cos(θt) ) cos(θ0) -
(
2nkBT V ln + γ 2K0λ
(( ) ) ) V 2K0λ
2
+1
1/2
(4)
Figure 8. Velocity field between 4.25 and 4.5 ns.
receding contact angle is here less than 80° for Csf ) Dsf ) 0.8. The velocity field is basically circular in the meniscus and strongly driven by the moving fiber. Close to the fiber surface and away from the three-phase line, the velocity is identical to the fiber velocity, indicative of a no-slip boundary condition during the displacement of the fiber. Regular and constant liquid recirculation is clearly shown in the whole meniscus. The first layer of liquid seems to be continuously entrained by the fiber, as can be seen in
with θt the dynamic contact angle, θ0 the equilibrium contact angle, γ the liquid-gas interfacial tension, n the number of adsorption sites per unit area, λ the typical length of each molecular displacement, K0 the jump frequency, kB the Boltzmann’s constant, T the temperature, and V the fiber velocity. To a first approximation n and λ are related by n ) λ-2. If we now plot the cosine of the contact angle versus the fiber speed, we get a nonlinear function as indicated in Figure 9. The results of fitting eq 4 to the data using LevenbergMarquad algorithm24 are given in Table 2. Assuming that λ is the distance between two solid atoms, which is exactly 3.5 × 10-10 m, we derive Table 3 from Table 2. (24) More´, J. J.; Garbow, B. S.; Hillstrom, K. E. User guide for MINPACK-1; Argonne National Laboratory Report ANL-80-74, 1980.
Liquid Coating of Moving Fiber
Langmuir, Vol. 20, No. 19, 2004 8389
Figure 10. Density profile as a function of the distance from the fiber for Csf ) Dsf ) 0.8. The inset gives the density profile close to the fiber. Table 2. Fitted Parameters Corresponding to Equation 5 for the Data of Figure 9 Csf ) Dsf
cos θ0
2nkBT/γ
1/(2K0λ) (m-1‚s)
quality of the fit (R2)
1.05 1.0 0.9 0.8
0.87 ( 0.02 0.68 ( 0.01 0.39 ( 0.01 0.17 ( 0.01
0.67 ( 0.04 0.75 ( 0.05 0.71 ( 0.06 0.67 ( 0.05
0.92 ( 0.13 0.54 ( 0.06 0.40 ( 0.05 0.27 ( 0.03
0.997 0.993 0.993 0.998
Table 3. Physical Parameters Derived from Table 2 Csf ) Dsf
θ0 (deg)
K0 (109 s-1)
1.05 1.0 0.9 0.8
29.54 ( 2.33 47.16 ( 0.78 67.05 ( 0.62 80.21 ( 0.58
1.55 ( 0.22 2.65 ( 0.30 3.57 ( 0.45 5.29 ( 0.60
One of the considerable advantage of the molecular dynamics is that the jump frequency K0 can be measured directly from the simulations. Following the procedure described by de Ruijter et al.,4 two frequencies can be distinguished, namely, those from jumps parallel to the solid fiber (displacement within the same liquid layer) and perpendicular to the solid fiber (displacement from one layer to another). To compute such frequencies, the density profile inside the liquid must first be examined to locate the liquid layers. Figure 10 gives the density of the liquid as a function of the distance from the fiber. It is obvious that the liquid shows strong layering close to the fiber. The parallel frequency, f|, then refers to displacement of one atom within the first layer and the perpendicular one, f⊥, to the displacement between the first and the second layer. One atom is considered to stay within a layer if its perpendicular travel is less than the width of the density peak assessed at its half-height. Five layers away from the solid fiber, the density becomes uniform and equal to 331.93 ( 4.19 kg‚m-3. This less than half the real value of hexadecane (770.0 kg‚m-3) but in agreement with previous studies.4 The unit length for the parallel frequency is λ, and the distance between the first and the second layer for the perpendicular one is 0.86λ. To account for a collective action of the atoms, we assumed that the average characteristic time is the time required for half atoms of a layer to travel one unit distance. The associated frequency is then the inverse of this time. The results of these computations are compared with the K0 values obtained by fitting eq 4 to the dynamic contact angle data (Table 4).
Figure 11. Plot of (ln K0) (squares) versus (1 + cos θ0) for the different values of Csf and Dsf. (ln f|) (triangles) is also drawn. The straight full line is the best linear fit of eq 5 applied to the data of Table 3 with ln(a) ) 24.24 ( 0.33 and b ) 1.58 ( 0.22. (R2 ) 0.95.) The dotted line represents eq 7 with γ/nkBT ) 2.86 ( 0.2, η ) 0.035 mPa‚s, and ν ) 358.35 Å3. Table 4. Comparison between the Frequencies Obtained from the Fit of Equation 4 and These Obtained Directly from the Measurement of the Displacements Observed in the Simulations Csf ) Dsf
K0 (109 s-1) (fit)
f| (109 s-1)
f⊥ (109 s-1)
1.05 1.0 0.9 0.8
1.55 ( 0.22 2.65 ( 0.30 3.57 ( 0.45 5.29 ( 0.60
1.65 ( 0.22 2.26 ( 0.32 3.5 ( 0.66 3.31 ( 0.8
17.32 ( 1.18 21.27 ( 1.24 24.55 ( 1.82 32.48 ( 2.4
From a schematic viewpoint, to advance the contact line of unit distance, the hole created by the displacement of one atom within the first layer (of constant density) has to be filled in by an atom coming from the second layer. The overall frequency is then a combination of both types of frequencies. However it turns out that here, the parallel frequency is the lowest frequency and is therefore the limiting factor. It is worth noting that in ref 4, the limiting factor proved to be the perpendicular frequency. This change is understandable if we consider that our solid fiber has weaker solid/solid interaction parameters than previously. In ref 4, Css ) 35.0 and Dss ) 5.0, whereas, here Css ) 1.0 and Dss ) 1.0. We are currently carrying out simulations with flat surfaces and the lower solid interaction parameters to study any possible influence of geometry. In any case, the agreement between this measured frequency and K0 is good, strongly supporting the evidence of the dissipation at the microscopic scale. Plotting (ln K0) versus (1 + cos θ0), gives us a noticeable linear relation as represented in Figure 11, indicating that K0 can be written in the form
K0 ) a exp(-b(1 + cos θ0))
(5)
This results is in broad agreement with the theoretical prediction of Blake and De Coninck.25 Following ref 25 and combining eqs 4 and 5, we can write
[
V ) 2λa exp(-b(1 + cos θ0)) sinh
]
γ(cos θ0 - cos θ) 2nkBT (6)
(25) Blake, T. D.; De Coninck, J. Adv. Colloid Interface Sci. 2002, 96, 21.
8390
Langmuir, Vol. 20, No. 19, 2004
Seveno et al.
K0 )
Figure 12. Predicted maximum speed of wetting as a function of the stationary contact angle calculated from eq 6 using the fitted parameters given in the caption to Figure 11. The maximum occurs at a stationary contact angle of 87° for a speed of 9.55 m‚s-1.
This equation, which describes the fiber velocity dependence of θ according to the molecular-kinetic theory, implies that there is a maximum speed of wetting (just before air entrainment) for a particular liquid/solid affinity (θ0 ) 87°, V ) 9.55 m‚s-1). The result of applying eq 6 using the data from the simulation is illustrated in Figure 12. It is worth noting that in ref 25, Blake and De Coninck further assume that the reversible work of adhesion is related to the surface component of the specific activation free energy of wetting (per unit area). For a one-to-one relationship, they derived
[
kB T γ exp (1 + cos θ0) ην nkBT
]
(7)
where ν is the volume of the “unit of flow” (molecular volume) and η is the liquid viscosity. For our simple liquid, the unit of flow is a single molecule, so that it equals the molecular volume. Equation 7 explicitly accounts for the effect of liquid viscosity. However, from Table 2, we get that γ/nkBT ) 2.86 ( 0.2, which is about twice value of the fitted parameter b. Thus, it seems that the 1:1 ratio between the work of adhesion and the surface component of the specific activation free energy of wetting (per unit area) is not verified. Nevertheless, a close relation between the work of adhesion and K0 is clearly established. Concluding Remarks. We have used molecular dynamics to simulate contact line motion of a liquid meniscus around a moving nanofiber. For several liquid/solid interactions, we have studied the shape of the meniscus versus time for a range of velocities. The associated contact angles have been estimated locally by fitting the liquid profiles using James’s equation. The corresponding flux lines describing the displacement of the liquid molecules inside the meniscus have also been computed. The analysis of the dynamic contact angle data is in agreement with the molecular-kinetic theory. Our data also support a close relationship between the jump frequency K0 and the energy of adhesion between the liquid and the fiber as predicted by Blake and De Coninck.25 This leads to the prediction that there is an optimal speed for the wetting of the fiber just before air entrainment. Acknowledgment. This research has been partially supported by the Structural European Funds and by the Re´gion Wallonne. LA049574Y
