Invited Feature Article pubs.acs.org/Langmuir

Molecular Dynamics Simulations: Insight into Molecular Phenomena at Interfaces Sepideh Razavi,† Joel Koplik,*,‡ and Ilona Kretzschmar*,† †

Department of Chemical Engineering and ‡Department of Physics and The Benjamin Levich Institute, City College of the City University of New York, New York, New York 10031, United States ABSTRACT: Molecular dynamics simulations, when aptly devised, can enhance our fundamental understanding of a system, set up a platform for testing theoretical predictions, and provide insight and a framework for further experimental studies. This feature article highlights the importance of molecular dynamics simulations in understanding interfacial phenomena using three case studies involving liquid−liquid and solid−liquid interfaces. After briefly reviewing molecular dynamics methods, we discuss velocity slip at a liquid−liquid interface, the coalescence of liquid drops in suspension and in free space, and the behavior of colloidal nanoparticles at a liquid−liquid interface. We emphasize the utility of simple intermolecular potentials and generic liquids. The case studies exemplify the significant insight gained through the molecular modeling approach regarding the interfacial phenomena studied. We conclude the highlight with a brief discussion illustrating potential shortcomings and pitfalls of molecular dynamics simulations.

1. INTRODUCTION Molecular dynamics (MD) methods play an important role in studying many-body problems by examining the motion of individual particles using a molecular-scale model of the system where all aspects of the behavior are open to inspection.1 The prominent feature of MD simulations is their access to the subtle details of molecular motion and structure where the speed of molecular events makes it difficult to attain this information experimentally.2 In light of experimental limitations, proper MD simulations can be devised to serve a twofold purpose: (i) an exploratory approach, where computer visualizations can help to relate obscure macroscopic behavior to a specific microscopic phenomenon, and (ii) a predictive approach, where a physical quantity can be consistently predicted for a potential future measurement.3 In view of these capabilities, MD simulations can apprise diverse research areas by complementing and enriching the traditional approaches of theory and experiment.4−7 One of the fields that benefits remarkably from the capacity of MD simulations is interface and colloid science. Interface and colloid science is a prominent interdisciplinary field with numerous applications in nanotechnology, the chemical industry, biotechnology, medical science, and many others.8 Very often, the behavior of these systems is dominated by the interfaces involved. Understanding some underlying aspects of these systems requires information on both the spatial organization and temporal dynamics of the interfacial boundaries that occur on the scale of molecules.9 Hence, MD simulations are well-suited and have proven to be an invaluable tool for studying system behavior across vast spatiotemporal domain−length scales of up to thousands of angstroms, with © XXXX American Chemical Society

atomic precision, and time scales of up to milliseconds, at femtosecond resolution.2 The aim of this feature article is to elucidate the unique contributions that MD simulations provide in the framework of interface and colloid science with a focus on the research activities of the authors. In this regard, we employ some of our recent work on a selected number of interesting problems that involve less heavily studied interfacial phenomena. In section 2, three representative interfacial problems and their importance in the field of colloid and interface science are introduced. These problems involve liquid−liquid and solid−liquid systems and a combination of both interfaces. The interested reader is directed to refs 10 and 11 for earlier comprehensive reviews of MD simulations treating phenomena at the vapor−liquid interface and examples of more recent work in refs 12−14. Section 3 briefly describes the fundamentals of MD simulations and the general procedure used to run the calculations. Next, in section 4, the simulation technique is applied to the three problems, and the specific findings accessed through the MD simulations are emphasized for each case. Some final comments on potential shortcomings and precautions when running MD simulations are given in section 5.

2. INTERFACIAL CASE STUDIES Interface and colloid science revolves around the study of vapor−liquid, liquid−liquid, and solid−liquid interfaces. Three representative examples are centered around the following Received: January 27, 2014 Revised: March 28, 2014

A

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

removing adjoining interfaces closer than some microscopic cutoff, and rupture is assumed to take place when a film thickness falls below some suitable value. In reality, on the molecular level a fluid interface has a finite thickness and may rupture when its molecules divide into spatially separated groups, and two interfaces can coalesce when their respective molecules occupy the same spatial region. MD simulations provide an ideal tool for studying interfacial rearrangements. Hence, the second problem described in section 4 examines the small-scale dynamics of droplet coalescence in which two liquid bodies merge and probes the importance of molecular interfacial rearrangements.26 The first two problems introduced here and studied in sections 4.1 and 4.2 solely involve liquid−liquid interfaces. Section 4.3 considers both liquid−liquid and solid−liquid interfaces in systems where colloidal particles are present at liquid−liquid interfaces.27,28 Colloidal particles of the proper wettability bind strongly to interfaces. The large binding energy, approximated via a continuum thermodynamic model,29 can be millions of times the thermal energy kBT for a micrometer-sized particle and is attributed to the reduction of the energetically unfavorable liquid−liquid interface. Hence, liquid−liquid interfaces can serve as promising templates for novel applications such as self-assembly and biphasic catalytic reactions.30,31 To develop these applications, it is important to advance our understanding of colloidal particles in the vicinity of interfaces. This is an emergent field with many open questions regarding the dynamics of binding,32 the impact of particle amphiphilicity and adsorption kinetics on the binding mechanism,27 the rotational and translational diffusion of particles at interfaces compared to their behavior in the bulk,33 and the influence of surface nanofeatures on the binding process.28 Despite some recently developed experimental techniques, it remains difficult to investigate these problems and probe all of the underlying characteristics of the system solely through experiments. For example, very recently the dynamics of a micrometer-sized colloidal particle binding to an oil−water interface was captured using digital holography microscopy.32 It was found that colloidal particles do not reach their equilibrium contact angle, described by Young’s equation, and it was speculated that the presence of surface heterogeneities is the reason for this behavior. Nevertheless, it was unclear whether the defects are the surface charge groups or some other features, such as asperities. Moreover, because of optical constraints, this technique can be applied only to a specific range of particle size and type of material. Thus, simulation methods can be regarded as a tool to elevate our understanding of such systems by identifying general features of energetics, dynamics and kinetics of the binding process and assessing the factors that affect interfacial adsorption. In all of the interfacial case studies presented here, information regarding interfacial boundaries is critical and simulations are needed to understand and inform experiments and provide insight into theoretical predictions.

three questions: Do two liquids experience slip at a liquid− liquid interface? How do two liquid interfaces coalesce? How does a solid particle interact with a liquid−liquid interface? The question of slip at a solid−liquid interface, whether the liquid’s velocity agrees with that of a bounding solid, has attracted enormous interest in the past decade15,16 because of its relevance to micro- and nanofluidic devices and flows. Because this boundary condition requires additional information beyond the equations of motion themselves, it is ostensibly a matter to be decided by experiment, but because measurements on the appropriate atomic scale are difficult to make, much of the relevant information has come from MD simulations.17 More recently, experiments involving newer techniques such as the surface force apparatus18,19 have begun to compete in resolution with the simulations, considerably enriching the information available. Numerous examples of slip have now been observed and quantified by means of a Navier boundary condition. Typically, this slip is associated with a weak interaction between solid and fluid molecules, particularly when the fluid is itself a gas or when surface irregularities trap bubbles or films on the surface. This phenomenology raises the question of slip at a liquid−liquid interface: normally it is assumed20 that the vector velocity is continuous across a fluid interface and that stress is continuous except for the surface tension and Marangoni effects. Underlying this continuity condition is the implicit assumption that molecules of the two fluids mix in the interfacial region and their mutual interactions would lead to a single velocity there. However, realistic liquid− liquid molecular interactions exhibit a similar diversity to the solid−liquid case, and one might question whether slip occurs when the attractive interaction between the molecules of the respective liquids are so weak that mixing is poor. Experiments at a liquid−liquid interface are difficult to conduct because thermal and hydrodynamic fluctuations occur there and the interface does not need to be a fixed planar surface. However, in the case of polymeric liquids there is both indirect and observational evidence that such slip occurs;21−23 furthermore, at low shear stress the velocity discontinuity is proportional to the stress, consistent with the Navier boundary condition widely used to describe solid−liquid slip. Moreover, MD simulations employing coarse-grained molecular models of polymer melts by Barsky and Robbins24 are consistent with the observations and the existence of slip at the interface between two melts. This liquid−liquid slip in polymers is often attributed to the incomplete entanglement of molecular chains of the different species, which suggests that weak attraction may play a role and that ordinary liquids may exhibit slip as well. With this motivation, as a first case study in section 4, we have investigated the possibility of having slip at liquid−liquid interfaces and have questioned the validity of the no-slip boundary condition at such interfaces.25 In the previous case, the interactions between two immiscible liquids were discussed. The second case study revolves around the interaction of two liquid drops initially separated by a surrounding immiscible liquid that are subsequently driven toward each other until they coalesce. In continuum fluid mechanics, an interface between fluids is a mathematical surface where boundary conditions are applied and surface tension exerts a stress. There is no natural provision in the Navier− Stokes equation for an interface to disappear or for a new interface to form, but in nature processes such as drop coalescence or breakup are quite common. In practice, changes in topology are inserted by hand: coalescence is modeled by

3. SIMULATION METHOD: AN OVERVIEW The vast diversity in the physical properties of materials results from the differences in the interaction potential energies of their atomic constituents, but the two essential common features are an attraction due to covalent or ionic bonds or dispersion forces at a separation of a few atomic diameters and a strong repulsion at short distances when electrons approach closely. To the extent that one is interested in generic B

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

Figure 1. (a) Lennard-Jones potential for cij = 1 in eq 1. (b) Phase diagram of the Lennard-Jones system. The numbers in the legend give the length of the simulation box in units of σ, and Exp. refers to experimental values. Adapted from ref 36.

fluctuations in position.40 Aside from simple liquid systems, the dispersion interactions encapsulated in LJ potentials can be critical in other situations (i.e., colloidal nanoparticles), where they form part of the DLVO interaction41 and may dominate when the electrostatic component is weak.30 In this case, the nanoparticle is constructed by carving a sphere with the desired radius out of an atomic lattice, and the motion of the nanoparticle is computed by rigid body dynamics, where all of the atoms that constitute the particle are fixed at their respective lattice sites. The mechanics of molecular dynamics (MD) simulations are well described in the monographs of Allen and Tildesley34 and Frenkel and Smit,42 for example. Typical simulations of flow, phase structure, mixing, and so forth operate in the classical regime, where gradients of the potential energy provide the forces on the atoms. Newton’s equations are then a coupled set of ordinary differential equations, which may be solved by predictor−corrector or Verlet-type algorithms, the former best for large forcing or rapid atomic rearrangements and the latter preferred when long-time stability is an issue. One important feature in many simulations is the use of periodic boundary conditions when it is desirable to avoid unwanted spatial heterogeneity. In this case, it must be noted that the size of the periodic box limits the wavelength of possible hydrodynamic instabilities and that the periodic images do interact with the system of interest. In most situations, a thermostat is used to control the temperature of the system, which may drift in equilibrium or rise because of viscous heating when flow occurs. Although numerous methods are available to adjust the temperature (usually proportional to the kinetic energy), the preferred methods produce the appropriate statistical ensemble in equilibrium and do not bias the velocity field when flow occurs (or, equivalently, conserve momentum). One method is to couple the atoms to Langevin dynamics, a fluctuating force plus frictional drag, in the neutral plane orthogonal to the flow if there is one, another is Nosé−Hoover dynamics, which ultimately involves extra forces proportional to the velocity fluctuation about a local average (a “profile-unbiased” thermostat),43 and a third recent alternative is a DPD thermostat that involves Langevin-like pairwise interactions.44 A key question in formulating MD simulations is the choice of the size of the system. A lower bound arises from the need to have enough atoms to obtain continuum behavior on average; a

properties rather than any specific material, it suffices to consider simple, uncharged, spherically symmetric atoms such as argon where the appropriate two-body potential energy is of the Lennard-Jones (LJ) form34,35 ⎡⎛ σ ⎞12 ⎛ σ ⎞6 ⎤ Vij = 4ε⎢⎜ ⎟ − cij⎜ ⎟ ⎥ ⎝r⎠ ⎦ ⎣⎝ r ⎠

(1)

Here, r is the distance between atomic centers, ε is the depth of the potential well, and σ is the approximate diameter of an atom, the separation at which the short-distance repulsion begins to increase. The constants cij are a generalization used below, which permits one to adjust the strength of attraction between any two coexisting atomic species i and j. Physically, the functional form of the attractive r−6 term originates in calculations of electric charge fluctuations, and the repulsive r−12 term is just a convenient parametrization. Because the potential decays rapidly as r increases, it is common to ignore it beyond a cutoff distance, commonly chosen as rc = 2.5σ. For argon, the parameter values that fit the experiment are ε/kB ≈ 120 K and σ ≈ 0.32 nm, and these, along with the mass of argon m = 40 a.u., provide the natural atomic time scale of τ = m(σ/ε)1/2 ≈ 2.2 ps. The basic potential (with all cij = 1) and the resulting phase diagram,36 which exhibits all of the relevant thermodynamic phases, are demonstrated in Figure 1. Even in generic material systems, additional interactions are often needed. For example, a monatomic LJ system has a high vapor pressure and a very diffuse vapor−liquid interface and may not be optimal for studying problems in drop dynamics. A simple remedy is to construct short linear chain molecules by binding adjacent atomic pairs using a harmonic or FENE potential,37 VH =

⎛ 1 2 1 r2 ⎞ kr VF = − kr0 2 log⎜1 − 2 ⎟ 2 2 r0 ⎠ ⎝

(2)

One finds that the resulting FENE dimer molecules have a vapor pressure similar to that of air at atmospheric pressure and room temperature, whereas tetrameric liquids are accompanied by only isolated vapor molecules and longer chains are not volatile.38,39 A similar device is commonly used to study coexisting LJ liquids and solids: solid atoms are tethered to fixed solid lattice sites using a harmonic potential, with the spring stiffness k chosen to produce the appropriate thermal C

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

common rule of thumb is that at least 10 constituent objects (atoms or molecules) per direction of space per chemical species are needed. The upper bound in principle is obtained as the size beyond which the quantity of interest does not vary with the simulation size, although in practice hardware considerations or accessibility and budgetary limitations dictate this. In the case of interacting vapor−liquid systems, an additional consideration is phase coexistence. The relative amounts of liquid and vapor in equilibrium are determined by the phase diagram in Figure 1b. So, for example, the number of atoms and the size of the simulation box needed to study a drop of a certain radius in equilibrium with vapor at some temperature can be anticipated. In the following sections, we discuss several simulations that are meant to mimic laboratory experiments: a system is first allowed to reach equilibrium in an uninteresting noninteracting state, and then a constraint is released and a new equilibrium is attained or a force is imposed to drive the system toward a new state.

4. APPLYING MD SIMULATIONS TO INTERFACIAL PROBLEMS 4.1. Liquid−Liquid Interfacial Slip. Within the context of LJ interactions, slip at a solid−liquid interface increases as the strength of the attractive term in eq 1 decreases.45 We have used MD simulations to test for the presence of liquid−liquid interfacial slip in simple liquids and to clarify the resulting boundary condition.25 We considered three pairs of liquids: (i) two monomers, (ii) a dimer and a tetramer chain (i.e., a (2, 4) system), and (iii) linear chains of length 4 and 16, with interactions corresponding to either complete immiscibility at the interface (c11 = c22 = 1 and c12 = 0 in eq 1) or partial interfacial miscibility following the Lorentz−Berthelot (LB) rule (c11 = 3/4, c22 = 5/4, and c12 = (c11c22)1/2), in both Couette and Poiseuille flows. The importance of simulating two different flows is that a genuine boundary condition on a partial differential (Navier−Stokes) equation should involve a material parameter that is dependent on the fluids involved but independent of the flow configuration. The simulations involved about 10 000 atoms arranged into contiguous slabs of the two liquids confined between flat and parallel atomistic walls with periodic boundary conditions in the directions parallel to the walls. Couette flow is generated by translating the walls, and Poiseuille flow is generated by a body force applied to each fluid atom. The results were analyzed by dividing the region between the walls into parallel slab-shaped sampling bins and recording the density, velocity, and stress as a function of the distance from one wall. The simulated Reynolds number based on the channel width ranged from 0.01 to 1.0, roughly in the creeping flow regime. Likewise, the fluid itself is in the Newtonian regime in these simulations since the Deborah number De = γ̇τ (the ratio of a characteristic molecular time of the fluid, the LJ time scale τ, to the characteristic time associated with the imposed strain, 1/γ̇) ranged from 0.001 to 0.1. The various fluid systems exhibit the same qualitative behavior, which is illustrated by the results for the dimer− tetramer system. The key issue is the degree of mixing between the two liquids, and the effect of miscible versus immiscible interactions can be seen visually by comparing snapshots of the fluid configuration in the two cases and quantitatively by comparing the density profiles. When there is an attraction between the two liquids, we see in Figure 2a that the two types

Figure 2. (Top panel) Interfacial region for the (2, 4) system with (a) miscible and (b) immiscible interactions. Frames a and b show 3D slabs centered on the interface, as viewed from a long distance along the interface, for the case of miscible and immiscible interactions, respectively. Molecules are represented by the line segments joining the atomic centers using red/blue lines for dimer/tetramer molecules. (Bottom panel) Couette flow density profiles as a function of coordinate y running normal to the interface, averaged over the other two directions and a 1000τ time interval for (c) miscible and (d) immiscible interactions. The (×) symbol refers to dimers and (*) refers to tetramers, whereas the continuous curve is the total liquid density. Reference 25. Copyright 2006 by the American Physical Society.

of molecules are well mixed in the center of the channel and in Figure 2c that the corresponding number density profile varies monotonically from one bulk value to the other. In contrast, when the liquids are immiscible, there is an obvious gap between their respective molecules in the snapshot (Figure 2b) accompanied by a distinct dip in the density profile (Figure 2d). The velocity profile, displayed in Figure 3a, is flat in the two bulk fluid regions for the miscible case (○), with different slopes corresponding to the different fluid viscosities, and has a smooth transition extending across the width of the interfacial region where the slope varies gradually. In the immiscible case (+), although the velocity still varies continuously it has a very rapid variation over a shorter distance, roughly the width of the density dip. From the point of view of a continuum description or that of a realistic experiment without angstrom resolution, one would describe the immiscible velocity profile as discontinuous with an apparent slip of about 0.15σ/τ. In both cases, the shear stress is constant across the channel with distinct numerical values. If the same molecular systems is subject to Poiseuille flow, then the snapshots and density profiles are approximately the D

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

Figure 3. Velocity profiles in (a) Couette and (b) Poiseuille flow for the (2, 4) system. Coordinate y runs normal to the interface, and the profiles average over the other two directions. Points labeled (○) and (+) refer to the miscible and immiscible interactions, respectively. Reference 25. Copyright 2006 by the American Physical Society.

Table 1. Numerical Results for Slipa

same. The velocity profiles (Figure 3b) are parabolic on either side of the interface with a transition that is smooth and gradual in the miscible case (○) but very abrupt for immiscible liquids (+). However, while in Couette flow the shear stress has a constant value throughout both liquids, for Poiseuille flow the shear stress profiles are piecewise linearstraight lines of different slope in the two bulk regions that join smoothly across the interface. Similar results are found in the two other molecular systems considered, and a boundary condition that fits all of the data is that the apparent velocity discontinuity Δu is proportional to the shear stress S at the position of the interface, Δu = αS. More precisely, we say that slip occurs when the extrapolations of the two bulk velocity curves fail to intersect in the interfacial region, and in this case, there is always a dip in the total density profile whose minimum we identify as the position of the interface. The results of the various simulations given in Table 1 show that slip coefficient α depends on the two liquids involved, but not the flow configuration or the value of the driving force. A natural question is whether the behavior described here for simplified liquid models carries over to real materials. There is no suitable laboratory data because of the difficulties mentioned earlier, but MD simulations using realistic potentials have used density profiles and show, for example, that a dip at the interface occurs in the water−octane system46 but not for water−carbon tetrachloride.47 In the absence of a direct measurement of the velocity slip at an interface, one could determine the slip coefficient α for any pair of liquids by the same flux-based indirect procedure used for polymeric liquids. The geometry is precisely that of the MD simulations, and one would compare the measured fluid flux to that of the straightforward analytic calculation of a two-layer Couette flow with the interfacial slip condition imposed and then back out the value of α. Because the slip coefficient is a material property, the experiment could be performed on systems of any convenient size. 4.2. Interface Coalescence. The first simulations of drop coalescence arose in the course of a study of the deformation of drops suspended in a second immiscible fluid undergoing shear,26 at the time an active research area in low Reynolds number fluid mechanics.48 It is observed that a spherical drop in a linear shear flow deforms into an ellipsoid oriented at a

system

flow

Δu

S

α

(1, 1) immiscible interface

0.1 C 0.2 C 0.01 P 0.1 C 0.01 P 0.1 C 0.2 C 0.01 P 0.02 P 0.1 C 0.2 C 0.01 P 0.02 P

0.031 0.050 0.0 0.0 0.0 0.070 0.13 0.062 0.14 0.032 0.069 0.15 0.24

0.019 0.020 0.0 0.016 −0.055 0.012 0.021 0.011 0.025 0.0099 0.022 0.049 0.75

2.6 2.5

(2, 4) miscible interface (2, 4) immiscible interface

(4, 16) immiscible interface

0.0 0.0 5.7 5.9 5.6 5.8 3.2 3.2 3.1 3.1

The notation is that (l1, l2) refers to a liquid made of flexible chains of length l1 in contact with a second liquid of chains of length l2, with interactions either of the immiscible or miscible (LB) type. 0.1 C means Couette flow with wall velocities of ±0.1, and 0.01 P means gravity-driven Poiseuille flow with an acceleration of 0.01. Δu, S, and α are the apparent slip, shear stress at the interface, and Navier coefficient, respectively. All entries are in MD units: (σ, τ, ε). Reference 25. Copyright 2006, American Physical Society. a

specific angle with respect to the imposed velocity, and MD calculations gave results consistent with experiments and continuum calculations. The initial simulations involved a monatomic 736-atom drop suspended in a second immiscible liquid undergoing Couette flow generated by oppositely translating solid walls and were converted to a coalescence study49 by modifying the initial condition to two drops (424 atoms each) placed either above or below the center line of the cell and driven together by the flow. The process is indicated in Figure 4 and begins with atoms of the respective drops thermally fluctuating (Figure 4b) until coming within interaction range and drawing together because of their attraction. These atoms in turn attract the atoms of their original drop toward themselves, and a connecting filament forms (Figure 4c). At this point, the two drops have formed a dumbbell (Figure 4d), which reduces its surface area by contracting and ultimately forms a single spherical drop (Figure 4e). E

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

complete coalescence event. To illustrate this point and also to address the questions of whether the drops above are large enough to be representative of continuum behavior and whether molecular bonding changes the picture, we show in Figure 5 the time evolution of the gap region between two larger tetramer drops. In this case, the liquid consisted of 4atom linear chains bonded by a FENE interaction, and the initial state involved two pre-equilibrated drops of 157 214 atoms separated by a narrow gap. More precisely, the drops were translated together so that the distance between the centers of their interfacial regions was 2.5σ. In this case, the gap region shows the drops fluctuating back and forth while connected by several pairs of interacting molecules for a time of about 60τ, following which more molecules enter the gap to form a connecting bridge that thickens and expands outward. After 125τ (not shown here), the density in the center of the merged region has reached the bulk density, and the two drops have merged into a peanutlike shape. If the simulation is repeated with different atomic realizations of the drops, we find that the duration of the initial stage of forming a tenuous bridge is quite random because it depends very strongly on the value of the initial separation and is furthermore dependent on thermal fluctuations. However, once the drops are well connected the remainder of the merger process proceeds in the same way, requiring 50−100τ to reach the bulk density in the center of the region of first contact. Complete merging occurs only at a later time, which is controlled by the size of the drop. In the calculations described above, the molecules were relatively small, which allowed the fluid in the drops to rearrange and mix quickly once coalescence began. In a related calculation,50 we focused on the properties of longer FENE chains of up to 100 monomers, and considered drops on a substrate that coalesced either spontaneously by proximity or else were forced together by the motion of the substrate. In the first case, initially hemispherical drops of radius 35σ required about 1000τ for the outer shape to equilibrate, and even after 2000τ, the molecules of the two drops had mixed only within a central region of about 20σ thickness. The initial stages of the coalescence process were much the same as observed for the atomic and tetramer casesthermal fluctuations intermittently bring individual molecules from the two drops into interaction range, and eventually one of these interactions would persist long enough for other molecules to be attracted and form a tendril of molecules joining the two drops. An additional observation is that this initial stage is accompanied by unusually large fluctuations in the position of molecules near the contact region, although this behavior damps out quickly as the drops merge. Molecular-level studies of drop coalescence are surprisingly infrequent in the literature, in comparison to studies of coalescence and healing processes at solid−vapor interfaces, but there are useful studies of water drop coalescence51,52 where, in view of applications and the availability of experimental data, a reliable molecular model of water has been employed. Additional interesting calculations involve the role of coalescence during vapor−liquid coexistence53 and in the phase separation of Lennard-Jones fluids.54 In continuum calculations, coalescence is often described as a singularity in the sense that an abrupt redefinition of the interface is required and also because the new interface has a high curvature (Figure 5 suggests that it is as small as the inverse of a few molecular diameters) and a correspondingly

Figure 4. (a) Starting configuration for the simulation of small-drop coalescence driven by Couette flow. The drop and wall atoms are indicated by points, and the suspending fluid is not shown. (b−e) Evolution of coalescing drops at times (b) 40τ, (c) 45τ, (d) 50τ, and (e) 110τ after wall motion commences. (See text for more details.) From ref 49. Reprinted with permission from AAAS.

The distinctive feature of this calculation is the use of an external fluid in shear flow to drive the original drops into close proximity so that the interatomic potential can act. Freefloating drops in vacuum or vapor tend to fluctuate in and out of interaction range of each other, making the simulations difficult to control in their initial stages, and in a simulation, the initial separation must be adjusted carefully to capture a F

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

Figure 5. Coalescence of larger tetramer drops in vacuum. (Left to right) The times are 0, 30τ, 60τ, and 90τ.

spherical nanoparticle of radius 3σ (composed of 88 atoms and carved out of a cubic lattice) is placed in the bulk of one liquid in the vicinity of the interface. In one study, the particle was chosen to have neutral wettability by adjusting the symmetric interaction coefficients with both liquids such that c1p = c2p = 0.5, where p denotes the particle. In another case study, a Janus particle was created by assigning different interaction coefficients to each face of the particle. In other words, the atoms in the lower half of the particle (J1) strongly prefer one liquid (i.e., liquid 1 (c1J1 = 1, c2J1 = 0)), whereas atoms in the upper half (J2) prefer the other liquid (i.e., liquid 2 (c2J2 = 1, c1J2 = 0)). The particle is quasi-statically pushed through the interface along the normal direction. At any normal distance y* from the interface, the information on the density profiles of each liquid is available during the translation process. Data on the region of interest (i.e., in the vicinity of the particle) is shown in Figure 6 (top panel). It was observed that the effective colloid−interface interaction altered by the particle wettability strongly influences the interface distortion, as depicted in the density profiles. Consequently, a Janus particle, initially dispersed in the bulk of one liquid, will carry a thin film of this same liquid after it is forced to pass through the interface and is pushed into the other liquid. Another valuable piece of information accessible in MD simulations is the average normal force ⟨Fy⟩ that the surrounding medium exerts on the particle at any given normal distance y*. Hence, because the temperature is fixed by a thermostat, the Helmholtz free energy difference Δ- for any given distance from the interface can be cast as

high surface tension force. On the molecular scale, however, a coalescence event is just the rearrangement of a modest number of molecules in the interfacial region along with a small inward displacement of the drops’ centers of mass, and nothing is singular. More quantitatively, one can examine the velocity and stress fields produced in the course of the MD simulation. Typical local values of stress during and after coalescence are similar to those for single drops in equilibrium, whereas the velocity field consists mostly of thermal fluctuations until the fluid begins to fill the gap between the drops, at which point the outward velocity becomes statistically significant. The subsequent rearrangement of the two connected drops into a single sphere is so slow that the velocity field is hidden by thermal fluctuations and the stress field is featureless. Lastly, having discussed the molecular nature of coalescence, one may ask the same question in the context of interface rupture, which occurs, for example, when a drop is torn apart in a strong shear flow48 or when a cylinder of liquid breaks up into spheres as a result of a (Rayleigh) surface tension instability.55 Here, MD simulations26 indicate a process roughly the reverse of coalescence. The interface evolves into a shape with a thinning neck separating two larger bodies of liquid, and the molecules in the neck gradually withdraw into one of the liquid bodies on either side of the neck until they are out of interaction range and the interface “breaks.” On the continuum level, the process may involve rapid rearrangement and appear to be singular, but as in coalescence, on smaller scales only a gradual evolution of the molecular configuration occurs, which does not produce any exceptionally large velocity or stress. 4.3. Particles at Interfaces. MD simulations are central to gaining insight into the inherent molecular-level behavior of colloidal particles and macromolecules at interfaces.56−59 We have employed simulation tools to examine the dynamics of colloidal particle translation through the interface of two immiscible liquids. Particles of various wettability have been studied with an emphasis on amphiphilic Janus particles. These particles, first highlighted by de Gennes,60 possess different surface chemistry on two sides and are currently leading-edge building blocks in the field of colloid and interface science because of their wide array of potential applications.61,62 In this set of simulations, the system was composed of about 20,000 atoms that form the four species: two immiscible liquids, a solid nanoparticle, and the confining solid walls.27 The liquid−liquid interface was created by setting the liquid self-interaction coefficients (eq 1) to c11 = c22 = 1, whereas the cross coefficient is set to c12 = 0.5. The solid walls that are parallel to the plane of the interface constrain the interface to the center of the simulation box, whereas periodic boundary conditions are applied in the two directions parallel to the interface. A

Δ-(y) = -(y) − -(yref ) = −

∫y

y

ref

⟨Fy⟩y * dy*

(3)

The juxtaposition of density profiles and free-energy landscapes reveals that as higher-energy states become accessible by driving the system through a series of nonequilibrium states the kinetic adsorption of a Janus particle exhibits a significant difference compared to a plain particle, and this disparity is embodied in the asymmetric free-energy landscape for the former (Figure 6, bottom panel). The observed phenomenon is highly relevant to experiments in the area of Pickering emulsions and interfacially templated materials because it points to the significance of high-energy states when amphiphilic particles are translated between immiscible liquids, a factor not pertinent to plain particles. This example illustrates the power of MD simulations as a tool for elucidating the microscopic details of a system invisible to most experimental tools. The reported findings are likely to inspire future experimental work on the stabilization of G

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

evaluating the rotational free-energy barriers of a system involving a nanoparticle to which an atomic ridge was attached. The simulation setup was very similar to the translational case discussed earlier: the two immiscible liquids were confined by a set of solid walls placed parallel to the plane of the interface, and the nanoparticle was positioned in the plane of the interface. A snapshot of the simulation box is shown in Figure 7a. Particles of size 3σ, 6σ, and 8.3σ were investigated, but here

Figure 7. (a) Snapshot of the simulation box. (b, top) Cylindrical nanoparticle of radius 6σ. (b, bottom) Same particle with the attached surface feature (blue atoms). (c) Free-energy landscape of the particle without satellite feature (Δ red line) and with the satellite feature (▼ blue line). The subtraction of the free energy computed for each case (dashed gray line) shows the net effect of the added satellite feature at θ = 0. Adapted with permission from ref 28. Copyright 2014, AIP Publishing LLC.

Figure 6. (Top panel) Density profiles of liquid 1 as a function of the 2D distance of a solvent atom from the axis normal to the interface passing through the particle’s center of mass (horizontal axis) and the height of the solvent atom from the particle’s center of mass (vertical axis) at different times during the process of particle adsorption to and desorption from the interface. (Top to bottom) For systems with (I) a plain particle and (II) a Janus particle. (Left to right) (a) A particle in the bulk of liquid 1, (b) a particle adsorbing onto the interface, (c) a particle sitting at the interface, (d) a particle desorbing from the interface, and (e) a particle in the bulk of liquid 2. (Bottom panel) System free energy as a function of the particle distance (yp) from the interface (y = 0). The top curve (red Δ) and bottom curve (purple ▽) correspond to the free-energy differences for plain and Janus particles, respectively. Reference 27. Reproduced by permission of the Royal Society of Chemistry.

we discuss only the findings for the particle with radius 6σ. To simulate the surface roughness, a ridge was attached to one side of the nanoparticle; see Figure 7b. Through successive quasistatic rotations of the particle out of the plane of the interface (around the z axis), we derived the system free energy for any specific orientation with a rotated angle θ*. These calculations are analogous to the translation case, and one needs only to substitute the force and distance in eq 3 with the torque and angle, respectively. The final results for the rotational freeenergy landscape are shown in Figure 7c for both particle structures displayed in Figure 7b. For a nanoparticle with an annexed satellite feature that locally extends the diameter, the rotational free energy exhibits two energy wells centered at θ = 0 and θ = π in which the surface ridge intersects the liquid− liquid interface. These wells are significantly deeper than those observed for the particle without the added feature, and in this example, the difference caused by the atomic ridge is as large as Δ- = |Δ- 2 − Δ- 1| ≃ 12−13ε. The presence of a larger energy barrier manifests itself in a prolonged lifetime of the

interfaces by utilizing the unique interfacial properties of Janus particles. In addition to chemical anisotropy (i.e., Janus character), another intriguing observation is associated with the intrinsic surface inhomogeneity (i.e., roughness) of colloidal particles and its consequences on the dynamic behavior of a particle straddling an interface. In a recent study,28 we focused on the rotation and angular behavior of cylindrical nanoparticles. Using simulation tools, we recognized that these ubiquitous surface nanofeatures indeed alter the binding dynamics and engender particle locking.28,63 This fact was corroborated by H

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

states located at θ = 0 and θ = π. We evaluated the average time that each nanoparticle spends in these two states before hopping to either θ = π/2 or θ = 3π/2 and recognized that at an elevated temperature of kBT/ε = 3 the satellite feature results in a 5-fold increase in the lifetime. Simulations were done at different temperatures, and it was concluded that the rotational locking caused by the attached feature can be tuned by adjusting the temperature (i.e., at lower temperatures, the angular locking is more severe and the proportionality is not linear). A very valuable feature of the MD simulations was the insight that they provided into the molecular structure of the liquid near the particle surface. We used the computer as a fine microscope to resolve the particle−liquid interfacial boundaries that are presented in the density (ρ) profiles of Figure 8a. A

These results are significant to a number of existing and emerging technological applications based on nanoparticles at liquid−liquid interfaces where it is desirable to tame Brownian effects. Because the roughness and the single protrusions studied are well within current synthesis capabilities,64,65 one can imagine science and technology that could be inspired by the observed rotational localization effects.

5. CONCLUDING REMARKS AND OUTLOOK The examples of MD simulations discussed above also illustrate some of the potential pitfalls and limitations of the technique. One obvious shortcoming is the scale of the calculation: no current simulation with atomic resolution can approach the physical size of a macroscopic system. To the extent that the potential energy is accurate for the materials in question, most MD simulations literally describe systems on the nanometer scale, and some justification is needed before the results can be used more generally. In channel flow simulations, such as those in section 4.1, it is relatively easy to address the size dependence by testing for any variation of the result as the dimensions of the channel increase. Interfacial parameters such as the slip length or surface tension depend only on the fluid structure near the interface, and the key requirement is that the system is wide enough for the liquids to relax to their bulk molecular structure on either side of it. Similarly, in “operational” computations of transport coefficients, where, for example, the viscosity is extracted from the velocity and stress fields obtained in a flow, the (mild) requirement is that the liquid exhibits a homogeneous region whose properties are insensitive to the overall dimensions. Beyond this, there is the issue of observing a robust signal in the presence of thermal fluctuations. Typical thermal velocities are hundreds of meters per second, and one is able to observe macroscopic flow velocities only because the statistical sample size is macroscopic as well. Simulations do not have the latter property, and it is often difficult to obtain accurate numerical results without very long averaging periods. In situations where the physical relaxation time is long, as in the free-energy calculations in section 4.3, the problem is exacerbated. Generally speaking, the stress tensor is the most difficult to resolve because it is directly related to the interatomic force, which is a rapidly varying function of atomic position. In practice, accurate stress results are easy to obtain only when there is a 1D variation; otherwise, they are often very noisy. The velocity is, so to speak, an integral of the force, and because integration is a smoothing operation, the velocity field is usually measurable in situations where there is a 1D or 2D variation. A density field involves atomic positions that are integrals of velocity and smoother still and it tends to be measurable with relatively high precision in almost all cases. The coalescence simulations, in contrast, illustrate a serious limitation of MD. In this problem, there is considerable interest in scaling laws for the evolution of the contact radius, for example, where laboratory experiments easily span many decades of time and radius.66 The larger simulations discussed in section 4.2 involve drops of radius 35σ, only about 1 order of magnitude greater than the interfacial thickness. Recent work by Pothier and Lewis67 examines the scaling behavior in MD simulations of 2D and 3D drop coalescence but under similar size limitations. Parallel computing platforms can now treat systems with a thousand times as many atoms, but this provides only an additional factor of 10 in radius. Furthermore, increases in the simulation size raise the concomitant issue of increased

Figure 8. (a) Number density ρ from MD calculations for θ = 0. (b) Mean-field potential um = −kBT log ρ from MD calculations for θ = 0. (Top panel) Nanoparticle without a ridge. (Bottom panel) Nanoparticle with the attached ridge. Note that because of the symmetry of the structure only a quarter of the solid particle is shown in these plots. Adapted with permission from ref 28. Copyright 2014, AIP Publishing LLC.

mean-field potential of um = −kBT log ρ was computed from the number density calculations to investigate the nanoparticle−liquid interaction near the boundary (Figure 8b). Besides the irregular nature of the particle surface, from the density profiles we readily observed the variations in the local number density on the particle surface and identified the regions for which the density is significantly higher. The meanfield potential um contours revealed that these are very low energy spots (um/ε < −1) to which fluid atoms are attracted and potentially remain trapped within a deep energy well and cause contact line pinning and other effects. Density fluctuations in the vicinity of the particle surface stem from the layering of the liquid atoms. These fluctuations decay further away from the solid surface. For the case with the added feature (Figure 8, bottom panel), a large region of low energy appears in the gap between the particle surface and the satellite feature. This gap can thus cause the pinning of the three-phase contact line at certain angular orientations, which can entail additional “locking” effects discussed in more detail in ref 28. I

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

Notes

simulation durations. The computation time per time step is at best linear in the number of atoms, and the physical time associated with the process often grows rapidly. For example, the development of viscous flows is controlled by the diffusion of vorticity, whose time scale is proportional to the square of the system size. Although we have stressed applications in which generic liquids based on simple potentials suffice, there are of course situations in which a realistic molecular model is required. Simulations of biological systems are an obvious example where very specific chemical processes and interactions occur, and molecular models incorporating the relevant structure and bonding68 are required. A situation where the use of both generic and realistic models complement each other occurs in surfactant systems. A “generic” surfactant molecular model may be based on a FENE chain where one end attracts and the other repels the solvent.69 We have used this technique ourselves to study the behavior of surface surfactant phases,70 and analogous models have successfully captured the terraced spreading71 of polymeric liquid drops72 but fail to reproduce the superspreading73 of trisiloxane surfactant drops.74 In the latter case, more realistic modeling of the trisiloxane molecule seems necessary, but even these simulations have not been successful so far,75 presumably because of limitations on the size of the simulated drops. Lastly, there are cases where an extrapolation from MD to macroscopic scales fails because the mechanism of the physical process is different on atomic and continuum scales. For instance, in the case of the coalescence phenomenon, although two bodies of the same liquid would be expected to coalesce when they are placed in contact, it has been observed76 that temperature differences or the relative velocity may suppress the coalescence indefinitely. For example, a small drop placed on the surface of a rotating pool of the same liquid can float there indefinitely because the thin layer of air separating the drop from the surface of the rotating pool is set into motion, which is a lubrication flow that develops a high enough pressure to overcome gravity and suspend the drop. An MD simulation of this process77 reproduced the phenomenon but through an entirely different mechanism. There, as in Figure 5, a tendril connecting the drop to the liquid pool forms and in the absence of rotation thickens and draws the drop downward to coalesce. The effect of the rotation is that liquid molecules evaporating from the pool retain their rotating velocity and when interacting with the drop atoms exert a torque on the drop, which causes it to rotate and break the incipient connecting tendril of atoms. This feature article aims to underscore the enormous potential that lies in exploiting simulation findings to present opportunities for experiments and theory. The topics addressed highlight only a few implementations of MD simulations, but we hope that they give a flavor of their potential utility in colloid and interface science. In essence, simulations are intended to facilitate the technological progress in the rapidly evolving field of colloid and interface science. Ultimately, as Frenkel already remarked: “simulations are not a good tool to summarize our understanding of nature; however, they can provide important insights.”3

■

The authors declare no competing financial interest. Biographies

Sepideh Razavi received her Masters’ degree in chemical engineering from Sharif University of Technology, Iran in 2007. From 2007 to 2009, she was a research assistant in Prof. A. Shojaei’s group at Sharif University, working on the rheological and mechanical behavior of recycling-based polymeric composites with application as composite railroad ties. In 2010, she joined the Chemical Engineering Department at the City College of New York and is currently pursuing her doctoral degree. Her research focuses on the behavior of colloidal Janus particles at interfaces using both experimental and simulation tools.

Joel Koplik received his Ph.D. from the University of California at Berkeley in 1974 in elementary particle theory. After postdoctoral positions at Columbia University, the Institute for Advanced Study, and the École Normale Supérieure, he joined Schlumberger-Doll Research in 1979 and switched to fluid mechanics. He moved to the City College of New York in 1989 as a professor of physics. His research interests are in molecular fluid mechanics, where he studies the dynamics of surfactants in solution, interfacial boundary conditions, nanoscale flows on patterned substrates, droplet coalescence and splashing, and systems of colloidal particles.

AUTHOR INFORMATION

Additional interests include fluid and particulate flow in geological

Corresponding Authors

fractures and random porous media and vortex dynamics in superfluid

*E-mail: [email protected]. *E-mail: [email protected].

turbulence. J

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

(11) Trokhymchuk, A.; Alejandre, J. Computer simulations of liquid/ vapor interface in Lennard-Jones fluids: some questions and answers. J. Chem. Phys. 1999, 111, 8510−8523. (12) Blas, F. J.; Bravo, A. I. M. V.; Miguez, J. M.; Pineiro, M. M.; MacDowell, L. G. Vapor-liquid interfacial properties of rigid-linear Lennard-Jones chains. J. Chem. Phys. 2012, 137, 084706−1-11. (13) Lukyanov, A. V.; Likhtman, A. E. Relaxation of surface tension in the free-surface boundary layer of simple Lennard-Jones liquids. J. Chem. Phys. 2013, 138, 034712-1−034712-19. (14) Neyt, J. C.; Wender, A.; Lachet, V.; Ghoufi, A.; Malfreyt, P. Molecular modeling of the liquid-vapor interfaces of a multicomponent mixture: prediction of the coexisting densities and surface tensions at different pressures and gas compositions. J. Chem. Phys. 2013, 139, 024701-1−024701-11. (15) Lauga, E.; Brenner, M. P.; Stone, H. A. Microfluidics: The NoSlip Boundary Condition. In Handbook of Experimental Fluid Mechanics; Foss, J.,Tropea, C., Yarin, A., Eds.; Springer: New York, 2005. (16) Vinogradova, O. I.; Belayev, A. V. Wetting, Roughness and Hydrodynamics Slip. In Nanoscale Liquid Interfaces; Ondarcuhu, T., Aimé, J.-P., Eds.; Pan Stanford Publishing: Singapore, 2013. (17) Li, Y.; Xu, J.; Li, D. Molecular dynamic simulation of nanoscale liquid flows. Microfluid. Nanofluid. 2010, 9, 1011−1031. (18) Neto, C.; Evans, D. R.; Bonaccurso, E.; Butt, H.-J.; Craig, V. S. J. Boundary slip in Newtonian liquids: a review of experimental studies. Rep. Prog. Phys. 2005, 68, 2859−2897. (19) Bouzigues, C. I.; Bocquest, L.; Charlaix, E.; Cottin-Bizonne, C.; Cross, B.; Joly, L.; Steinberger, A.; Ybert, C.; Tabeling, P. Using surface force apparatus, diffusion and velocimetry to measure slip lengths. Philos. Trans. R. Soc. A 2008, 366, 1455−1468. (20) Batchelor, G. K. An Introduction to Fluid Dynamics; Cambridge University Press: Cambridge, U.K., 1967. (21) Lam, Y. C.; Jiang, L.; Yue, C. Y.; Tam, K. C.; Li, L.; Hu, X. Interfacial slip between polymer melts studied by confocal microscopy and theological measurements. J. Rheol. 2003, 47, 795−807. (22) Park, H. E.; Lee, P. C.; Macosco, C. W. Polymer-polymer interfacial slip by direct visualization and by stress reduction. J. Rheol. 2010, 54, 1207−1218. (23) Zartman, G. D.; Wang, S.-Q. A particle tracking velocimetric study of interfacial slip at polymer−polymer interfaces. Macromolecules 2011, 44, 9814−9820. (24) Barsky, S.; Robbins, M. O. Molecular dynamics study of slip at the interface between immiscible polymers. Phys. Rev. E 2001, 63, 021801-1−021801-7. (25) Koplik, J.; Banavar, J. R. Slip, immiscibility and boundary conditions at the liquid-liquid interface. Phys. Rev. Lett. 2006, 96, 044505-1−044505-5. (26) Koplik, J.; Banavar, J. R. Molecular dynamics of interface rupture. Phys. Fluids A 1993, 5, 521−536. (27) Razavi, S.; Koplik, J.; Kretzschmar, I. The effect of capillary bridging on the Janus particle stability at the interface of two immiscible liquids. Soft Matter 2013, 9, 4585−4589. (28) Razavi, S.; Kretzschmar, I.; Koplik, J.; Colosqui, C. E. Nanoparticles at liquid interfaces: rotational dynamics and angular locking. J. Chem. Phys. 2014, 140, 014904-1−014904-12. (29) Pieranski, P. Two-dimensional interfacial colloidal crystals. Phys. Rev. Lett. 1980, 45, 569−572. (30) Glotzer, S. C.; Solomon, M. J.; Kotov, N. A. Self-assembly: From nanoscale to microscale colloids. AIChE J. 2004, 50, 2978−2985. (31) Crossley, S.; Faria, J.; Shen, M.; Resasco, D. E. Solid nanoparticles that catalyze biofuel upgrade reactions at the water/oil interface. Science 2010, 327, 68−72. (32) Kaz, D. M.; McGorty, R.; Mani, M.; Brenner, M. P.; Manoharan, V. N. Physical ageing of the contact line on colloidal particles at liquid interfaces. Nat. Mater. 2012, 11, 138−142. (33) Lin, Y.; Boker, A.; Skaff, H.; Cookson, D.; Dinsmore, A. D.; Emrick, T.; Russell, T. P. Nanoparticle assembly at fluid interfaces: Structure and dynamics. Langmuir 2005, 21, 191−194.

Ilona Kretzschmar received her Ph.D. from the Department of Chemistry at the Technical University of Berlin in 1999. During her graduate research with Prof. H. Schwarz, she studied reactions of metal cations with organic molecules in the gas phase, employing Fourier transform ion-cyclotron resonance and guided-ion-beam mass spectrometry. From 2000 to 2002, she was a Feodor-Lynen postdoctoral fellow at Harvard University, working on hydrocarbon radical rearrangement reactions on clean and modified metal surfaces. In 2002, she joined Prof. M. A. Reed’s molecular electronics group in the Department of Electrical Engineering at Yale University as research associate to familiarize herself with the field of molecular electronics. Currently, she is a professor of chemical engineering at the City College of New York. Her research interests focus on Janus and patchy particles including their preparation, molecular and field-directed assembly, and application.

■

ACKNOWLEDGMENTS The work highlighted here was partially supported by the City University of New York High Performance Computing Center under NSF grants CNS-0855217 and CNS-0958379. Support from the National Science Foundation through awards NSF CBET-1067501 and NSF CBET-1264550 is acknowledged. Much of J.K.’s research was supported in the past by the late, lamented NASA microgravity program. We thank Hiroshi Watanabe for providing a colored version of Figure 1b.

■

REFERENCES

(1) Haile, J. M. Molecular Dynamics Simulation: Elementary Methods; Wiley Press: New York, 1997. (2) Borhani, D. W.; Shaw, D. E. The future of molecular dynamics simulations in drug discovery. J. Comput. Aided Mol. Des. 2012, 26, 15−26. (3) Frenkel, D. Simulations: the dark side. Eur. Phys. J. Plus 2013, 128, 1−21. (4) Cheatham, T. E.; Case, D. A. Twenty-five years of nucleic acid simulations. Biopolymers 2013, 99, 969−977. (5) Friedman, R.; Boye, K.; Flatmark, K. Molecular modelling and simulations in cancer research. Biochim. Biophys. Acta 2013, 1836, 1− 14. (6) Zahn, D. Molecular dynamics simulation of ionic conductors: perspectives and limitations. J. Mol. Model. 2011, 17, 1531−1535. (7) Lindahl, E. Membrane proteins: molecular dynamics simulations. Curr. Opin. Struct. Biol. 2008, 18, 425−431. (8) Lyklema, J. Fundamentals of Interface and Colloid Science; Academic Press: London, 1995. (9) Feller, S. E. Molecular dynamics simulations of lipid bilayers. Curr. Opin. Colloid Interface Sci. 2000, 5, 217−223. (10) Mecke, M.; Winkelmann, J.; Fischer, J. Molecular dynamics simulation of the liquid-vapor interface: the Lennard-Jones fluid. J. Chem. Phys. 1997, 107, 9264−9270. K

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Invited Feature Article

(34) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, U.K., 1987. (35) Parseghian, V. A. van der Waals Forces; Cambridge University Press: New York, 2005. (36) Watanabe, H.; Ito, N.; Hu, C.-K. Phase diagram and universality of the Lennard-Jones gas-liquid system. J. Chem. Phys. 2012, 136, 204102-1−204102-9. (37) Kremer, K.; Grest, G. S. Molecular dynamics simulation for polymers in the presence of a heat bath. Phys. Rev. A 1986, 36, 3628− 3631. (38) Cheng, S. F.; Lechman, J. B.; Plimpton, S. J.; Grest, G. S. Evaporation of Lennard-Jones fluids. J. Chem. Phys. 2011, 134, 2247041−224704-13. (39) Koplik, J.; Zhang, R. Nanodrop impact on solid surfaces. Phys. Fluids 2013, 25, 022003−1-12. (40) Thompson, P. A.; Robins, M. O. Simulations of contact-line motion: slip and the dynamic contact angle. Phys. Rev. Lett. 1989, 63, 766−769. (41) Hunter, R. J. Foundations of Colloid Science, 2nd ed.; Oxford University Press: New York, 2001. (42) Frenkel, D.; Smit, B. Understanding Molecular Simulation, 2nd ed.; Academic Press: San Diego, 2002. (43) Evans, D. J.; Morriss, G. Statistical Mechanics of Nonequilibrium Liquids, 2nd ed.; Cambridge University Press: New York, 2008. (44) Soddemann, T.; Dünweg, B.; Kremer, K. Dissipative particle dynamics: a useful thermostat for equilibrium and nonequilibrium molecular dynamics simulations. Phys. Rev. E 2003, 68, 046702-1− 046702-11. (45) Thompson, P. A.; Troian, S. A general boundary condition for flow at solid surfaces. Nature 1997, 389, 360−362. (46) Patel, H. A.; Naumann, E. B.; Garde, S. Molecular structure and hydrophobic solvation thermodynamics at an octane-water interface. J. Chem. Phys. 2003, 119, 9199−9206. (47) Senapathi, S.; Berkowitz, M. L. Computer simulation study of the interface width of the liquid/liquid interface. Phys. Rev. Lett. 2001, 87, 176101-1−176101-4. (48) Rallison, J. M. The deformation of small viscous drops and bubbles in shear flow. Annu. Rev. Fluid Mech. 1984, 16, 45−65. (49) Koplik, J.; Banavar, J. R. Molecular structure of the coalescence of liquid interfaces. Science 1992, 257, 1664−1666. (50) Pal, S.; Koplik, J.; Banavar, J. R. Dynamics of nanoscale droplets. Phys. Rev. E 2002, 65, 021504-1−021504-14. (51) Svanberg, M.; Ming, L.; Markovic, N.; Pettersson, J. B. C. Collision dynamics of large water clusters. J. Chem. Phys. 1998, 108, 5888−5897. (52) Zhao, L.; Choi, P. Molecular dynamics simulation of the coalescence of nanometer-sized water droplets in n-heptane. J. Chem. Phys. 2004, 120, 1935−1942. (53) Roy, S.; Das, S. K. Dynamics and growth of droplets close to the two-phase coexistence curve in fluids. Soft Matter 2013, 9, 4178−4187. (54) Peters, G. H.; Eggebrecht, J. Obervation of droplet growth and coalescence in phase-separating Lennard-Jones fluids. J. Phys. Chem. 1995, 99, 12335−12340. (55) Lord Rayleigh. On the instability of a viscous liquid under capillary force. Philos. Mag. 1892, 34, 145−150. (56) Ranatunga, R. J. K. U.; Nguyen, C. T.; Wilson, B. A.; Shinoda, W.; Nielsen, S. O. Molecular dynamics study of nanoparticles and nonionic surfactant at an oil-water interface. Soft Matter 2011, 7, 6942− 6952. (57) Luo, M.; Mazyar, O. A.; Zhu, Q.; Vaughn, M. W.; Hase, W. L.; Dai, L. L. Molecular dynamics simulation of nanoparticle self-assembly at a liquid-liquid interface. Langmuir 2006, 22, 6385−6390. (58) Cheung, D. L.; Bon, S. A. F. interaction of nanoparticles with ideal liquid-liquid interfaces. Phys. Rev. Lett. 2009, 102, 066103−1-4. (59) Cheung, D. L.; Carboneb, P. How stable are amphiphilic dendrimers at the liquid−liquid interface? Soft Matter 2013, 9, 6841− 6850. (60) de Gennes, P.-G. Soft matter (Nobel Lecture). Angew. Chem., Int. Ed. 1992, 31, 842−845.

(61) Kumar, A.; Park, B. J.; Tu, F.; Lee, D. Amphiphilic Janus particles at fluid interfaces. Soft Matter 2013, 9, 6604−6617. (62) Walther, A.; Müller, A. H. E. Janus particles: synthesis, selfassembly, physical properties, and applications. Chem. Rev. 2013, 113, 5194−5261. (63) Colosqui, C. E.; Morris, J. F.; Koplik, J. Colloidal adsorption at fluid interfaces: regime crossover from fast relaxation to physical aging. Phys. Rev. Lett. 2013, 111, 028302−1-5. (64) Glotzer, S. C.; Solomon, M. J. Anisotropy of building blocks and their assembly into complex structures. Nat. Mater. 2007, 6, 557−562. (65) Xia, Y.; Xia, X.; Wang, Y.; Xie, S. Shape-controlled synthesis of metal nanocrystals. MRS Bull. 2013, 38, 335−344. (66) Eggers, J.; Lister, J. R.; Stone, H. A. Coalescence of liquid drops. J. Fluid Mech. 1999, 401, 293−310. (67) Pothier, J.-C.; Lewis, L. J. Molecular dynamics study of the viscous to inertial crossover in nanodroplet coalescence. Phys. Rev. B 2012, 85, 115447−1-10. (68) Schlick, T. Molecular Modeling and Simulation; Springer: New York, 2002. (69) Schmid, F. Systems Involving Surfactants. In Computational Methods in Surface and Colloid Science; Borówko, M., Ed.; Surface Science Series; Marcel Dekker: New York, 2000; Vol. 89, pp 631−683 (70) Tomassone, M. S.; Couzis, A.; Maldarelli, C.; Banavar, J. R.; Koplik, J. Molecular dynamics simulation of gaseous-liquid phase transitions of soluble and insoluble surfactants at a fluid interface. J. Chem. Phys. 2001, 115, 8634−8642. (71) Czabat, A. M.; Heslot, F.; Fraysse, N. Molecular layering in the spreading of wetting liquid drops. Nature 1989, 338, 640−642. (72) d’Ortona, U.; De Coninck, J.; Koplik, J.; Banavar, J. R. Terraced spreading mechanisms for chain molecules. Phys. Rev. E 1996, 53, 562−569. (73) Venzmer, J. Superspreading-20 years of physicochemical research. Curr. Opin. Colloid Interface Sci. 2011, 16, 307−313. (74) Shen, Y.; Couzis, A.; Koplik, J.; Maldarelli, C.; Tomassone, M. S. Molecular dynamics study of the influence of surfactant structure on the spreading of droplets on solid surfaces. Langmuir 2005, 21, 12160−12170. (75) Halverson, J.; Maldarelli, C.; Couzis, A.; Koplik, J. Wetting of hydrophobic substrates by nanodroplets of aqueous trisiloxane and alkyl polyethoxylate surfactant solutions. Chem. Eng. Sci. 2009, 64, 4657−4667. (76) Nietzel, G. P.; Dell’Aversana, P. Noncoalescence and nonwetting behavior of liquids. Annu. Rev. Fluid Mech. 2002, 34, 267−89. (77) Dell’Aversana, P.; Koplik, J.; Banavar, J. R. Suppression of coalescence by shear and temperature gradients. Phys. Fluids 1996, 8, 15−28.

L

dx.doi.org/10.1021/la500376z | Langmuir XXXX, XXX, XXX−XXX