Home | About | Journals | Submit | Contact Us | Français |

**|**Scientific Reports**|**PMC5589961

Formats

Article sections

Authors

Related links

Sci Rep. 2017; 7: 10880.

Published online 2017 September 7. doi: 10.1038/s41598-017-11350-6

PMCID: PMC5589961

0000 0001 0694 4940grid.438526.eDepartment of Mechanical Engineering,Virginia Polytechnic Institute and State University, Blacksburg, VA 24061 USA

Jiangtao Cheng, Email: ude.tv@tjgnehc.

Received 2017 May 10; Accepted 2017 August 23.

Copyright © The Author(s) 2017

In this paper, we report molecular kinetic analyses of water spreading on hydrophobic surfaces via molecular dynamics simulation. The hydrophobic surfaces are composed of amorphous polytetrafluoroethylene (PTFE) with a static contact angle of ~112.4° for water. On the basis of the molecular kinetic theory (MKT), the influences of both viscous damping and solid-liquid retarding were analyzed in evaluating contact line friction, which characterizes the frictional force on the contact line. The unit displacement length on PTFE was estimated to be ~0.621 nm and is ~4 times as long as the bond length of C-C backbone. The static friction coefficient was found to be ~10^{−3} Pa·s, which is on the same order of magnitude as the dynamic viscosity of water, and increases with the droplet size. A nondimensional number defined by the ratio of the standard deviation of wetting velocity to the characteristic wetting velocity was put forward to signify the strength of the inherent contact line fluctuation and unveil the mechanism of enhanced energy dissipation in nanoscale, whereas such effect would become insignificant in macroscale. Moreover, regarding a liquid droplet on hydrophobic or superhydrophobic surfaces, an approximate solution to the base radius development was derived by an asymptotic expansion approach.

The dynamic wetting of liquids on solid surfaces is a ubiquitous phenomenon in nature and has many applications in industry. Familiar examples associated with natural processes are capillary suction in plants and self-cleaning of lotus leaves, i.e., lotus effect; and a wide range of practically important applications including surface coating^{1}, phase change heat transfer enhancement^{2, 3} and surface treatment-assisted flotation^{4}. The early attempts on dynamic wetting research can be traced back to the study of capillary tubes about two hundred years ago. During the last two decades, the long-lasting interest in dynamic wetting has been further stimulated due to the substantial progress in micro/nano-fluidics^{1} and microelectromechanical systems (MEMS) processing^{5, 6} that enables advanced surface treatments and micro\nano-structure characterization^{7, 8}. According to the continuum theory of fluid mechanics, dynamic wetting is a thermodynamic process governed by the capillary force, viscosity and gravity^{9}. However, the hydrodynamic analysis of the dynamic wetting yields unphysical predictions regarding the motion of three-phase (solid, liquid and vapor) contact line. When a contact line is moving at a spreading velocity *u*_{c}, the classical non-slip boundary condition renders the velocity gradient at the three-phase contact line unbounded and hence leads to a shear stress singularity^{10, 11}. In a manner analogous to the boundary layer theory, some researchers truncated the velocity profile to a slip length *L*_{s}where the continuum concept breaks down and the Navier-Stokes equations can then be solved by relaxing the no-slip condition in the vicinity of contact line^{12–14}. This approach is essentially a mathematical treatment to the microscopic phenomenon within the continuum domain and doesn’t provide any physical interpretation to the origin of interfacial slip. Nonetheless, there is little consensus on the order of magnitude of *L*
_{s}
^{10, 15} and the influence of intrinsic properties of the solid surface on wetting is still not well understood. Therefore, the dynamic wetting needs to be further addressed in the sub-continuum domain, *i.e*., at the molecular level.

In addition to the viscous dissipation and adhesion energy that are manifest in the macroscopic processes^{16}, the excess surface free energy can be dissipated in the form of contact line friction, which becomes more prominent in the molecular level or in nanoscale, as a result of the oscillations of contact line molecules from their equilibrium positions^{17}. The contact line friction^{15, 18, 19} is an important concept in the molecular kinetic theory (MKT), a conceptualized model that is based on a series of phenomenological parameters such as the equilibrium frequency *K*
_{0} and the unit displacement length λ^{20}. The basic idea underlying MKT is that the contact line motion is essentially a rate process controlled by the corresponding energy barriers^{20}. The continuous and macroscopic displacement of a contact line results from the collective manner of discretized forward or backward jumps of fluid molecules within the contact zone on a solid surface. In this scheme, the spontaneous spreading of contact line stems from the tilted energy barrier that favors the unidirectional displacement of contact line molecules. Hence MKT provides an efficient approach to tackling the interfacial slip paradox in dynamic wetting.

In spite of the tremendous research efforts in studying the dynamic wetting process at the sub-continuum level, there continues to be a lack of understanding on the physical meanings of molecular kinetic parameters and the dominant factors governing the dynamic wetting process still remain elusive. de Ruijter *et al*.^{21} applied the MKT to the spontaneous spreading of liquid droplets on various solid materials and found the unit displacement length λ was around 1 nm. The MKT was also applied in processes involving dewetting^{22}, wetting on chemically heterogeneous surfaces^{23}, forced reactive wetting^{24}, and even to solid-liquid-liquid systems^{25, 26} containing ionic liquids^{27}. In these studies, MKT-based predictions of some macroscopic parameters yielded satisfactory agreement with the experimental values. Besides, the MKT was coupled with molecular dynamics (MD) simulations and was reported to be able to depict the dynamic behaviors of nano-droplets^{17, 21, 28–30}. However, a recent review^{31} pointed out that the predictive power of the MKT is limited due to the insufficient knowledge on the physical interpretation of the molecular kinetic parameters, which are hard to be related to experimentally accessible quantities or material properties. Therefore, a meticulous analysis of the molecular kinetics in dynamic wetting process is needed to understand the physical meanings of MKT parameters and the fundamental mechanisms governing wetting dynamics.

Another topic of vast scientific interest is the contact zone formation, a multiscale phenomenon spanning from the molecular scale to the macroscopic scale. It is widely accepted that the macroscopic (apparent) contact angle is accompanied by distinct and velocity-dependent microscopic contact angles^{14, 32}. de Gennes *et al*.^{33} proposed that the contact zone could be divided into four subdomains, i.e., the molecular region, the proximal region, the central region and the distal region. Each region falls within different length scales and exhibits distinct curvatures controlled by *in situ* forces. In a recent experimental study, Chen and coworkers^{34, 35} reported a convex nanobending of contact line in the proximal region and concluded that the microscopic contact angle in MKT actually evolves beyond the molecular region of only several molecules thick. A follow-up study by Lukyanov and Likhtman^{32} demonstrated that the microscale contact angle forms as a result of the nonlinear distribution of the frictional force at the solid surface and hence the Young’s equation needs to be modified in order to be applied at nanoscale.

Via MD simulation in this study, we aim to provide an exhaustive analysis on the molecular kinetics of the dynamic wetting of water especially on hydrophobic surfaces, on which only limited wetting studies have been carried out. For this purpose, the solid material was chosen to be the amorphous polytetrafluorethylene (PTFE) because PTFE is the main component of Teflon^{®}, which is one of the most widely used hydrophobic coating materials in industry. In particular, adoption of the real chemical structure of PTFE can effectively eliminate the possible artifacts that would otherwise be introduced by artificially-assembled solid structures. From an *ab initio* perspective, the dynamic wetting of water on hydrophobic surfaces is scrutinized at the molecular level in this work. Besides, the size effect of nano-droplets on the molecular kinetic parameters were investigated and the inherent contact line fluctuations were examined. Finally, an asymptotic solution to the dynamic wetting on hydrophobic or superhydrophobic surfaces was put forward. Our results will help advance the understanding of dynamic wetting on hydrophobic surfaces at the molecular level and may assist active control of wetting-related applications.

Most of previous molecular dynamics studies used the artificially-assembled lattice structures as the solid phase. In general, the wettability of such virtual lattice structure was tuned by varying the depth of the Lennard-Jones potential with other parameters fixed. This approach would incur inconsistency between the lattice constant and the Lennard-Jones potential and the resultant uncertainties would inextricably interfere the recognition of the physical meanings of MKT parameters. To eliminate such effects in wetting analyses, a real solid-liquid system must be simulated and the confined layer method was used in this work to construct the smooth PTFE surface (the details are given in the Supplementary Information). However, the correctness of the force field parameters for as-formed PTFE must be carefully examined. To this end, the glass transition temperature and the variation of the specific volume of PTFE were compared with the experimental data respectively. The detailed procedures to determine the glass transition temperature *T*
_{g} and the cubic thermal expansion coefficient *α*
_{v} for rubbery and glassy states can be found in the Supplementary Information. As shown in Table 1, our MD results are in good agreement with the experimental values of PTFE. Therefore, thus-formed PTFE, which owns comparable properties with the real PTFE material (Teflon^{®}), was adopted as the substrate surface in our study. The least square fit of the MD results is presented in Fig. 1(a). The glass transition temperature, defined as the intersection of the two straight lines representing the rubbery state and the glassy state, was found to be 116.5 °C, which is in excellent agreement with the reported glass transition temperature of 115 °C for Teflon^{®}
^{36}.

(**a**) Specific volume development of as-formed amorphous PTFE with respect to temperature. (**b**) Distribution of C and F atoms within the amorphous PTFE layer in lateral (*x* and *y*) directions. (**c**) Snapshots (front view) of a water droplet (*d* = 30 nm) spreading **...**

The isotropic properties of a surface would also tremendously affect the dynamic spreading of a liquid on it. When a droplet spreads on an anisotropic surface, such as crystalline materials, its leading edge may be inclined to advance along certain directions, giving rise to irregular contact on the surface. The anisotropic effects of crystalline polymer surfaces on dynamic wetting can be found in a previous study^{37}. The isotropy of the amorphous PTFE surfaces in our study was scrutinized in two ways. Firstly, the lateral distributions of C and F atoms in lateral (*x* and *y*) directions were checked through calculating their number density distributions^{38}. For an amorphous polymer, the orientations of molecular chains are homogenized and the component atoms exhibit an isotropic distribution. From Fig. 1(b), it is apparent that the distributions of C and F atoms in both *x* and *y* directions share the identical profile, evidencing the exceptional isotropy of the as-formed amorphous PTFE layer. Secondly, the droplet shape and profile during wetting were also inspected. The front view snapshots of the morphing water droplet on the amorphous PTFE surface are shown in Fig. 1(c) where the droplet exhibits an approximately symmetric meniscus shape. Obviously, thus-formed PTFE is rather smooth and homogeneous, reflecting the effectiveness of the confined layer method in smoothing out the surface irregularities. Note that there are some water molecules dispersed into the vacuum environment due to evaporation. Figure 1(d) shows the snapshots of the contact area with some irregularities at the perimeter edge due to the inherent fluctuations of the contact line, which will be discussed in detail in the following section. As shown in Fig. 1(c,d), the ideally meniscus-shaped droplet profile and the almost uniform wetting perimeter demonstrate the exceptional isotropy of the amorphous PTFE surface generated by confined layer method.

Water droplets with initial diameters *d*=15 nm, 20nm, 25nm and 30nm were simulated to investigate the possible effects of droplet size. Each case was simulated 5 times with randomly generated Maxwell velocity distribution at 300K and the dynamic quantities including contact angle and droplet base radius were evaluated as the averaged values. To analyze the dynamic wetting behavior of water on hydrophobic PTFE surfaces, MKT is adopted in this work. The basic idea underlying the MKT is that contact line motion is actually a rate process controlled by the corresponding energy barriers among the adsorption sites on a solid surface^{11, 39}. As illustrated in Fig. 2(a), the macroscopic three-phase contact line motion can be viewed as the statistical results of water particles’ adsorption to and detachment from the adsorption sites^{20, 21, 28, 39}. The adsorption and detachment of molecules are controlled by the outward/forward frequency *K*
^{+} and inward/backward frequency *K*
^{-}, respectively:

1

$${K}^{+}=\frac{{k}_{B}T}{h}\mathrm{exp}\left(-\frac{\mathrm{\Delta}G}{{N}_{A}{k}_{B}T}\right)\mathrm{exp}\left(\frac{w{\lambda}^{2}}{2{k}_{B}T}\right)$$

2

$${K}^{-}=\frac{{k}_{B}T}{h}\mathrm{exp}\left(-\frac{\mathrm{\Delta}G}{{N}_{A}{k}_{B}T}\right)\mathrm{exp}\left(-\frac{w{\lambda}^{2}}{2{k}_{B}T}\right)$$

3

4

Note that λ is the unit displacement length defined in MKT and reflects the statistically averaged distance that liquid molecules can move in a single step. At equilibrium state where the dynamic contact angle *θ* reaches the static contact angle *θ*
_{0}, both the adsorption frequency and detachment frequency are equal to a constant value *K*
_{0}:

$${K}^{+}={K}^{-}={K}_{0}=\frac{{k}_{B}T}{h}\mathrm{exp}\left(-\frac{\Delta \mathrm{G}}{{N}_{A}{k}_{B}T}\right)$$

5

The equilibrium state frequency *K*
_{0} is primarily controlled by the total energy barrier *ΔG*. Intrinsically, the MKT assumes that *K*
_{0} is significantly smaller than the frequency associated with the bulk thermal motion. This assumption holds true for most cases since the attractive solid-liquid interactions may restrict the random motion of contact line molecules. Herein, we took the *ΔG* as the total resistance against the contact line and divided it into two main components^{11} as illustrated in Fig. 2(b). One part, *ΔG*
_{w}, is ascribed to the solid-liquid retarding as a result of the work of adhesion between the solid surface and the liquid:

Δ*G*_{w} = *γ*_{LV}(1 + cos *θ*_{0})*N*_{A}*λ*^{2}

6

The other part, Δ*Gv*, arises from the viscous damping between liquid molecules within the contact zone and those in the bulk. Therefore the wetting equation of MKT can be written in its final form as

$${u}_{c}=\frac{2{k}_{B}T\lambda}{h}\mathrm{exp}\left(-\frac{\Delta {G}_{v}}{{N}_{A}{k}_{B}T}\right)\mathrm{exp}\left(-\frac{{\gamma}_{LV}{\lambda}^{2}\left(1+\phantom{\rule{.25em}{0ex}}\mathrm{cos}\phantom{\rule{.10em}{0ex}}{\theta}_{0}\right)}{{k}_{B}T}\right)\mathrm{sinh}\left(\frac{{\gamma}_{LV}{\lambda}^{2}}{2{k}_{B}T}\left(\mathrm{cos}\phantom{\rule{.10em}{0ex}}{\theta}_{0}-\phantom{\rule{.25em}{0ex}}\mathrm{cos}\phantom{\rule{.10em}{0ex}}\theta \right)\right)$$

7

As reflected in equation 7, the wetting velocity *u*_{c} pertains to the dynamic contact angle *θ* via parameters λ, *ΔG*
_{v} and *θ*
_{0}. To determine these parameters, nonlinear regression with 95% confidence interval was used to fit our MD results to the MKT. It was once suggested by some researchers that the relatively high speed range (i.e., the early stage of wetting in this study) should be excluded in applying MKT. However, recent studies^{31, 40} found these suggestions are actually misguiding and the contact line friction is actually dominant Figure 3(a) confirmed these findings since the wetting data obtained from MD simulation shows excellent agreement with the MKT predictions throughout the whole range of dynamic contact angle *θ*. The fitted values of λ, *ΔG*
_{v} and *θ*
_{0} for water droplets with different diameters are shown in Fig. 3(b–d).

(**a**) Comparison of contact line velocity of our MD simulation with the MKT predictions; (**b**) Fitted unit displacement length *λ* and static contact angle *θ*_{0} for different droplet sizes; (**c**) Viscous damping *ΔG*
_{v} and equilibrium frequency *K*
_{0} with respect **...**

Figure 3(b) unveils a weak dependence of static contact angle *θ*
_{0} on the droplet size as *θ*
_{0} of water on PTFE surfaces exhibit a constant value of ~112.4°, which is in line with the experimental value. The landscape of the interatomic potentials on a solid surface can be depicted as an array of potential energy wells (or barriers) as shown in Fig. 2(a). The adsorption sites are located inside each well, where water molecules are inclined to jump into in an advance movement or jump out in a retreat movement. The unit displacement length λ is the statistical result of the distance that liquid molecules within contact zone can move in a single step under the influence of both solid-liquid retarding and viscous damping. Therefore, given the prescribed solid-liquid configuration, λ should remain constant regardless of the droplet size. The independence of λ on water droplet size is confirmed in Fig. 3(b) and λ was found to be ~0.621±0.014 nm. In this study, we suggest that λ be correlated to the backbone configuration of the PTFE chains. Hence we compared the unit displacement length λ to the backbone bond (C-C) length, which is 0.1529 nm, and found that they are on the same order of magnitude and more specifically, λ is ~4 times as long as the backbone bond length.

Figure 3(c) also shows an apparent increase of the molar free energy *ΔG* with the droplet size. From Fig. 2(b), the motion of contact line molecules is resisted by the viscous damping *ΔG*
_{v} and the solid-liquid retarding *ΔG*
_{w}. The solid-liquid retarding refers to the mechanical resistance of the solid surface against the contact line motion in the form of molecular friction. Indeed, *ΔG*
_{w} exhibits a constant value of ~10.35±0.47 kJ/mol, which can be explicitly calculated from equation ^{6}. The viscous damping *ΔG*
_{v} is essentially related to the viscous dissipation within the droplet. In prior studies^{11, 17, 31}, *ΔG*
_{v} is simply attributed to the viscosity as

$$\mathrm{\Delta}{G}_{v}={N}_{A}{k}_{B}T\phantom{\rule{.25em}{0ex}}\mathrm{ln}\phantom{\rule{.10em}{0ex}}\frac{\eta {\nu}_{L}}{h}$$

8

Equation ^{8} infers a constant *ΔG*
_{v} for droplets with different diameters. Nonetheless, *ΔG*
_{v} in Fig. 3(c) shows a systematic positive correlation with the droplet size, indicating *ΔG*
_{v} is not solely related to liquid viscosity *η*. Here we put forward an alternative analysis to clarify the contradiction. *ΔG*
_{v} is regarded as the cohesion forces exerted on the contact line molecules by all other water molecules, either at the vicinity of contact line or in the bulk liquid. In this scenario, *ΔG*
_{v} should be a function of the strength of intermolecular potential and the fraction of contact line molecules

Δ*G*_{v} = *f*(*P*_{ij}, *f*_{c})

9

The meaning of *f*
_{c} is illustrated in Fig. 4. *P*_{ij} is the strength of intermolecular potential and is supposed to be related to the bulk viscosity. Apparently, increasing the intermolecular potential *P*_{ij} is equivalent to enhancing the viscosity and gives rise to a stronger viscous damping. This trend coincides with equation ^{8}. However, a larger *f*
_{c} implies that the motion of a single water molecule at the contact line is damped by fewer water molecules, thus resulting in a smaller *ΔG*
_{v}, so we have

$$\frac{\partial \mathrm{\Delta}{G}_{v}}{\partial {P}_{ij}}>0\phantom{\rule{.25em}{0ex}}\mathrm{and}\phantom{\rule{.25em}{0ex}}\frac{\partial \mathrm{\Delta}{G}_{v}}{\partial {f}_{c}}<0$$

10

For the current study, *P*_{ij} is related to the viscosity *η* and is a constant for TIP4P water model, thus *ΔG*
_{v} should be mainly controlled by *f*
_{c}. As shown in Fig. 4, the contact zone is considered to span *λ* inward and upward from the contact line. The number *n*
_{c} of water molecules within this region can be approximated by

$${n}_{c}=\frac{2\pi \rho R{\lambda}^{2}}{M}$$

11

The total volume of water droplet is scaled as

12

Eventually *f*
_{c} can be calculated from its definition as

$${f}_{c}~\frac{\frac{2\pi \rho R{\lambda}^{2}}{M}}{\frac{\rho {R}^{3}}{M}}=\frac{2\pi {\lambda}^{2}}{{R}^{2}}\phantom{\rule{.25em}{0ex}}$$

13

Apparently a larger droplet has a larger droplet base, giving rise to $\frac{dR}{dd}>0$, therefore

$$\frac{d{f}_{c}}{dd}=\frac{d{f}_{c}}{dR}\frac{dR}{dd}=-\phantom{\rule{-.25em}{0ex}}\frac{2\pi {\lambda}^{2}}{{R}^{3}}\frac{dR}{dd}<0$$

14

Thus, from the chain rule, we have

$$\frac{\partial \mathrm{\Delta}{G}_{v}}{\partial d}=\frac{\partial \mathrm{\Delta}{G}_{v}}{\partial {f}_{c}}\frac{d{f}_{c}}{dd}>0$$

15

Equation ^{15} provides a qualitative prediction for the magnitude of viscous damping and is capable of clarifying the observed increase of viscous damping with respect to droplet diameter in Fig. 3(c). Note that in equation ^{5}, the equilibrium frequency *K*
_{0} is an exponential function of *ΔG*. Thus a slight increase of *ΔG*
_{v} with droplet size may eventually induce a significant decrease of *K*
_{0} as shown in Fig. 3(c).

The singular phenomena in the interfacial flow are usually tackled either in the viscous regime or in the inertial regime^{41, 42}. Recently a pioneering study heralded the existence of a distinct scheme at the onset of drop coalescence^{43}. Similarly, the MKT implies the contact line friction, which stems from the oscillations of contact line molecules from their equilibrium positions, as a nonnegligible mechanism of energy dissipation. From Rayleigh dissipation function^{44}, the friction coefficient *ξ* in the MKT framework can be calculated as:

$$\xi =\frac{{\gamma}_{LV}\left(\mathrm{cos}\phantom{\rule{.10em}{0ex}}{\theta}_{0}-\phantom{\rule{.25em}{0ex}}\mathrm{cos}\phantom{\rule{.10em}{0ex}}\theta \right)}{2{K}_{0}\lambda \phantom{\rule{.25em}{0ex}}\mathrm{sinh}\left(\frac{{\gamma}_{LV}{\lambda}^{2}\left(\mathrm{cos}\phantom{\rule{.10em}{0ex}}{\theta}_{0}-\phantom{\rule{.25em}{0ex}}\mathrm{cos}\phantom{\rule{.10em}{0ex}}\theta \right)}{2{k}_{B}T}\right)}$$

16

Also if *γ*_{LV}*λ*^{2}(cos *θ*_{0} − cos *θ*) ≪ 2*k*_{B}*T*, equations 7 and ^{16} can be linearized with

$${u}_{c}=\frac{{K}_{0}{\gamma}_{LV}{\lambda}^{3}\left(\mathrm{cos}\phantom{\rule{.10em}{0ex}}{\theta}_{0}-\phantom{\rule{.25em}{0ex}}\mathrm{cos}\phantom{\rule{.10em}{0ex}}\theta \right)}{{k}_{B}T}$$

17

and a static friction coefficient *ξ*_{0} can be obtained as:

$${\xi}_{0}=\frac{{k}_{B}T}{{K}_{0}{\lambda}^{3}}$$

18

On one hand, *ξ*
_{0} in equation ^{18} measures the dissipation rate of surface free energy in the form of contact line friction; on the other hand, it can be viewed as the static resistance that contact line experiences at the equilibrium state. Figure 3(d) shows that the static friction coefficient of water on a smooth PTFE surface is on the same order of magnitude with the dynamic viscosity *η*. This parameter is central to many energy analyses in wetting-related phenomena. For example, in the dropwise condensation process, the accurate prediction of droplet growth and motion entails the knowledge of the magnitude of contact line friction^{16}. The self-cleaning of lotus leaves by virtue of roll off of condensate or rain droplets, i.e., the so-called lotus effect^{45}, partially results from the relatively small *ξ*
_{0} on the superhydrophobic surfaces.

According to MKT, water molecules in the contact zone are assumed to move inward and outward with an identical frequency *K*
_{0} at the equilibrium state. Therefore, it can be viewed as a stochastic process in which each water molecule is forced to move or jump 2*K*
_{0}
*τ* times in a period of *τ* and has equal probability to either advance or retract in each step of movement. From this perspective, the motion of water particles at the equilibrium state under the MKT framework can be deemed as a Bernoulli trial *ε*, as shown in Fig. 5(a) and the probability of each outcome is equal to 0.5. Hence the probability *p*(*k*) in an event where water particles move a distance of (2*k*−2*K*
_{0}
*τ*)*λ, k* = 0, 1, 2, … 2*K*
_{0}
*τ*, can be calculated as:

$$p\left(k\right)=\left(\begin{array}{c}2{K}_{0}\tau \\ k\end{array}\right){\left(\frac{1}{2}\right)}^{2{K}_{0}\tau}$$

19

where $\left(\begin{array}{c}2{K}_{0}\tau \\ k\end{array}\right)$ represents the binomial coefficient. Apparently, the mathematical expectation *E*(*ε*) of this Bernoulli trial is equal to 0, which coincides with the fact that the mean wetting velocity is 0 at the equilibrium state. The standard deviation *σ* from the mathematical expectation is given by

$$\sigma =\sqrt{E\left({\epsilon}^{2}\right)-{\left(E\left(\epsilon \right)\right)}^{2}}=\sqrt{\sum _{k=0}^{2{K}_{0}\tau}\left(\begin{array}{c}2{K}_{0}\tau \\ k\end{array}\right){\left(\frac{1}{2}\right)}^{2{K}_{0}\tau}{\left(2k-2{K}_{0}\tau \right)}^{2}{\lambda}^{2}}=\sqrt{2{K}_{0}\tau}\lambda $$

20

Here *σ* characterizes the deviations of water molecules from their equilibrium locations in a period of *τ* and is subsequently transformed into the standard deviation of wetting velocity, *σ*
_{u}, as:

$${\sigma}_{u}=\frac{\sigma}{\tau}=\sqrt{\frac{2{K}_{0}}{\tau}}\lambda $$

21

For a certain system, *τ* represents the characteristic time in the wetting process and can be calculated as the characteristic length *d* divided by the characteristic wetting velocity ${u}^{\u204e}=\frac{{K}_{0}{\gamma}_{LV}{\lambda}^{3}}{{k}_{B}T}$, which is derived from equation ^{17}, so we have

$$\tau =\frac{d}{{u}^{\u204e}}=\frac{{k}_{B}Td}{{K}_{0}{\gamma}_{LV}{\lambda}^{3}}$$

22

Equation ^{22} implies a linear relationship between the time constant *τ* and the droplet size *d*. To characterize the strength of the inherent fluctuations, *σ*_{u} needs to be nondimensionalized and a nondimensional number *Z* is obtained by combining equations ^{21} and ^{22}:

$$Z=\frac{{\sigma}_{u}}{{u}^{\u204e}}=\sqrt{\frac{2{k}_{B}T}{{\gamma}_{LV}\lambda}\frac{1}{d}}=\sqrt{\frac{1.15\times {10}^{-19}}{\lambda}\frac{1}{d}}$$

23

(**a**) MKT at the equilibrium state can be viewed as a Bernoulli trial. (**b**) The fluctuation strength *Z* decays with respect to the characteristic length *d*.

For a certain type of liquid on a given solid surface, the unit displacement length *λ* has been confirmed to be constant. Thus *Z* is supposed to be only affected by the droplet size and the fluctuation of wetting velocity becomes less significant for a larger droplet. More specifically, λ is 0.621 nm for the current case, so we have $Z=\sqrt{\frac{0.185}{d(\mathrm{nm})}}.$ Fig. 5(b) indicates the contact line fluctuations can only be neglected when the characteristic length exceeds micron scale. In a recent study, Wang *et al*.^{18} proved that the oscillation of microscopic contact line suggests a mechanism of enhanced energy dissipation. Subsequently, the prominent contact line fluctuation in the submicron scale is indicative of intensified dissipation of surface free energy in the nanoscale. Such effects need to be accounted for in wetting studies and imply an extra term *F* in modifying the Young’s equation at nanoscale:

24

This equation was also suggested by a more recent study^{32} regarding the contact angle at nanoscale, which ascribes the local contact angle variation at nanoscale to the fluctuations of microscopic force. More specifically, the addition of *F* in calculating contact angle at nanoscale would give rise to the convex nanobending of contact line^{34} for large-scale droplets.

The base radius *R*, defined as the radius of a circle with the same area as the (irregular) contact area of a water droplet, is also investigated in this work. In previous studies^{10, 39}, a scaling law of $R~{t}^{\frac{1}{7}}$ was put forward to describe the development of base radius in spontaneous wetting process. The scaling approach is an effective method especially when rigorous analytical solution is unachievable for complex non-equilibrium phenomena. Nonetheless, the derivation of such a scaling law (Supplementary Information) is based on the assumption of a small contact angle, *i.e*., on a hydrophilic surface. However, the PTFE surface in this study is hydrophobic and the static contact angle *θ*
_{0} is larger than 110°. We thus resorted to an alternative approach to derive a heuristic scaling law of contact line spreading on hydrophobic surfaces. In the dynamic wetting process, the droplet is assumed to follow a meniscus shape, which is given by:

$${R}^{3}=\frac{3V}{\pi}\frac{{\mathrm{sin}}^{3}\theta}{2-3\phantom{\rule{.25em}{0ex}}\mathrm{cos}\phantom{\rule{.10em}{0ex}}\theta +{\mathrm{cos}}^{3}\theta}$$

25

Equations ^{17} and ^{25} can be expanded as the Taylor series at *θ* = *π* due to the obtuse contact angle throughout the wetting process and the simplicity of the Taylor series of the trigonometric functions at *θ* = *π*. Therefore, equations ^{17} and ^{25} can be rewritten as

$$\frac{dR}{dt}={u}_{c}~{C}_{1}-{\left(\theta -\pi \right)}^{2}$$

26

$$R\approx -\sqrt[3]{\frac{3V}{4\pi}}\left(\theta -\pi \right)~-\phantom{\rule{-.25em}{0ex}}{C}_{2}(\theta -\pi )$$

27

where *C*
_{1} and *C*
_{2} are constants and larger than 0. By combining equations ^{26} and ^{27}, the base radius development can be summarized as

$$R\approx \frac{R(t=0)}{\mathrm{tanh}\phantom{\rule{.10em}{0ex}}a}\phantom{\rule{.25em}{0ex}}\mathrm{tanh}\left(a+\frac{t}{\phantom{\rule{.25em}{0ex}}\tau}\right)$$

28

where *a* and *τ* are the coefficients obtained by fitting the MD results and *τ* can be taken as the characteristic time in the wetting process. Figure 6(a) shows the comparison of our MD data with equation ^{28}. The proposed *R* ~ tanh(*t*) relationship demonstrates excellent agreement with the MD results. Also as shown in Fig. 6(b), the linear relationship between the time constant τ and droplet size *d* is in accord with equation ^{22}.

(**a**) Base radius development of a water droplet on a hydrophobic PTFE surface with time. (**b**) Linear relationship between time constant *τ* and droplet diameter *d*.

The contact angle in the studied wetting process ranges from 160° to 110° with a divergence from the truncation angle (*θ* = 180°) in the Taylor expansion, from which the slight deviation of equation ^{28} from the MD data may origin. However, it can be expected that the *R*~ tanh(*t*) relationship would yield a more accurate prediction of the base radius evolution on a superhydrophobic surface with *θ* > 150°.

In this study, the dynamic wetting of water on hydrophobic surfaces is explored using large-scale molecular dynamics simulations. The dynamic wetting process is elucidated in the molecular kinetic framework and the molecular kinetic parameters are interpreted from an *ab initio* perspective. Our results indicate that the contact line friction originates from both the solid-liquid retarding and viscous damping. The solid-liquid retarding arises from the work of adhesion, while the viscous damping cannot be simply related to the viscosity. The contact line friction coefficient of water on smooth PTFE surfaces is on the same order of magnitude as the dynamic viscosity of water, indicating the contact line friction cannot be neglected in the analysis of processes involving moving contact line, such as condensate droplet growth and droplet coalescence. Based on these observations, we formulated a nondimensional number to characterize the inherent contact line fluctuations that are responsible for the nanobending of contact line. Finally, an asymptotic solution was derived to describe the moving contact line on a hydrophobic surface. It is believed that these findings can aid wetting-related studies and can be applied to various engineering applications.

All the simulations in this study were carried out on the single-precision molecular dynamics package Gromacs 5.1.2^{46}. The total potential energy between two atoms *i* and *j* separated by *r*
_{ij} is the sum of Lennard-Jones potential and the Coulombic pairwise potential:

$${E}_{ij}=4{\epsilon}_{ij}\left[{\left(\frac{{\sigma}_{ij}}{{r}_{ij}}\right)}^{12}-{\left(\frac{{\sigma}_{ij}}{{r}_{ij}}\right)}^{6}\right]+\frac{{k}_{e}{q}_{i}{q}_{j}}{{r}_{ij}}$$

29

where *ε*_{ij} is the depth of the potential well, *σ* is the zero-crossing distance, *k*_{e} is the Coulomb constant and *q* is the charge of each atom. The time step in each simulation was chosen to be 1 femtosecond (fs). The short range forces were cut-off at 1.2 nm with the Verlet scheme. The particle mesh Ewald method was adopted to handle the long range Coulombic interactions. LINear Constraint Solver (LINCS) algorithm^{46} were incorporated to enforce all the bonding constraints applied on each atom. A constant temperature of 300 K was maintained in both NPT (constant number *N*, pressure *P* and temperature *T*) and NVT (constant number *N*, volume *V* and temperature *T*) ensembles through a velocity-rescaling scheme.

To construct the amorphous PTFE surface, a single chain molecule consisting of 30 monomers was firstly created, which gives the exact chemical formula of CF_{3}-(CF_{2})_{28}-CF_{3}. The force field parameters are based on the Optimized Potentials for Liquid Simulations-all atom (OPLSAA) force filed^{47, 48}, which has shown excellent capability in simulating dynamic wetting polyethylene and polyvinyl chloride^{49}. The detailed preparation of PTFE surfaces can be found in the Supplementary Information.

The initial structures of water droplets with diameters of *d* = 15 nm, 20 nm, 25 nm and 30 nm were prepared in the open source software Packmol^{50}. The amount of the water molecules that are encapsulated into the prescribed spherical space should be consistent with the water density of 998.2 kg/m^{3}, for example, 270930 water molecules for *d* = 25 nm. The TIP4P^{51} water model was used to account for the topology information of water molecules, since it shows great agreement with the liquid water properties^{52}. The water sphere was then placed 0.2 nm away from the confined PTFE surface to allow interactions between the solid and the liquid atoms to take effect. The dynamic wetting process was then simulated with temperature maintained at 300 K through an NVT ensemble. The measurement of the instantaneous contact angle followed the algorithm reported by previous MD studies^{53, 54}.

This work is financially supported by NSF CBET under grant number 1550299 and Institute for Critical Technology and Applied Science (ICTAS) at Virginia Tech. The authors also acknowledge Advanced Research Computing at Virginia Tech (http://www.arc.vt.edu) for providing computational resources and technical support that have contributed to the results reported within this paper.

Author Contributions

J.C. and L.Z. designed the project. L.Z. carried out the simulation and data analysis. J.C. and L.Z. prepared and reviewed the manuscript.

The authors declare that they have no competing interests.

**Electronic supplementary material**

**Supplementary information** accompanies this paper at doi:10.1038/s41598-017-11350-6

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

1. Giordano N, Cheng JT. Microfluid mechanics: progress and opportunities. J Phys-Condens Mat. 2001;13:R271–R295. doi: 10.1088/0953-8984/13/15/201. [Cross Ref]

2. Miljkovic N, Enright R, Wang EN. Effect of Droplet Morphology on GrowthDynamics and Heat Transfer during Condensation on Superhydrophobic NanostructuredSurfaces. Acs Nano. 2012;6:1776–1785. doi: 10.1021/nn205052a. [PubMed] [Cross Ref]

3. Cheng, J. T., Vandadi, A. & Chen, C. L. Condensation Heat Transfer on Two-Tier Superhydrophobic Surfaces. Proceedings of the Asme International Mechanical Engineering Congress and Exposition–2012, *Vol. 7*, *Pts a-D*, 2649–2653 (2013).

4. Jensen KE, et al. Wetting and phase separation in soft adhesion. P Natl Acad Sci USA. 2015;112:14490–14494. doi: 10.1073/pnas.1514378112. [PubMed] [Cross Ref]

5. Mugele F, Baret J-C. Electrowetting: from basics to applications. Journal of Physics:Condensed Matter. 2005;17:R705–R774.

6. Zhao Y-P, Wang Y. Fundamentals and Applications of Electrowetting. Reviews of Adhesion and Adhesives. 2013;1:114–174. doi: 10.7569/RAA.2013.097304. [Cross Ref]

7. Chow TS. Wetting of rough surfaces. J Phys-Condens Mat. 1998;10:L445–L451. doi: 10.1088/0953-8984/10/27/001. [Cross Ref]

8. Quere D. Wetting and roughness. Annu Rev Mater Res. 2008;38:71–99. doi: 10.1146/annurev.matsci.38.060407.132434. [Cross Ref]

9. Grewal, H. S., Kim, H. N., Cho, I. J. & Yoon, E. S. Role of Viscous Dissipative Processes on the Wetting of Textured Surfaces. Sci Rep-Uk **5** (2015). [PMC free article] [PubMed]

10. Blake TD. The physics of moving wetting lines. J Colloid Interface Sci. 2006;299:1–13. doi: 10.1016/j.jcis.2006.03.051. [PubMed] [Cross Ref]

11. Snoeijer JH, Andreotti B. Moving Contact Lines: Scales, Regimes, and Dynamical Transitions. Annu Rev Fluid Mech. 2013;45:269–292. doi: 10.1146/annurev-fluid-011212-140734. [Cross Ref]

12. Ralston J, Popescu M, Sedev R. Dynamics of wetting from an experimental point of view. Annu Rev Mater Res. 2008;38:23–43. doi: 10.1146/annurev.matsci.38.060407.130231. [Cross Ref]

13. Voinov OV. Hydrodynamics of wetting. Fluid dynamics. 1977;11:714–721. doi: 10.1007/BF01012963. [Cross Ref]

14. Cox RG. The Dynamics of the Spreading of Liquids on a Solid-Surface .1. Viscous-Flow. J Fluid Mech. 1986;168:169–194. doi: 10.1017/S0022112086000332. [Cross Ref]

15. Li H, Paneru M, Sedev R, Ralston J. Dynamic electrowetting and dewetting of ionic liquids at a hydrophobic solid-liquid interface. Langmuir. 2013;29:2631–2639. doi: 10.1021/la304088t. [PubMed] [Cross Ref]

16. Gokhale SJ, Plawsky JL, Wayner PC. Experimental investigation of contact angle, curvature, and contact line motion in dropwise condensation and evaporation. J Colloid Interf Sci. 2003;259:354–366. doi: 10.1016/S0021-9797(02)00213-8. [PubMed] [Cross Ref]

17. Blake TD, De Coninck J. Dynamics of wetting and Kramers’ theory. Eur Phys J-Spec Top. 2011;197:249–264. doi: 10.1140/epjst/e2011-01467-2. [Cross Ref]

18. Wang, J. Y. *et al*. Surface structure determines dynamic wetting. *Sci Rep-Uk***5** (2015). [PMC free article] [PubMed]

19. Wang F-C, Zhao Y-P. Slip boundary conditions based on molecular kinetic theory: The critical shear stress and the energy dissipation at the liquid–solid interface. Soft Matter. 2011;7:8628. doi: 10.1039/c1sm05543g. [Cross Ref]

20. Blake, T. D. & Haynes, J. M. Kinetics of Liquid/Liquid Displacement. *J Colloid Interf Sci***30**, 421-& (1969).

21. de Ruijter MJ, De Coninck J, Blake TD, Clarke A, Rankin A. Contact angle relaxation during the spreading of partially wetting drops. Langmuir. 1997;13:7293–7298. doi: 10.1021/la970825v. [Cross Ref]

22. Fetzer R, Ralston J. Dynamic Dewetting Regimes Explored. J Phys Chem C. 2009;113:8888–8894. doi: 10.1021/jp901719d. [Cross Ref]

23. Ray S, Sedev R, Priest C, Ralston J. Influence of the Work of Adhesion on the Dynamic Wetting of Chemically Heterogeneous Surfaces. Langmuir. 2008;24:13007–13012. doi: 10.1021/la802264d. [PubMed] [Cross Ref]

24. Blake TD. Forced wetting of a reactive surface. Adv Colloid Interfac. 2012;179:22–28. doi: 10.1016/j.cis.2012.06.002. [PubMed] [Cross Ref]

25. Seveno D, Blake TD, Goossens S, De Coninck J. Predicting the Wetting Dynamics of a Two-Liquid System. Langmuir. 2011;27:14958–14967. doi: 10.1021/la2034998. [PubMed] [Cross Ref]

26. Ramiasa M, Ralston J, Fetzer R, Sedev R. Contact Line Friction in Liquid–Liquid Displacement on Hydrophobic Surfaces. The Journal of Physical Chemistry C. 2011;115:24975–24986. doi: 10.1021/jp209140a. [Cross Ref]

27. Li H, Sedev R, Ralston J. Dynamic wetting of a fluoropolymer surface by ionic liquids. Physical Chemistry Chemical Physics. 2011;13:3952–3959. doi: 10.1039/c0cp02035d. [PubMed] [Cross Ref]

28. Blake TD, Clarke A, DeConinck J, deRuijter MJ. Contact angle relaxation during droplet spreading: Comparison between molecular kinetic theory and molecular dynamics. Langmuir. 1997;13:2164–2166. doi: 10.1021/la962004g. [Cross Ref]

29. Wang F-C, Zhao Y-P. Contact angle hysteresis at the nanoscale: a molecular dynamics simulation study. Colloid and Polymer Science. 2012;291:307–315. doi: 10.1007/s00396-012-2747-2. [Cross Ref]

30. De Coninck J, Blake TD. Wetting and molecular dynamics simulations of simple liquids. Annu Rev Mater Res. 2008;38:1–22. doi: 10.1146/annurev.matsci.38.060407.130339. [Cross Ref]

31. Sedev R. The molecular-kinetic approach to wetting dynamics: Achievements and limitations. Adv Colloid Interfac. 2015;222:661–669. doi: 10.1016/j.cis.2014.09.008. [PubMed] [Cross Ref]

32. Lukyanov AV, Likhtman AE. Dynamic Contact Angle at the Nanoscale: A Unified View. Acs Nano. 2016;10:6045–6053. doi: 10.1021/acsnano.6b01630. [PubMed] [Cross Ref]

33. Degennes PG, Hua X, Levinson P. Dynamics of Wetting - Local Contact Angles. J Fluid Mech. 1990;212:55–63. doi: 10.1017/S0022112090001859. [Cross Ref]

34. Chen L, Yu JP, Wang H. Convex Nanobending at a Moving Contact Line: The Missing Mesoscopic Link in Dynamic Wetting. Acs Nano. 2014;8:11493–11498. doi: 10.1021/nn5046486. [PubMed] [Cross Ref]

35. Deng YJ, Chen L, Yu JP, Wang H. Nanoscopic morphology of equilibrium thin water film near the contact line. Int J Heat Mass Tran. 2015;91:1114–1118. doi: 10.1016/j.ijheatmasstransfer.2015.08.057. [Cross Ref]

36. Teflon, P. & Resin, P. F. Properties handbook. *DuPont Fluoroproducts, Washington* (1996).

37. Fan CF, Cağin T. Wetting of crystalline polymer surfaces: A molecular dynamics simulation. The Journal of Chemical Physics. 1995;103:9053. doi: 10.1063/1.470016. [Cross Ref]

38. Giorgino T. Computing 1-D atomic densities in macromolecular simulations: The density profile tool for VMD. Comput Phys Commun. 2014;185:317–322. doi: 10.1016/j.cpc.2013.08.022. [Cross Ref]

39. de Ruijter MJ, Blake TD, De Coninck J. Dynamic wetting studied by molecular modeling simulations of droplet spreading. Langmuir. 1999;15:7836–7847. doi: 10.1021/la990171l. [Cross Ref]

40. Paneru M, Priest C, Sedev R, Ralston J. Static and Dynamic Electrowetting of an Ionic Liquid in a Solid/Liquid/Liquid System. J Am Chem Soc. 2010;132:8301–8308. doi: 10.1021/ja9106397. [PubMed] [Cross Ref]

41. Thete, S. S., Anthony, C., Basaran, O. A. & Doshi, P. Self-similar rupture of thin free films of power-law fluids. *Phys Rev E***92** (2015). [PubMed]

42. Castrejon-Pita JR, et al. Plethora of transitions during breakup of liquid filaments. P Natl Acad Sci USA. 2015;112:4582–4587. doi: 10.1073/pnas.1418541112. [PubMed] [Cross Ref]

43. Paulsen JD, et al. The inexorable resistance of inertia determines the initial regime of drop coalescence. P Natl Acad Sci USA. 2012;109:6857–6861. doi: 10.1073/pnas.1120775109. [PubMed] [Cross Ref]

44. Goldstein, H. *Classical mechanics*. (Pearson Education India, 1965).

45. Cheng J, Vandadi A, Chen C-L. Condensation heat transfer on two-tier superhydrophobic surfaces. Appl Phys Lett. 2012;101:131909. doi: 10.1063/1.4756800. [Cross Ref]

46. Bekker, H. *et al*. Gromacs: A parallel computer for molecular dynamics simulations. *Physics computing***92**, 252–256 (1993).

47. Jorgensen WL, Maxwell DS, TiradoRives J. Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids. J Am Chem Soc. 1996;118:11225–11236. doi: 10.1021/ja9621760. [Cross Ref]

48. Kaminski GA, Friesner RA, Tirado-Rives J, Jorgensen WL. Evaluation and reparametrization of the OPLS-AA force field for proteins via comparison with accurate quantum chemical calculations on peptides. J Phys Chem B. 2001;105:6474–6487. doi: 10.1021/jp003919d. [Cross Ref]

49. Hirvi, J. T. & Pakkanen, T. A. Molecular dynamics simulations of water droplets on polymer surfaces. Journal of Chemical Physics **125** (2006). [PubMed]

50. Martinez L, Andrade R, Birgin EG, Martinez JM. PACKMOL: A Package for Building Initial Configurations for Molecular Dynamics Simulations. J Comput Chem. 2009;30:2157–2164. doi: 10.1002/jcc.21224. [PubMed] [Cross Ref]

51. Jorgensen, W. L., Chandrasekhar, J., Madura, J. D., Impey, R. W. & Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. *Journal of Chemical Physics***79**, 926–935 (1983).

52. Vega C, De ME. Surface tension of the most popular models of water by using the test-area simulation method. The Journal of chemical physics. 2007;126:154707. doi: 10.1063/1.2715577. [PubMed] [Cross Ref]

53. Webb EB, Grest GS, Heine DR, Hoyt JJ. Dissolutive wetting of Ag on Cu: A molecular dynamics simulation study. Acta Mater. 2005;53:3163–3177. doi: 10.1016/j.actamat.2005.03.021. [Cross Ref]

54. Sun, Y. & Webb, E. B. The atomistic mechanism of high temperature contact line advancement: results from molecular dynamics simulations. *J Phys-Condens Mat***21** (2009).

Articles from Scientific Reports are provided here courtesy of **Nature Publishing Group**

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's national science library in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |