An All-Atom Molecular Dynamics Study of Water-Dodecane Interface

Article pubs.acs.org/JPCC

Cite This: J. Phys. Chem. C 2018, 122, 687−693

All-Atom Molecular Dynamics Study of Water−Dodecane Interface in the Presence of Octanol Baofu Qiao*,† and Wei Jiang*,‡ †

Chemical Sciences and Engineering Division, Argonne National Laboratory, Argonne, Illinois 60439, United States Leadership Computing Facility, Argonne National Laboratory, Argonne, Illinois 60439, United States

‡

S Supporting Information *

ABSTRACT: In this work, classical atomistic molecular dynamics (MD) simulations were performed to investigate the structural properties of octanol at water/dodecane interface. Multiple combinations of simulation conditions and force fields were employed in order to obtain high fidelity interfacial structure. The analysis of noncovalent interactions shows that the convergence of interfacial octanol distribution requires a simulation time scale larger than 100 ns. Finite size effect diminishes with increasing dimension of simulation box both laterally and vertically. In the largest simulated systems with vanishing finite size effect, the additive CHARMM force field shows similar interfacial aggregation behavior with the DRUDE polarizable model. Self-assembly of octanol molecules at water/oil interface is observed, expected to reduce the free energy barrier between water phase and oil phase. This work provides an MD simulation benchmark for structural study of complex liquid/liquid interface with cosurfactant reducing interfacial tension between two immiscible phases.

■

INTRODUCTION Liquid−liquid interfaces are of paramount interest due to their abundance and significant roles in fundamental research and technical applications, e.g., assembly of nanoparticles, 1 functionalized carbon materials,2 colloidal self-assembly,3 and ion transport,4,5 to name a few. Because of the sharp change in the dielectric permittivities of two immiscible liquids across the interface, the structures are greatly different from those in bulk solutions.6−8 Experimentally X-ray and neutron scattering techniques9 have been extensively employed in investigating the assembly of amphiphilic molecules at liquid−liquid interfaces, along with some other experimental techniques including sum frequency generation spectroscopies,10−12 nuclear magnetic resonance,5 and so on. However, the molecular structures are still in debate due to the limited resolutions of experimental instruments. Complementary to experimental techniques, all-atom computer simulations have been proven valuable in a broad spectrum of aqueous13 and organic systems.14,15 Nevertheless the liquid−liquid interfaces are only limitedly investigated with computer simulations,6,7,16−24 among which the interfacial tension25−27 was focused on. Furthermore, due to the high computing cost, systematic convergence analysis of all-atom simulation was generally avoided. Alternatively, coarse-grained model has been employed to avoid the demanding need of high performance computing resources, at the expense of loss of quantitative or fine structure information.28 In the present work, we aim to study structural properties of liquid−liquid interface with full capability of all-atom brutalforce molecular dynamics (MD) simulations, namely running long single-trajectory without sampling enhancement ap© 2017 American Chemical Society

proaches. In specific, ternary solution octanol/water/dodecane is employed for its generality in microemulsions, for which short-chain alcohol molecules are frequently added as cosurfactants for fine-tuning of microemulsions.5 Increasing number of studies on alcohol molecules have been reported recently on the novel structures of the so-called “surfactantfree” microemulsions in bulk solutions.29−31 Nevertheless our quantitative knowledge in terms of their influences at the liquid−liquid interfaces is still limited. The present work is thus poised to be a benchmark of atomistic simulation in studying the weak self-assembly of alcohol at liquid−liquid interface. We expect that the results will pave the way toward more complex multiphase liquid−liquid interface simulations for various engineering applications. From point of view of MD simulation, finite size effect and induced polarization due to heterogeneous interfacial structure have been long-standing issues for liquid−liquid interface simulations. In this report, a variety of systems were studied to calibrate the effects of finite size simulation box in liquid−liquid interface simulations, as well as the influences of atomistic (additive CHARMM and polarizable DRUDE) force fields. In the last two decades, polarizable models have demonstrated advantage over the traditional additive force fields in simulating interface systems.32,33 A simulation with the DRUDE model is thus performed to catch significant polarization effects that might be missed in additive force fields. It needs to be emphasized that this report focuses on interfacial structure investigation. A Received: November 7, 2017 Revised: December 15, 2017 Published: December 19, 2017 687

DOI: 10.1021/acs.jpcc.7b10997 J. Phys. Chem. C 2018, 122, 687−693

Article

The Journal of Physical Chemistry C Table 1. Compositions and Dimensions of the Simulated Systemsa system CHARMM-s1 CHARMM-s2 CHARMM-s3 CHARMM-L DRUDE-L

no. of octanols 300 300 300 2400 2400

(1.0) (1.0) (1.0) (1.4) (1.4)

no. of dodecanes 460 460 460 3600 3600

(1.6) (1.6) (1.6) (2.0) (2.0)

no. of waters 8000 8000 8000 32000 32000

(27.6) (27.6) (27.6) (18.1) (18.1)

LXY/nmb

LZ/nmb

4.7 6.1 7.2 10 10

∼21.8 ∼13.0 ∼9.3 ∼29.4 ∼29.4

a

The overall concentrations of the components (M) are listed in the parentheses. bEquilibrium simulation box edge lengths in XY plane (LXY) and Z dimension (LZ).

oil separately coupled with the Nosé−Hoover algorithm.41 The semiisotropic Parrinello−Rahman algorithm42 was applied to couple the pressure in the Z dimension (reference pressure PZ = 1 bar, compressibility of 4.5 × 10−5 bar−1 and relaxation time of 4 ps) while keeping the XY plane fixed with a compressibility of 0. Three-dimensional periodic boundary conditions (PBC) were used. The short-range Coulomb interactions were calculated up to 1.2 nm with the longrange interactions calculated using the particle mesh Ewald (PME) method.43,44 The short-range van der Waals (vdW) interactions using the Lennard-Jones 12−6 potential were calculated up to 1.2 nm, along with the long-range dispersion correction for energy and pressure. All chemical covalent bonds were constrained by means of the LINCS algorithm,45,46 so that a leapfrog integration time of 2 fs was stably supported. Each simulation was performed for 200 ns with a saving frequency of 10 ps per frame to collect the simulation trajectory. Influences of atomistic force fields were investigated with another two simulations with significantly larger simulation box and around 0.3 million atoms (CHARMM-L and DRUDE-L in Table 1). Specifically the polarizable DRUDE force field47−49 was compared with the additive CHARMM potential.36 The DRUDE force field parameters of alkane and octanol molecule have been reported by Ansimov et al.50 and those of the SWM4-NDP water model by Lamoureux et al.51 Given the large sizes of the two “L” simulations, the greatly scalable package NAMD (version 2.10)52 was employed in double precision on Mira of Argonne Leadership Computing Facility. For each of the two systems, 4096 CPU cores were employed in performing the simulation. In the simulations of CHARMM-L and DRUDE-L, the nonbonded pair lists were searched up to the distance of 13.5 Å. The nonbonded switching function was turned on between 10 and 12 Å. The short-range Coulomb interactions were

systematic kinetic study of interfacial self-assembly process, which demands a sophisticated sampling algorithm, such as local sampling enhancement,34 transition path sampling,35 etc., is our ongoing work.

■

SIMULATION METHODOLOGY Two sets of simulations were carried out to investigate the effects of finite size simulation box and force fields in the water/dodecane interface simulations. For the former, we simulated three systems, which are named as systems CHARMM-s1, CHARMM-s2, and CHARMM-s3 in Table 1. The composition ratio of octanol/water/dodecane molecules is fixed, with the simulation box dimension varied. Each simulation box suffixed with s1, s2, or s3 contains around 0.05 million atoms. In these three simulations, the CHARMM General Force Field (CGenFF),36 which has already been implemented in the MD package GROMACS,37 was employed for all the molecules with no modifications. The CHARMM TIP3P water model was used with the structures constrained using the SETTLE algorithm.38 GROMACS version 4.5.539 was employed for all the finite size effect investigations. The force field parameters are also presented in the Supporting Information for completion. The initial structures were generated with the package Packmol.40 The octanol molecules were randomly distributed in the whole simulation box to eliminate any possible biased effect on the self-assembly process. Initially water and dodecane molecules formed separate water phase and oil phase, respectively (Figure 1a). The energy of the initial structure was first minimized using the steepest descent algorithm, followed by 1 ns pre-equilibration under the condition of constant pressure and temperature (NPT). Subsequent production runs were performed under NPT condition, with the temperature (298 K) of octanol, water and

Figure 1. (a) Snapshots of the octanol (highlighted) self-assembly at different simulation times and (b) noncovalent interactions (Coulomb and vdW) between octanol and different molecules in the system CHARMM-s1. 688

DOI: 10.1021/acs.jpcc.7b10997 J. Phys. Chem. C 2018, 122, 687−693

Article

The Journal of Physical Chemistry C

Figure 2. Last snapshots (200 ns) in systems (a) CHARMM-s1, (b) CHARMM-s2, and (c) CHARMM-s3. (d) Density profiles of the octanol oxygen atoms as a function of their coordinates along the Z dimension. z = 0 denotes the center of the oil phase on the basis of the dodecane carbon atoms. The density profiles of water and dodecane are provided in Figure S2.

with the snapshots in Figure 1a. Within the first 10 ns simulation all the energy profiles changed sharply, and gradually converged later on. Impressively, observable changes occurred at simulation time of ∼110−120 ns, due to the “translocation” of those octanol-bearing micelles from the water phase to the oil (Figure 1a). Moreover, these calculations indicate that the vdW interactions dominate over the Coulomb interactions between octanol molecules, and between octanol and dodecane, whereas the Coulomb interactions dominate between octanol and water. These simulated results are in accordance with the fact that the Coulomb interactions play a more crucial role when polar molecules (water) are involved, whereas the vdW interactions are more significant for the nonpolar (dodecane) and weak polar (octanol) molecules. To investigate finite size effect, the systematically increased lateral areas (X × Y) for the octanol disbribution at the water/ oil interface are adopted, as shown in Table 1 and Figure 2a−c. To quantify the finite size effect, the density profile of octanol oxygens in the Z dimension is calcualted on the basis of the simulation trajectories of 150−200 ns (Figure 2d). It can be seen that varying box lengths in XY plane greatly affects the density profiles. First the peak height gradually decreases with increasing box length in XY plane. A striking finding is that in the system CHARMM-s1, a second aggregation layer at |z| ≈ 3.4 nm is stably formed beneath the primary water/dodecane interfacial layer at |z| ≈ 5.5 nm. This secondary peak becomes weaker, however still observable, in the system CHARMM-s2, and disppears in the system CHARMM-s3 with the largest lateral area in XY plane. 2. Effects of Atomistic Force Fields. The above solid evidence of how finite size effect induces misguiding interfacial behavior endorses the necessity of constructing large simulation box. In addition to the elimination of finite size effect, we ran two parallel simulations using the pairwise additive CHARMM and the polarizable DRUDE potentials. The initial and final structures are illustrated in Figure 3a and b, respectively, for the CHARMM-L system. To quantify the convergence of the simulations, we calculated the intermolecular interactions between octanol and different types of molecules in the CHARMM simulation. Figure 3c shows that the simulation reached equilibration at a time scale of 100 ns. Similar to Figure 1b, the electrostatic interactions play a predominant role in the interactions between water molecules and octanol, while the vdW interactions dominate the

calculated up to a distance of 12 Å, with the long-range interactions calculated using the PME algorithm43,44 under the three-dimensional PBC conditions. The six-order interpolation was utilized with the tolerance of 10−6 and the grid spacing of 1 Å for the Ewald algorithm. The isothermal−isobaric ensemble was employed. The Langevin dynamics was used at the reference temperature of 298.15 K and the damping coefficient of 5 ps−1. The integration time step was 2 fs with all the covalent bond lengths constrained. The impulse multiple time step algorithm53 was utilized, where the nonbonded interactions were calculated every one time step along with every two time steps for the long-range electrostatic interactions. The Nosé−Hoover Langevin piston semiisotropic pressure coupling54,55 was employed to reach the optimal system density with PZ = 1 atm (damping time scale 50 fs and oscillation period 100 fs) while fixing the simulation box edge length in X and Y dimensions. This simulation was performed for 200 ns. Because of the large storage requirement, this simulation trajectory was saved at a frequency of 1 frame per ns. In the simulation on the system DRUDE-L, an integration time step of 1 fs was employed, and the nonbonded interactions were calculated every time step. The other simulation control parameters were the same as those in the CHARMM-L simulation.

■

RESULTS AND DISCUSSION 1. Finite Size Effects. In the small size simulations, as illustrated in Figure 1a with the system CHARMM-s1, the octanol molecules rapidly aggregated from random distribution into large micelles in aqueous phase in a few ns time scale. As expected, the weak polar headgroups of octanol molecules extend outward and are stabilized by H-bonds with surrounding water environment. These micelles sojourn in water phase for 100 ns or so before they diffuse to the organic phase. These findings are consistent with the low solubility of octanol in water (540 mg/L at 298 K56). In the oil phase, the octanol molecules form polar-nonpolar self-assembly through H-bond interactions with neighboring octanol and water molecules and nonpolar interactions with dodecane molecule (Figure S1).57 Figure 1b exhibits the short-range noncovalent interactions between octanol and different molecules as a function of simulation time. In these calculations all the Coulomb and vdW interactions are calculated up to the atom− atom distance of 1.2 nm. The energy profiles are consistent 689

DOI: 10.1021/acs.jpcc.7b10997 J. Phys. Chem. C 2018, 122, 687−693

Article

The Journal of Physical Chemistry C

3. Structures at Liquid/Liquid Interface. All the simulations evidence that the convergence of the interfacial adsorption of octanol molecules occurs at a time scale of around 100 ns. Compared with many biomolecule simulations of order of microseconds on GPUs, 100 ns time scale on CPU seems not long enough. However, it needs to be stressed that the current simulation involves over three hundreds of thousands atoms with full range of interactions. For such systems current serial computation on GPU does not show any time-to-solution advantage over massively parallel computing on large CPU cluster. Double-precision and large pairwise cutoff distance on CPU, which guarantees energy accuracy and conservation, is also an important factor to consider. On the other side, the simulated problem here is dominated by weak electrostatic interactions and vdW interactions. Compared to complex biophysical processes such as protein folding/ unfolding, this problem does not possess any comparable complicated conformation sampling and kinetics and therefore shows no necessity of microseconds or longer simulation. Note that to speed up the aggregation behavior, simulated annealing technique has been suggested recently by us5,14,57−60 and others.61 With the simple sampling acceleration technique, simulated cosurfactant/surfactant self-assembly is observed occurring at time scale of even smaller than 100 ns.5 Retrieved from the trajectory of 160−200 ns, the strengths of the interactions in the whole systems roughly follow the order: vdW (octanol-dodecane) > vdW (octanol−octanol) ≈ Coulomb (octanol−water) > vdW (octanol−water) > Coulomb (octanol-dodecane) > Coulomb (octanol−octanol). The vdW interactions dominate between molecules of octanols, and between octanol and dodecane, which are ascribed to the weak polar feature of octanol and the nonpolar feature of dodecane molecules. With a high solubility of octanol in oil phase, strong vdW interactions are presented between octanol and dodecane and independent of the system size. As expected, the electrostatic interactions play a critical role in the interactions between octanol and water. From the point of view of free energy, the octanol layer plays a role of reducing free energy barrier between water and oil phases, driven by electrostatic interactions with water and vdW interactions with dodecane. The two interactions tend to align octanols at the interface. As a result of such an ordered orientation, the presence of positive Coulomb interactions between octanol molecules is observed, demonstrated in Figure 1b and Figure 3c. The ordered octanol layer indicates decrease of entropy, which is effectively compensated by the more negative interaction energy with water phase and oil phase to reduce the overall free energy. At any rate, the ordering of octanols is the critical factor for surface tension reduction. In this regard, we calculated the second order orientational parameters É ÄÅ Å3 1 1 ÑÑ O(z) = ∑ ÅÅÅÅÅ cos2 α(z) − ÑÑÑÑÑ N (z ) 2Ö (1) Ç2

Figure 3. Snapshots of (a) the initial and (b) the final (200 ns) structures of the system CHARMM-L. All octanol molecules are dissolved in the dodecane phase. (c) Convergence of interaction energies between octanol and different types of molecules in the system CHARMM-L.

interactions between octanol and octanol, and between octanol and dodecane. In addition to the noncovalent interactions in Figure 3c, we calculated the average lateral area per alcohol as a function of simulation time (Figure S3), which also evidence the convergence of the simulation at around 100 ns. The last simulation snapshots for the CHARMM and DRUDE force fields are presented in Figure 4a and b,

Figure 4. Last snapshots of (a) system CAHRMM-L and (b) system DRUDE-L. The number density profiles of the octanol oxygen atoms are inserted as blue and yellow curves, respectively. The standard deviations are presented as the error bars of the density profiles.

respectively. The number density profiles of octanol oxygen atoms with different force fields were calculated and provided in the insets of Figure 4, parts a and b. These calculations were performed with the trajectory generated during the simulation period of 160−200 ns. As shown in Figure 4, the calculated number density profiles resemble each other with no observable distribution difference. Indeed, as octanol molecule only has a weak polar headgroup interacting with water phase, it is natural that polarization effect plays a negligible role in the current study.

where α(z) denotes the angle between the Z direction and the vector pointing from octanol oxygen to the terminal carbon on the same molecule (Figure 5b), N(z) is the number of octanol oxygen atoms distributed in a slab between (z − 0.05) nm and (z + 0.05) nm with a bin width of 0.1 nm. Within this context a value of unit of O(z) stands for that octanol chains exactly parallel the Z direction, O(z) = −0.5 for perfect vertical orientation to Z direction, and O(z) = 0 for random 690

DOI: 10.1021/acs.jpcc.7b10997 J. Phys. Chem. C 2018, 122, 687−693

Article

The Journal of Physical Chemistry C

simulation results for octanol adsorption at water/oil interface. It is shown that well converged interfacial structures require simulation time scale of 100 ns or so with all-atom MD simulation. Effective elimination of finite size effect demands large system dimension and long MD trajectory with massively parallel computing, a challenging task for the majority of computational researchers. The simulated interfacial aggregation behavior with additive CHARMM potential and more accurate DRUDE polarizable model does not show considerable difference, due to the weak polar nature of octanol molecules. The calculated interfacial area of octanol at a water−dodecane interface, ∼0.30 nm2 per octanol molecule, is in reasonable agreement with the experiment value of 0.24 nm2 at water/dodecane interface. As a summary, by rigorously quantifying the convergence of the liquid−liquid simulations with high fidelity all-atom force fields and vanishing finite size effect, we expect that this work will provide an MD simulation benchmark for complex liquid/liquid interface with cosurfactants reducing interfacial tension between two immiscible phases.

Figure 5. (a) Second order orientational parameter of octanol molecule and number density profile of octanol oxygen in the system CHARMM-L. Error bars (standard deviations) are included for the orientational parameters. (b) Definition of the angle α in eq 1, where Z represents the normal vector of the water/dodecane interface.

orientation. The calculated result is shown in Figure 5a for the system CHARMM-L. Strong orientation of the octanol molecules was demonstrated at the interfacial regimes of 9 < |z| < 11 nm, where the two primary orientation peaks of octanol molecules can be clearly observed. The two orientation peaks indicate a tilt angle of approximately 35° of the octanol molecules at the interface. It needs to be noted that due to the limited amount of octanols in the simulation box and self-assembly of octanols at water/oil interface, the overall concentration of octanol in oil is different from the “effective” concentration of octanols in the bulk dodecane phase. The amount of octanol molecules at water/dodecane interface is in equilibration with that in the bulk oil phase. That is, the distribution of interfacial octanol molecules is affected by both the lateral area and the “effective” concentration of octanol molecules in the oil phase. The effective bulk concentration of octanol in the simulations is estimated from the number of octanols distributed in the central region of the simulation box and the corresponding volume of the central region. Along z direction, the central region is defined as a slab between the upper position and the lower position of the minimum density (right adjacent to the primary distribution peaks) of octanol oxygen. In this context, the thickness of the central region is 18 nm in the CHARMML system (−9 < z < 9 nm, inset of Figure 4). The effective bulk concentrations are thus estimated as 1.1, 0.3, 0.03, 1.6, and 1.7 M for systems CHARMM-s1, CHARMM-s2, CHARMM-s3, CHARMM-L, and DURDE-L systems, respectively. The interfacial area per octanol molecule is calculated as 0.30 and 0.35 nm 2 for CHARMM-L and DRUDE-L systems, respectively (0.27, 0.28, and 0.35 nm2 for CHARMM-s1, CHARMM-s2, and CHARMM-s3, respectively). These values are in reasonable agreement with the experimental value of 0.24 nm2 per octanol at the dodecane/water interface.62 They also agree well with the reported values of 0.18−0.27 nm2 per octanol at water−air interfaces by a variety of experimental approaches63 (around 0.3 nm2 per hexadecanol at water−air interfaces64). From the density profiles of octanol oxygen atoms, the calculated full width at half-maximum are ∼0.6 nm, indicating a sharp interface in the systems investigated.65

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcc.7b10997.

■

Octanol aggregates in oil phase, density profiles in the small systems, the interfacial area per octanol as a function of the simulation time, and CHARMM and DRUDE force field parameters (PDF)

AUTHOR INFORMATION

Corresponding Authors

*(B.Q.) E-mail: [email protected]. *(W.J.) E-mail: [email protected]. ORCID

Baofu Qiao: 0000-0001-8870-5985 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work is supported by the U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences, Division of Chemical Sciences, Biosciences and Geosciences, under Contract DE-AC02-06CH11357. The computing resources, provided on Blues, a high-performance computing cluster, operated by the Laboratory Computing Resource Center at Argonne National Laboratory, are gratefully acknowledged. This research used resources of the Argonne Leadership Computing Facility, which is a DOE Office of Science User Facility supported under Contract DE-AC02-06CH11357.

■

REFERENCES

(1) Boker, A.; He, J.; Emrick, T.; Russell, T. P. Self-Assembly of Nanoparticles at Interfaces. Soft Matter 2007, 3, 1231−1248. (2) Shao, J. J.; Lv, W.; Yang, Q. H. Self-Assembly of Graphene Oxide at Interfaces. Adv. Mater. 2014, 26, 5586−5612. (3) Isa, L.; Kumar, K.; Müller, M.; Grolig, J.; Textor, M.; Reimhult, E. Particle Lithography from Colloidal Self-Assembly at Liquid− Liquid Interfaces. ACS Nano 2010, 4, 5665−5670. (4) Scoppola, E.; Watkins, E. B.; Campbell, R. A.; Konovalov, O.; Girard, L.; Dufrêche, J.-F.; Ferru, G.; Fragneto, G.; Diat, O. Solvent

■

CONCLUSIONS By systematically increasing the atomistic simulation box size and comparing the additive CHARMM force field and the polarizable DRUDE potential, we have obtained converged 691

DOI: 10.1021/acs.jpcc.7b10997 J. Phys. Chem. C 2018, 122, 687−693

Article

The Journal of Physical Chemistry C Extraction: Structure of the Liquid−Liquid Interface Containing a Diamide Ligand. Angew. Chem., Int. Ed. 2016, 55, 9326−9330. (5) Qiao, B.; Muntean, J. V.; Olvera de la Cruz, M.; Ellis, R. J. Ion Transport Mechanisms in Liquid−Liquid Interface. Langmuir 2017, 33, 6135−6142. (6) Benjamin, I. Molecular Structure and Dynamics at Liquid-Liquid Interfaces. Annu. Rev. Phys. Chem. 1997, 48, 407−451. (7) Luo, G.; Malkova, S.; Yoon, J.; Schultz, D. G.; Lin, B.; Meron, M.; Benjamin, I.; Vanýsek, P.; Schlossman, M. L. Ion Distributions Near a Liquid-Liquid Interface. Science 2006, 311, 216−218. (8) Ariga, K.; Nakanishi, T.; Hill, J. P.; Shirai, M.; Okuno, M.; Abe, T.; Kikuchi, J.-i. Tunable pK of Amino Acid Residues at the Air− Water Interface Gives an L-zyme (Langmuir Enzyme). J. Am. Chem. Soc. 2005, 127, 12074−12080. (9) Schlossman, M. L. Liquid−Liquid Interfaces: Studied by X-Ray and Neutron Scattering. Curr. Opin. Colloid Interface Sci. 2002, 7, 235−243. (10) Shen, Y. R. Surface Properties Probed by Second-Harmonic and Sum-Frequency Generation. Nature 1989, 337, 519−525. (11) Shen, Y. R. Basic Theory of Surface Sum-Frequency Generation. J. Phys. Chem. C 2012, 116, 15505−15509. (12) Johnson, C. M.; Baldelli, S. Vibrational Sum Frequency Spectroscopy Studies of the Influence of Solutes and Phospholipids at Vapor/Water Interfaces Relevant to Biological and Environmental Systems. Chem. Rev. 2014, 114, 8416−8446. (13) Leung, C.-Y.; et al. Molecular Crystallization Controlled by pH Regulates Mesoscopic Membrane Morphology. ACS Nano 2012, 6, 10901−10909. (14) Qiao, B.; Demars, T.; Olvera de la Cruz, M.; Ellis, R. J. How Hydrogen Bonds Affect the Growth of Reverse Micelles around Coordinating Metal Ions. J. Phys. Chem. Lett. 2014, 5, 1440−1444. (15) Ferru, G.; Gomes Rodrigues, D.; Berthon, L.; Diat, O.; Bauduin, P.; Guilbaud, P. Elucidation of the Structure of Organic Solutions in Solvent Extraction by Combining Molecular Dynamics and X-ray Scattering. Angew. Chem., Int. Ed. 2014, 53, 5346−5350. (16) Ito, M.; Cosgrove, T. Atomistic Simulation of Short-Chain Polymers at the Liquid/Liquid Interface. Colloids Surf., A 1994, 86, 125−131. (17) Luo, M.; Mazyar, O. A.; Zhu, Q.; Vaughn, M. W.; Hase, W. L.; Dai, L. L. Molecular Dynamics Simulation of Nanoparticle SelfAssembly at a Liquid−Liquid Interface. Langmuir 2006, 22, 6385− 6390. (18) Hantal, G.; Terleczky, P.; Horvai, G.; Nyulászi, L.; Jedlovszky, P. Molecular Level Properties of the Water−Dichloromethane Liquid/Liquid Interface, as Seen from Molecular Dynamics Simulation and Identification of Truly Interfacial Molecules Analysis. J. Phys. Chem. C 2009, 113, 19263−19276. (19) Vácha, R.; Rick, S. W.; Jungwirth, P.; de Beer, A. G. F.; de Aguiar, H. B.; Samson, J.-S.; Roke, S. The Orientation and Charge of Water at the Hydrophobic Oil Droplet−Water Interface. J. Am. Chem. Soc. 2011, 133, 10204−10210. (20) Vazdar, M.; Pluhařová, E.; Mason, P. E.; Vácha, R.; Jungwirth, P. Ions at Hydrophobic Aqueous Interfaces: Molecular Dynamics with Effective Polarization. J. Phys. Chem. Lett. 2012, 3, 2087−2091. (21) Harwood, D. B.; Peters, C. J.; Siepmann, J. I. A Monte Carlo Simulation Study of the Liquid−Liquid Equilibria for Binary Dodecane/Ethanol and Ternary Dodecane/Ethanol/Water Mixtures. Fluid Phase Equilib. 2016, 407, 269−279. (22) Fukuto, M.; Ocko, B. M.; Bonthuis, D. J.; Netz, R. R.; Steinruck, H. G.; Pontoni, D.; Kuzmenko, I.; Haddad, J.; Deutsch, M. Nanoscale Structure of the Oil-Water Interface. Phys. Rev. Lett. 2016, 117, 256102. (23) Zhang, J.; Dong, Z.; Zhang, Y.; Wang, M.; Yan, Y. Effects of the Methane Content on the Water−Oil Interface: Insights from the Molecular Level. Energy Fuels 2017, 31, 7026−7032. (24) Tan, J. S. J.; Zhang, L.; Lim, F. C. H.; Cheong, D. W. Interfacial Properties and Monolayer Collapse of Alkyl Benzenesulfonate Surfactant Monolayers at the Decane−Water Interface from Molecular Dynamics Simulations. Langmuir 2017, 33, 4461−4476.

(25) Neyt, J.-C.; Wender, A.; Lachet, V.; Ghoufi, A.; Malfreyt, P. Quantitative Predictions of the Interfacial Tensions of Liquid−Liquid Interfaces through Atomistic and Coarse Grained Models. J. Chem. Theory Comput. 2014, 10, 1887−1899. (26) Ghoufi, A.; Malfreyt, P.; Tildesley, D. J. Computer Modelling of the Surface Tension of the Gas-Liquid and Liquid-Liquid Interface. Chem. Soc. Rev. 2016, 45, 1387−1409. (27) Sresht, V.; Lewandowski, E. P.; Blankschtein, D.; Jusufi, A. Combined Molecular Dynamics Simulation−Molecular-Thermodynamic Theory Framework for Predicting Surface Tensions. Langmuir 2017, 33, 8319−8329. (28) Katiyar, P.; Singh, J. K. A Coarse-Grain Molecular Dynamics Study of Oil−Water Interfaces in the Presence of Silica Nanoparticles and Nonionic Surfactants. J. Chem. Phys. 2017, 146, 204702. (29) Schoettl, S.; Marcus, J.; Diat, O.; Touraud, D.; Kunz, W.; Zemb, T.; Horinek, D. Emergence of Surfactant-Free Micelles from Ternary Solutions. Chem. Sci. 2014, 5, 2949−2954. (30) Zemb, T.; Kunz, W. Weak Aggregation: State of the Art, Expectations and Open Questions. Curr. Opin. Colloid Interface Sci. 2016, 22, 113−119. (31) Zemb, T. N.; et al. How to Explain Microemulsions Formed by Solvent Mixtures Without Conventional Surfactants. Proc. Natl. Acad. Sci. U. S. A. 2016, 113, 4260−4265. (32) Yan, T.; Li, S.; Jiang, W.; Gao, X.; Xiang, B.; Voth, G. A. Structure of the Liquid−Vacuum Interface of Room-Temperature Ionic Liquids: A Molecular Dynamics Study. J. Phys. Chem. B 2006, 110, 1800−1806. (33) Jiang, W.; Hardy, D. J.; Phillips, J. C.; MacKerell, A. D.; Schulten, K.; Roux, B. High-Performance Scalable Molecular Dynamics Simulations of a Polarizable Force Field Based on Classical Drude Oscillators in NAMD. J. Phys. Chem. Lett. 2011, 2, 87−92. (34) Jo, S.; Jiang, W. A Generic Implementation of Replica Exchange with Solute Tempering (REST2) Algorithm in NAMD for Complex Biophysical Simulations. Comput. Phys. Commun. 2015, 197, 304− 311. (35) Bolhuis, P. G.; Chandler, D.; Dellago, C.; Geissler, P. L. Transition Path Sampling: Throwing Ropes Over Rough Mountain Passes, in the Dark. Annu. Rev. Phys. Chem. 2002, 53, 291−318. (36) Vanommeslaeghe, K.; et al. CHARMM General Force Field: A Force Field for Drug-Like Molecules Compatible with the CHARMM All-Atom Additive Biological Force Fields. J. Comput. Chem. 2010, 31, 671−690. (37) Bjelkmar, P.; Larsson, P.; Cuendet, M. A.; Hess, B.; Lindahl, E. Implementation of the CHARMM Force Field in GROMACS: Analysis of Protein Stability Effects from Correction Maps, Virtual Interaction Sites, and Water Models. J. Chem. Theory Comput. 2010, 6, 459−466. (38) Miyamoto, S.; Kollman, P. A. SETTLE: An Analytical Version of the SHAKE and RATTLE Algorithm for Rigid Water Models. J. Comput. Chem. 1992, 13, 952−962. (39) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 435−447. (40) Martínez, L.; Andrade, R.; Birgin, E. G.; Martínez, J. M. PACKMOL: A Package for Building Initial Configurations for Molecular Dynamics Simulations. J. Comput. Chem. 2009, 30, 2157−2164. (41) Cheng, A.; Merz, K. M. Application of the Nosé−Hoover Chain Algorithm to the Study of Protein Dynamics. J. Phys. Chem. 1996, 100, 1927−1937. (42) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. Appl. Phys. 1981, 52, 7182−7190. (43) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An N· log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. 692

DOI: 10.1021/acs.jpcc.7b10997 J. Phys. Chem. C 2018, 122, 687−693

Article

The Journal of Physical Chemistry C

(64) Menger, F. M.; Shi, L.; Rizvi, S. A. A. Re-evaluating the Gibbs Analysis of Surface Tension at the Air/Water Interface. J. Am. Chem. Soc. 2009, 131, 10380−10381. (65) Luo, G.; Malkova, S.; Pingali, S. V.; Schultz, D. G.; Lin, B.; Meron, M.; Graber, T. J.; Gebhardt, J.; Vanysek, P.; Schlossman, M. L. The Width of the Water/2-Heptanone Liquid−Liquid Interface. Electrochem. Commun. 2005, 7, 627−630.

(44) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577−93. (45) Hess, B. P-LINCS: A Parallel Linear Constraint Solver for Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 116−122. (46) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. LINCS: A Linear Constraint Solver for Molecular Simulations. J. Comput. Chem. 1997, 18, 1463−1472. (47) Lamoureux, G.; Harder, E.; Vorobyov, I. V.; Roux, B.; MacKerell, A. D., Jr A Polarizable Model of Water for Molecular Dynamics Simulations of Biomolecules. Chem. Phys. Lett. 2006, 418, 245−249. (48) Schröder, C.; Steinhauser, O. Simulating Polarizable Molecular Ionic Liquids with Drude Oscillators. J. Chem. Phys. 2010, 133, 154511−154513. (49) Huang, J.; Lopes, P. E. M.; Roux, B.; MacKerell, A. D. Recent Advances in Polarizable Force Fields for Macromolecules: Microsecond Simulations of Proteins Using the Classical Drude Oscillator Model. J. Phys. Chem. Lett. 2014, 5, 3144−3150. (50) Anisimov, V. M.; Vorobyov, I. V.; Roux, B.; MacKerell, A. D., Jr. Polarizable Empirical Force Field for the Primary and Secondary Alcohol Series Based on the Classical Drude Model. J. Chem. Theory Comput. 2007, 3, 1927−1946. (51) Lamoureux, G.; MacKerell, A. D., Jr.; Roux, B. A Simple Polarizable Model of Water Based on Classical Drude Oscillators. J. Chem. Phys. 2003, 119, 5185−5197. (52) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kale, L.; Schulten, K. Scalable Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781− 1802. (53) Bishop, T. C.; Skeel, R. D.; Schulten, K. Difficulties with Multiple Time Stepping and Fast Multipole Algorithm in Molecular Dynamics. J. Comput. Chem. 1997, 18, 1785−1791. (54) Martyna, G. J.; Tobias, D. J.; Klein, M. L. Constant Pressure Molecular Dynamics Algorithms. J. Chem. Phys. 1994, 101, 4177− 4189. (55) Feller, S. E.; Zhang, Y.; Pastor, R. W.; Brooks, B. R. Constant Pressure Molecular Dynamics Simulation: The Langevin Piston Method. J. Chem. Phys. 1995, 103, 4613−4621. (56) Yaws, C. L.; Hopper, J. R.; Sheth, S. D.; Han, M.; Pike, R. W. Solubility and Henry’s Law Constant for Alcohols in Water. Waste Manage. 1998, 17, 541−547. (57) Knight, A. W.; Qiao, B.; Chiarizia, R.; Ferru, G.; Forbes, T.; Ellis, R. J.; Soderholm, L. Subtle Effects of Aliphatic Alcohol Structure on Water Extraction and Solute Aggregation in Biphasic Water/nDodecane. Langmuir 2017, 33, 3776−3786. (58) Qiao, B.; Ferru, G.; Olvera de la Cruz, M.; Ellis, R. J. Molecular Origins of Mesoscale Ordering in a Metalloamphiphile Phase. ACS Cent. Sci. 2015, 1, 493−503. (59) Bera, M. K.; Qiao, B.; Seifert, S.; Burton-Pye, B. P.; Olvera de la Cruz, M.; Antonio, M. R. Aggregation of Heteropolyanions in Aqueous Solutions Exhibiting Short-Range Attractions and LongRange Repulsions. J. Phys. Chem. C 2016, 120, 1317−1327. (60) Qiao, B.; Ferru, G.; Ellis, R. J. Complexation Enhancement Drives Water-to-Oil Ion Transport: A Simulation Study. Chem. - Eur. J. 2017, 23, 427−436. (61) Vo, Q. N.; Hawkins, C. A.; Dang, L. X.; Nilsson, M.; Nguyen, H. D. Computational Study of Molecular Structure and SelfAssociation of Tri-n-butyl Phosphates in n-Dodecane. J. Phys. Chem. B 2015, 119, 1588−1597. (62) Aveyard, R.; Briscoe, B. J. Adsorption of n-Octanol at the nDodecane/Water Interface. Trans. Faraday Soc. 1970, 66, 2911− 2916. (63) Kinoshita, K.; Parra, E.; Needham, D. New Sensitive MicroMeasurements of Dynamic Surface Tension and Diffusion Coefficients: Validated and Tested for the Adsorption of 1-Octanol at A Microscopic Air-Water Interface and its Dissolution into Water. J. Colloid Interface Sci. 2017, 488, 166−179. 693

DOI: 10.1021/acs.jpcc.7b10997 J. Phys. Chem. C 2018, 122, 687−693