Ultra-Coarse-Grained Models Allow for an Accurate and Transferable

Subscriber access provided by UNIV OF NEW ENGLAND ARMIDALE

Article

Ultra-Coarse-Grained Models Allow for an Accurate and Transferable Treatment of Interfacial Systems Jaehyeok Jin, and Gregory A. Voth J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.7b01173 • Publication Date (Web): 26 Feb 2018 Downloaded from http://pubs.acs.org on March 13, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Theory and Computation is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Ultra-Coarse-Grained Models Allow for an Accurate and Transferable Treatment of Interfacial Systems Jaehyeok Jin and Gregory A. Voth1 Department of Chemistry, James Franck Institute, and Institute for Biophysical Dynamics The University of Chicago, Chicago, Illinois, 60637 USA Abstract Interfacial systems are fundamentally important in many processes. However, constructing coarse-grained (CG) models for such systems is a significant challenge due to their inhomogeneous nature. This problem is made worse due to the generally non-transferable nature of the interactions in CG models across different phases. In this paper, we address these challenges by systematically constructing ultra-coarse-grained (UCG) models for interfaces, in which the CG sites are allowed to have internal states. We find that a multiscale coarse-grained (MS-CG) representation of a single CG site model fails to identify the directionality of a molecule and is unable to reproduce the correct phase coexistence for aspherical molecules. In contrast with conventional MS-CG models, the UCG methodology allows chemical and environmental changes to be captured by modulating the interactions between internal states. In this work, we design the internal states to depend on local particle density in order to distinguish different phases in liquid/vapor or liquid/liquid interfaces. These UCG models are able to capture phase coexistence and recapitulate structures, notably at state points in which the MS-CG method yields poor results. Interestingly, effective pairwise forces and potentials from the UCG models are almost identical to those of the bulk liquids that correspond to each phase, indicating that the UCG approach can provides transferable interactions. This approach is expected to be applicable to other systems that exhibit phase coexistence, and also to complex macromolecular systems, by modulating interactions based on local density or other order parameters in order to unravel the complex nature underlying heterogeneous system boundaries.

1

To whom correspondence should be addressed; E-mail: [email protected]

1

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1! Introduction Interfacial phenomena have attracted significant interest over many decades due to their prevalence in nature and importance in various scientific fields.1 Practical applications range from industrial separation processes, including liquid extraction2 and charge transfer,3 to multiphasic catalytic reactions4, 5 and absorption processes.6 In addition to the practical aspects, further studies in complex biomolecular systems,7 such as micelles, lipids, and membranes betters the understanding of liquid/vapor and liquid/liquid interfaces that are of fundamental importance in many scientific fields. Computer simulations can offer valuable insights for such chemical and physical phenomena.8-13 However, studying complex liquids and biological processes using conventional atomistic molecular dynamics (MD) simulations is computationally expensive and sometimes it is not feasible to access the physically relevant time and length scales. In this regard, several coarse-graining (CG) approaches14-21 have been developed by averaging over atomistic details in some way and then creating CG models that can reproduce such essential properties. There are two different ways of constructing coarse-grained force fields (FFs): the “bottom-up” and “top-down” approaches. In the top-down approach, the main objective is to reproduce experimental thermodynamic properties.22-24 Certain top-down approaches such as the MARTINI model21 can be used to model various liquid state molecules as well as more complex lipid bilayers and even biomolecules in some instances. However, these empirical CG models do not have strong connections with certain fine-grained (FG) models, indicating that the models lack detail at a molecular level and do not clearly correspond to a physical system at the atomistic resolution.25 In the bottom-up approach, the CG interactions are extracted and constructed from reference all-atom systems. Employing iterative Boltzmann inversion (IBI),26 inverse Monte Carlo (IMC),27 force matching (FM),28-31 and relative entropy minimization32-35 provides a systematic approach to construct the effective coarse-grained interactions. However, the CG FF

2


Page 2 of 47

Page 3 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


from bottom-up approaches is usually not transferable; developing models that retain both accuracy and transferability remains a challenge.19, 25, 36 In other words, one might obtain a CG interaction of liquid/vapor or liquid/liquid interfaces by employing bottom-up coarse-graining strategies, but the resultant CG interaction may not be transferable to bulk liquid states. Namely, transferable CG models of liquid/vapor or liquid/liquid interfaces are needed to account for different molecular environments in interfacial systems. An important feature of interfacial systems, such as liquid/vapor or liquid/liquid interfaces, is that two different phases in the system are identifiable from the density profile along the axis normal to the surface (slab axis). The liquid/vapor interface has a liquid phase with high density, a vapor phase with low density, and a phase boundary where the density changes monotonically between the two phases. Therefore, instead of trying to describe an interfacial system using just one CG site or “bead” per molecule and a pairwise interaction between these CG molecules (which can work well for a single phase system), one might envision that the CG site has two different possible states: a liquid state which is a locally more dense state in the liquid phase and a locally less dense state in the vapor phase or at the phase boundary. Also, for a liquid/liquid interface, the density profile might be decomposed into two distinct density profiles. In turn, one could classify the state of each CG bead based on its local density. These internal CG states could be determined from the particle’s coordinates relative to the density profile, but such an approach based only on the Eulerian coordinate fails to account for the endogenous nature of changes in the interaction between different phases or environments across the interface. Instead, the local density of the CG particle could be used as a collective variable or order parameter to modulate the interactions between CG particles.37-41 However, this necessitates that the CG model considers changes in structures or chemical environments within the CG beads at its target resolution. The Ultra-Coarse-Graining (UCG) methodology has recently been developed to include both configurational variables and discrete state variables within CG beads as part of a systematic bottom-up coarse-graining approach.42-44 The first paper of the UCG series42 3



introduced the general framework of UCG theory and the statistical procedure for deriving UCG methods related to force matching or relative entropy methods. The second43 and third44 papers focused on algorithms and implementation of UCG models with different assumed dynamics of the internal CG states. The second article introduced the UCG model with internal states where states are in the infrequent transition limit. The third article (denoted as Paper III hereafter), the most relevant for this manuscript, developed the theory and implementation for the opposite limit where the internal UCG states are always in the quasi-equilibrium limit. The UCG with rapid local equilibrium, or UCG-RLE model showed good performance in reproducing cooperative solvophobic association by modulating the solvent interactions based on the local coordination number. In the UCG-RLE framework, each CG site has internal states that rapidly equilibrate, and thus the quasi-equilibrium distribution expressed as the substate probability governs the state changes and further provides the effective CG interactions. Similar to the previously reported CG force fields that depend directly on collective variables or order parameters,37-41 the substate probability can be designed to be dependent on local environmental factors, such as the local density. In our previous work (Paper III), the UCG-RLE methodology was applied to the liquid/vapor phase boundaries as a proof-of-concept. By designing the probability of internal states as a function of the local coordination number, we have shown that UCG-RLE can modulate the phase-dependent substates. However, the application in Paper III had two major drawbacks: 1) the study was restricted to idealized systems governed by Lennard-Jones interactions and did not provide a systematic way of modulating interactions for realistic liquid/vapor systems, and 2) the main objective of employing UCG-RLE was to show that the phase coexistence profile can be suppressed or enhanced by modulating the interaction between internal states. Therefore, interaction parameters used in the UCG-RLE system are not directly obtained from realistic LJ liquid/vapor interfaces but designed to impose internal states where a more dense state has a higher ! and a less dense state has the coexistence ! value. It was a purely model study. 4


Page 4 of 47

Page 5 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60




Page 6 of 47

While several examples regarding the construction of CG models and simulations of CG FFs in liquid/vapor interfaces have been reported, the liquids considered have largely been limited to water,45, 46 methanol,47 or ionic liquids.48-51 However, a wide variety of subtle intermolecular interactions (e.g., hydrogen bonding, dipole-dipole, or induced dipole-induced dipole interactions) can govern the behavior of liquid molecules and our intention here is to extend our UCG approach to capture such subtleties. A recent theoretical investigation at all-atom resolution has suggested that this broad range of intermolecular interactions may affect the interfacial properties including surface tension or tangential pressure.52 From the molecules studied in ref. 52

, we also added chloroform and ethanol in order to investigate the importance of symmetry

within the molecule, which will be discussed later. The chosen set of molecules depicted in Figure 1 contains hydrogen bonding, dipolar aprotic, and van der Waals interactions. For each molecule in Figure 1, we first constructed a liquid system composed of 1000 molecules in a cubic box under periodic boundary conditions. Initial structures of each system were constructed using the VMD53 and Packmol54 packages. The OPLS-AA force field55, 56 was used to model the system with Ewald summation to include long-range electrostatics.57 All calculations were performed with the LAMMPS MD simulation package.58-60 After energyminimization, the system was annealed with a Nosé-Hoover thermostat61,

62

to a final

temperature T = 300 K. We then equilibrated the system for 1 ns with a Nosé-Hoover thermostat and Andersen barostat63 under constant NPT conditions at a pressure P = 1 atm with a damping constant of 1 ps. Next, a liquid/vapor interface was created by inserting a vacuum slab with a length of Lz (e.g., the equilibrated length of the simulation domain) in the z direction. A bulk liquid system was also prepared from the equilibrated system in a cubic box after constant NPT dynamics in order to check transferability of the CG interactions. Finally, constant NVT dynamics were performed for both liquid/vapor interface and bulk liquid systems at T = 300 K with a Nosé-Hoover thermostat for 5 ns using a 0.1 ps damping constant.

6


Page 7 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60




Page 8 of 47

constant NVT MD simulations. (c) Final configuration after NVT runs from the generated liquid/liquid structure. An initial configuration of the interface was constructed from Figure 2-(b) by integrating phases. Figure 2(a) shows the separately prepared bulk liquid systems. The overall system is composed of 168 carbon tetrachloride molecules, thus 84 molecules in each of the two simulation domains, and 400 methanol molecules. After the system depicted in Figure 2(b) is propagated using constant NVT MD simulations, we constructed the liquid/liquid interface by integrating the three domains into a single system. Then, energy minimization was performed to relax the jointed phases. For the relaxed interface, the system was equilibrated for 1 ns with constant NVT dynamics at T = 300 K, and followed by constant NPT dynamics at P = 1 atm with same damping constants from the liquid/vapor interface. From the equilibrated system, all-atom trajectories for analysis were collected every 0.25 ps by performing a final NVT simulation for 5 ns at 300 K.

2.2! Multiscale Coarse-Graining (MS-CG) Procedure The multiscale coarse-graining method (MS-CG) proposed by Izvekov and Voth provides a systematic framework of constructing a CG model28-30, 68 based on the following relation.

exp !U CG (R N ) / k BT " # dr n exp ! u(r n ) / k BT $ % M RN (r n ) ! R N

(

)

(

) (

)

(1)

This equation indicates that the ideal equilibrium distribution of the CG site or “beads” should be identical to that projected from its reference atomistic coordinate distribution by virtue of a set of n mapping functions. In Equation 1, u(r ) denotes the interaction energy as a function of the fine-

grained (usually the all-atom resolution) coordinates, U CG (R N ) denotes the effective CG interaction energy followed by the operator M RN : r n ! R N that maps the all-atom configurations r n on to CG configurations R N by integrating out some degrees of freedoms. (It should be

noted that the notation here is shorthand for a product of delta functions, ! (M RN (r n ) " R I ) , one I

8


Page 9 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


for each mapping to the CG site I.) From the MS-CG theory,28-30, 68 it has been shown that minimizing the following residual in eq 2 as a variational least-square problem provides an efficient way of obtaining the effective CG interaction U CG (R N ) by variationally minimizing the 2 residual ! , given by

1 ! = 3N 2

N

#

f I (r n ) " FI (M RN (r n ))

2

(2)

I =1

where f I (r n ) denotes the fine-grained (atomistic) forces mapped onto CG bead I, while

(

)

FI M RN (r n ) = ! "U CG (R N ) "R NI denotes the effective CG force acting on the CG bead I. I

In this work, CG pair potentials were constructed from the force matching procedure as described above. We used a linear combination of sixth-order B-splines to describe the pair potential. The resolution of B-splines was chosen as 0.20 Å. The outer cutoff for each pair potential was 10.0 Å and the inner cutoffs were chosen from the minimum value of each pair. Then, MS-CG simulations of both liquid/vapor and liquid/liquid interfaces were run for 5 ns for each system following constant NVT dynamics at T = 300 K with a Nosé-Hoover thermostat. The initial structure for the CG simulation was from the CG coordinates of the last frame of the preceding all-atom simulations.

2.3! Ultra-Coarse-Graining with Rapid Local Equilibrium The derivation and implementation of UCG-RLE are detailed in Paper III.44 In this section, a brief description of the main concepts in Paper III will be provided. To introduce the internal states, or substates, of N CG beads, we define ! as all possible configurations from the given states. Therefore, ! = {sg } = {s1 ,s2 ,!,sN } , where g indexes the specific substates and sI denotes the state of CG bead I after assuming that every CG bead has the same number of states, which is a reasonable assumption because of the homogeneous composition of the systems in the current study.

9



Page 10 of 47

Furthermore, in order to perform an efficient computation over all possible substates of the CG beads, we assume that the substate assignments ! are statistically independent conditionals on the coordinates of CG beads. In this case, one can express the mixed free energy of the UCG UCG (R N ) , or the effective CG potential, from the internal states in rapid quasiconfiguration U mix

N equilibrium using the conditional probability p(! | R ) of each ! given R N , since the

substate probability is based on UCG configurational variables R N . Here, in the same way as in the previous papers from the UCG series, the unique pair potential between CG bead pairs

U (! ,R N ) with a substate assignment of ! is composed of pair interactions between the state sI of particle I and sJ of particle J with an interparticle distance of RIJ , i.e. ,

"

U (! ,R N ) =

U 2,s ,s (RIJ ) . I

J

I ,J neighs

(

)

UCG U mix (R N ) = ! ! p(sI | R N ) k BT ln p(sI | R N ) ! I

sI

(3)

! ! p(s

! +

N

I

N

| R ) p(sJ | R ) U 2,s ,s (RIJ ) I

I ,J neighs sI ,sJ

J

(It should be noted that there is a sign error on the term !k BT ln p in eqs 5-11 and 14 of Paper III, which is corrected here.) Therefore, the final UCG force expression can be obtained from the derivative of eq 3, as follows UCG !U mix (R N ) = " " k BT ln p(sI | R N ) !p(sI | R N ) !

(

I

)

sI

! +

" " p(s

I

| R N ) p(sJ | R N )!U 2,s ,s (RIJ ) ! I

I ,J neighs sI ,sJ

! +

" "U I ,J neighs sI ,sJ

2,sI ,sJ

J

(RIJ ) #$ p(sI | R N )!p(sJ | R N ) + p(sJ | R N )!p(sI | R N ) %&

It is worth noting that the form of Equation 4 does not confine the pairwise force fields

U 2,s ,s (RIJ ) to any certain form as long as there are well-defined substate probabilities (with I J 10


(4)

Page 11 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


their derivatives) and pair interactions between the substates of CG beads. However, in the previous UCG papers, we have shown that MS-CG theory can efficiently and effectively implement and parameterize pairwise, per-particle per-substate force fields by minimizing the following residual, similar to Equation 2, given by

1 1 * U " n , M (rtn ) + k BT ln p+ " n | M (rtn ) , M + *u(rtn ) ) ) Ts '( T Ns n s t

! 2 #$U (" ,R N ) %& = lim

( (

)

(

))

(

)

2

(5)

In the present work, the pre-sampled CG trajectory with each CG bead’s coordinates and forces from Equation 4 was duplicated by N s = 20 times. For each duplication, the internal states of each CG bead were assigned differently based on its local density employing a LAMMPS58-60 rerun. A detailed description of designing the UCG model will be discussed in Sections 3.2 and 3.5.

3! Results 3.1! Density Profile from All-atom and MS-CG Simulations for Liquid/Vapor Interfaces The density profile of the liquid/vapor interface has been regarded as a fundamental property of the interfacial structures both from theory and experiment.69 From the all-atom reference trajectories, CG trajectories were prepared using the center of mass of each molecule (a one-site CG model) and the slab density profiles ! (z) were calculated from the all-atom reference trajectories. From the CG reference trajectories, MS-CG force fields were obtained by force matching, and slab density profiles were calculated from subsequent simulations. Figure 3 compares the density profiles along the slab axis ( z -axis) from all-atom and MS-CG simulations for methanol, ethanol, acetone, acetonitrile, neopentane, carbon tetrachloride, chloroform, and benzene.

11




Page 12 of 47

Page 13 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the interface are needed to understand and capture the interactions compared to the bulk liquid.71 The calculation of sum-frequency vibrational spectra for the acetone interface also resulted in a structural orientation similar to that of methanol, which demonstrated that this ordering at the phase boundary can originate not only from hydrogen bonding, but also from polar ordering.72 This behavior also supports the notion that the agreement between all-atom and MS-CG simulations was only observed with carbon tetrachloride and neopentane because the spherical symmetry of the CG beads matches the all-atom symmetry. A detailed analysis will be followed below to support our hypothesis. In comparison to a general tangent hyperbolic shape from the slab density profile, the MS-CG simulation of methanol shows a sharper density profile without a steep regime at the liquid phase (see Figure 3(a)). Also, a non-zero binned density within the sparse region around 50 – 75 Å and the density gradient at the phase boundary indicate that the MS-CG model might not be able to correctly reproduce the interactions that occur at the phase boundary. For the ethanol molecule, we expect less agreement since ethanol has one more methyl group than methanol, thus the more linear structure of ethanol than that of methanol deviates even more from the spherical symmetry of a single CG bead. As predicted, the density profile of ethanol using the MS-CG method shows only a homogeneous fluid with a normalized density of 1.0 rather than two phases in coexistence; note the similarity between the all-atom density profiles of the methanol and ethanol systems. A similar homogeneous density profile is also observed in Figure 3(d) for acetonitrile, which has a linear structure formed from a strong triple C ! N bond, indicating that the CG model is not able to capture any interfacial interactions. Acetone has a higher symmetry than methanol or acetonitrile due to its C2 symmetry axis and therefore better agreement between the mapped all-atom and MS-CG density is expected. Figure 3(c) supports this prediction, but the MS-CG model does not fully reproduce the width and binned density of the liquid and vapor phases. For example, the density profile from the mapped all-atom system clearly distinguishes the sparse regimes at nearly 55 – 95 Å with zero binned density from the

13



dense regimes. On the other hand, in the MS-CG model, the lowest binned density is about 0.45 and the overall shape of the density profile is not reproduced. An interesting behavior is seen by comparing carbon tetrachloride to chloroform. Carbon tetrachloride has a highly symmetric (tetrahedral) structure, as described above, and thus shows almost identical density profiles between the mapped all-atom and MS-CG models as described above and shown in Figure 3(f). In comparison, chloroform has one hydrogen atom in place of a chlorine atom and no longer exhibits high symmetry. Due to its asymmetry, the chloroform system exhibits a MS-CG density profile that slightly differs from that of the mapped all-atom profile as seen by the reduced normalized density at the high-density regime and shallower density gradient at the phase boundary. However, since only a few structural symmetries are broken, it does not demonstrate a drastic change such as in the ethanol or acetonitrile cases. Finally, the density profile of benzene confirms this story. Even if benzene is a highly symmetric molecule with a point group of D6h, it has planar rather than spherical symmetry. This planar symmetry along with conjugated ! electrons enables benzene to form strong intermolecular

! " ! stacking interactions, which yields strong directionality along the C6 axis. Next, we investigate the aforementioned liquid molecules, except neopentane and carbon tetrachloride which are already well described by just MS-CG, using the UCG formalism (in a state of local quasi-equilibrium of CG sites, i.e., UCG-RLE) to better understand and improve the interactions at liquid/vapor interfaces.

3.2! UCG Model Construction for Liquid/Vapor Interface Figure 4 depicts a schematic diagram of how the inner and outer regions can be defined on the methanol interface system as an example. To distinguish the two distinct regions, we used the density profile from Figure 4, which shows the two phases separated by a hyperbolic tangent density profile at the phase boundary. The inner region (the red dashed part in the left side of Figure 4) is defined as the region where the CG particles are distributed with an approximately

14


Page 14 of 47

Page 15 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60




Separating the interface into two distinct regions provides a reliable way to calculate the radial distribution function (RDF or g(r)) on the interface with correct asymptotic behavior. Unlike in bulk liquids, the RDF defined within the interface encounters the following problem. By inserting a slab, the spherical symmetry is broken and the system is no longer homogeneous.69 Also, the number density becomes smaller due to the arbitrarily large slab length, resulting in a problem that the normalized g(r) at large r does not become unity. We are aware of some structure-based CG approaches to overcome this problem by introducing a slab RDF from the bulk RDF and the density profile.47 However, this structure-based CG approach mainly focuses on constructing CG potentials from analytical potentials that are parameterized using inverse Boltzmann iterations (IBI), and therefore does not provide a clear correspondence between the slab and bulk liquid. Since the inner region exhibits a uniform density along the slab axis, the structural characteristics of the inner region are expected to be similar to that of bulk liquid. To confirm the structural similarity between the inner region and the bulk liquid, we evaluated the g(r) of the inner region by applying periodic boundary conditions. We chose the following six molecules that showed different density profiles with MS-CG compared to the all-atom simulations: methanol, ethanol, acetone, acetonitrile, chloroform, and benzene. Inner regions for each interface were simply identified from the all-atom density profile depicted in Figure 3. The results are shown in Figure 5 with naïvely calculated g(r) profiles for the slab and bulk liquid systems.

16


Page 16 of 47

Page 17 of 47

!"#$%&'( !'

6728/99!:385 ;0028!:385 !:385

!&#$

6728/99!:385 ;0028!:385 !:385

!&#$

!& !%#$

!"#$

!"#$ !(

!)

!*

!%"

!%&

!"#$ !" !&

!(

+,-./012!345

5617.88!9274 ://17!9274 ;

4-5506!7889-,:; 4 ?=>

!" !%

-./01234!567

!&#

!%#

!"#$

!"#$%&'%&''(#%#(!&)'

!"#$%&'%&''(#%#(!&)'

!'#&

!(#$

*+"#'&"

!(#'

!"#

!"#$%&'%&''(#%#(!&)'

!(

!"#$%&'%&''(#%#(!&)'

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


!$ !(#' !(#& !(#% !(#$ !( !"#' !"#& !"#% !"#$ !"

4-5506!7889-,:; 4 ?=>

!(#

!%

!&

!'

!(" !($ !(% !(& !('

)*+,-./0!123

Figure 8: Radial distribution function g(r) calculated in the inner region (see Figure 4 for detailed definition) from the all-atom (red line), MS-CG (blue line), and UCG simulations (green line) for liquid/vapor interfaces composed of (a) methanol, (b) ethanol, (c) acetone, (d) acetonitrile, (e) chloroform, and (f) benzene molecules. In Figure 8(a), the RDF of methanol obtained from the UCG simulation captures the peak at the first maximum while the larger peak predicted from the MS-CG model represents an overly structured liquid. In order to calculate the inner RDF of the MS-CG simulation, which has a non-uniform density at the liquid-like part of the slab structure, we used the same interval defined in the UCG model. Figure 8(c) especially supports the notion that the UCG model

25



captures the inner RDF from the all-atom reference rather than the MS-CG model for acetone molecule; the acetone MS-CG model does not locate the first maximum and inflection point at the right distances. However, the g(r) from the UCG calculation can reproduce these peaks at exact distances. Similarly, the MS-CG model of acetonitrile represents the ordered structures in the long-range regime as seen from peaks at 3.5 and 7.0 Å in Figure 8(d). In contrast, the UCG model is capable of reproducing major peaks at 3.5, 5.0, 7.0 and 9.0 Å, which validates that the UCG model and reference model have almost identical structures within the liquid phase. This almost exact agreement suggests that this density-based UCG approach can even enhance the quality of CG simulation from a strongly-maintained linear molecule. The reason the UCG acetonitrile model performs better than the UCG ethanol model is attributed to the much smaller size of acetonitrile than of ethanol (-CN < -CH2OH), even though acetonitrile is more linear than ethanol due to the rigidity of the C ! N bond. Close agreement between the all-atom and UCG inner g(r) functions might imply that the interactions from the inner region are transferable to bulk liquids from Henderson’s theorem75 if one could extract the interactions in the inner region from the interface system properly. However, Henderson’s theorem is not directly applicable to slab systems because it only holds in homogeneous systems. Also, the actual structural correlation between CG pairs at the inner region might slightly deviate from the inner RDF, since the inner RDFs are obtained from periodic images of solely the inner regions, thus neglecting the contributions from CG pairs between the inner and outer regions. In the current UCG model, we envisage that a locally more dense state would correspond to a CG particle in the inner region. To check the transferability of these CG interactions, we calculated the effective CG pair force !"U (R) obtained from force matching. We chose to calculate the effective pair force because in the liquid/vapor system the interaction potential U (R ) at a large distance can depend on the interfacial structure. Since the inner and outer regions are defined in the interface system by excluding the vapor phase, CG beads in the vapor phase might affect the long-range interactions in the inner and the outer regions. Even minor differences between long-range forces

26


Page 26 of 47

Page 27 of 47

will accumulate (through integration) into non-negligible differences in short-range interaction potentials. Therefore, we chose the pair force !"U (R) rather than the pair potential U (R ) to eliminate the effect of interfacial interactions of the long-range on short-range interactions. The resulting effective pair interaction forces between different interfacial systems are illustrated in Figure 9. A slight difference of pair forces between the MS-CG and UCG is seen at long distances (> 8 Å), but the difference is negligible, smaller than 0.02 kcal/mol/Å in the methanol case, and show the same asymptotic behavior., i.e. the long-range forces in the interface and bulk cubic systems are in reasonable agreement.

!"#$%&'(

)#$%&'(

%" %#

&(

%$!"#$%&!'#(&')$*+,-"+,.(

%$ %& !$ !#

#( ;/1!>/1 ;/1 ;/1 ;/1 ;?

!"#$ +,-#$%$#&'("%&)$*

!"# !"#$ !""#$%$#&'("%&)$*

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 47

!&#( !&#' !"#% !"#$

!'

!(

!$

!)

!*

-./01234!567

!+

!,

899:10;< =>:?@ A?@

!'#$ !' !&#( !&#' !"#% !"#$ !"

!" !&

!%#

!'#%

!'

!)

!$

!*

!(

!+

!%

!,

-./01234!567

!'

!)

!$

!*

!(

!+

!%

!,

-./01234!567

Figure 14: Radial distribution function g(r) in the liquid/liquid interface from the all-atom (a red line), MS-CG (a blue line), and UCG simulations (a green line) for (a) methanol-methanol pairs in the inner region, (b) methanol-carbon tetrachloride pairs in the outer region, and (c) carbon tetrachloride-carbon tetrachloride pairs in the outer region.

36


Page 37 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


In Figure 14, methanol-methanol, methanol-carbon tetrachloride, and carbon tetrachloridecarbon tetrachloride pairs are chosen to provide the structural characteristics of each pair in its most populous region. For example, carbon tetrachloride is sparsely distributed in the inner region, thus carbon tetrachloride-carbon tetrachloride pairs in the inner region are not well sampled. Improvements can be seen by the agreement between UCG and all-atom RDFs. For methanol-methanol pairs in the inner region, the MS-CG model gives an over-structured liquid with higher g(r) values at the first peak, while the UCG model corrects this intensity as well as the correct long-range structure at the second minimum of the RDF. Similarly, the g(r) of methanol-carbon tetrachloride pairs in the outer region largely deviates in the MS-CG model, but the UCG model is capable of successfully capturing the structure at all-atom resolution. The final carbon tetrachloride-carbon tetrachloride pairs in the outer region confirm that the UCG model gives the correct local structures in the system, despite its strong phase-separation. As discussed earlier, we also constructed the UCG model by using a homogeneous density as an order parameter. However, the resultant density profiles from homogeneous UCG model give uniform liquids with invariant density along the slab axis. From the detailed discussion given in Supporting Information (see Section S3), we found out that the homogeneous local density is not suitable for imposing internal states for this system. Furthermore, we checked the transferability of the UCG interactions by comparing them to interactions from bulk liquid states. Distinct from the liquid/vapor interface where the system is comprised of a single homogeneous molecule, the liquid/liquid interface contains two different molecules, which are distributed throughout both inner and outer regions. Thus, the bulk system that corresponds to the liquid/liquid interface is not simply a pure homogeneous liquid, but a mixture of two liquid molecules. Since the inner and outer regions have two molecules with different density, we constructed two bulk systems (Bulk 1 and Bulk 2) according to the ratio of methanol number density to carbon tetrachloride number density obtained from the all-atom simulation. In other words, Bulk 1 was designed to correspond to the outer region, and Bulk 2 was set to resemble the chemical environments in the inner region. In practice, we constructed

37



the Bulk 1 system using 110 methanol molecules and 60 carbon tetrachloride molecules with a number ratio of 1.83 in a cubic box under periodic boundary condition. Similarly, the Bulk 2 system represents a number ratio of 13.0 that was composed of 520 methanol molecules and 40 carbon tetrachloride molecules in a cubic box. From the initial structure prepared by the Packmol54 package, we followed the procedure applied to the bulk system for the liquid/vapor interface systems. We minimized the energy of the initial configurations, then annealed the system to T = 300 K using constant NVT dynamics with a Nosé-Hoover thermostat61, 62. After a constant NPT equilibration with a Nosé-Hoover thermostat and Andersen barostat63, the final NVT MD simulation was performed for 5 ns. To verify the final bulk structures, we calculated the local cross-density histogram from the Bulk 1 and Bulk 2 systems and compared them to that of the outer and inner regions, respectively. Normalized density histograms described in Supporting Information Figure S7 and S8 imply that the Bulk 1 and Bulk 2 systems are able to represent the cross-density from the outer and inner regions of the interface. To examine the transferability of UCG interactions to describe bulk phase interactions, we first calculated the effective CG pair force !"U (R) as discussed in the previous liquid/vapor interface section. Even though the effective pair force between the methanol CG pairs show transferability with a smooth force profile (see Figure S9), the pair force between carbon tetrachloride molecules obtained through force matching show large fluctuations. In particular, this fluctuation was mainly observed in the far-far interaction that corresponds to the inner region. This statistical noise can be understood from the average number of carbon tetrachloride in the inner region during the simulation: only 18 molecules were present in the inner region out of 168 molecules. Thus, in order to reduce the statistical noise, we calculated the effective CG potentials. In the case of the liquid/liquid interface, the interaction potential U (R ) is not expected to be significantly affected by the interfacial structure. This is mainly due to the structural characteristics of the methanol/carbon tetrachloride interface. First, we divided the whole system

38


Page 38 of 47

Page 39 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60




between methanol far from carbon tetrachloride is the same as the methanol pair interaction in Bulk 2. As expected, an opposite correspondence between the UCG statewise interactions and bulk liquids is obtained. Due to the different chemical environments described in Figure 11, the MS-CG potentials from the Bulk 1 liquid is in a close agreement with carbon tetrachloride far from methanol, while those from Bulk 2 are close to the potential between more dense UCG states.

Table 4: Euclidian distance between the pair interactions (RMSD) from the MS-CG and UCG models for methanol-methanol and chloroform-chloroform pair interactions. Types of UCG models: “Near” is a more dense state and “Far” is a less dense state. In order to avoid numerical instability at the hard-core region of the PMF, reported distances are calculated in a following range: U (r) < 5 kcal/mol. (a) Methanol-Methanol Inter-PMF Distance (kcal/mol) MS-CG Model Type Bulk 1 Bulk 2 Interface 0.2660 1.2144 1.0514 Near-Near UCG 1.5074 0.1086 0.3956 Far-Far (b) Chloroform-Chloroform Inter-PMF Distance (kcal/mol) MS-CG Model Type Bulk 1 Bulk 2 Interface Near-Near 1.0299 0.1892 0.8023 UCG Far-Far 0.4074 0.9387 0.0367 Again, we calculated the average Euclidian distance between the PMFs, listed in Table 4, to examine the transferability in a quantitative manner. Comparable to Figure 15, Table 4(a) confirms an agreement between the locally more dense UCG state and Bulk 1 interactions (closer than Bulk 2 by 4.5 times), and between the locally less dense UCG state and Bulk 2 interactions (closer than Bulk 1 by 15 times). In the case of the chloroform-chloroform interaction, the plain MS-CG interaction resembles the UCG interactions between locally less dense states rather than that of the more dense states by 22 times, suggesting that a naïve treatment of the liquid/liquid interface fails to capture a locally more dense environment. Therefore, from Figure 9 and 15 in conjunction with Table 3 and 4, we have demonstrated how

40


Page 40 of 47

Page 41 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the UCG methodology can be applied to complex interfacial systems and result in transferable CG interactions.

4! Conclusions In this paper, we have developed ultra-coarse-grained (UCG) models for liquid/vapor and liquid/liquid interfaces by implementing the recently proposed UCG-RLE theory.44 This work provides a more rigorous procedure for applying the UCG methodology in realistic inhomogeneous systems instead of the idealized examples in the third paper of the UCG series. From the density profile along the slab axis, we identified that the inner and outer regions experience different chemical environments that can be classified by local densities. To modulate the interactions, the substate probabilities were designed to be dependent on the local density with a density threshold that distinguishes locally more dense and locally less dense states. The UCG simulations were efficiently implemented along with pre-existing MS-CG theory via the LAMMPS simulation engine. We constructed both MS-CG and UCG models of various molecules with different intermolecular interactions and geometries. In turn, the present work improves the quality of CG simulations at interfaces compared to the conventional MS-CG approach. Though the MS-CG approach provides a systematic coarse-graining methodology for reproducing important properties in many systems, we have shown that single-bead MS-CG models fail to capture the proper interfacial interactions when the target molecule has aspherical symmetry. Even in the case of near-spherical molecules, the MS-CG model predicts a collapse in the liquid/liquid interface where spherical carbon tetrachloride interacts with methanol. In contrast, we demonstrate that the resulting UCG models dramatically reproduce the density profiles from the all-atom reference while preserving local structures. For the liquid/vapor interfaces, the UCG interactions between CG particles in the more dense state were almost identical to the interactions between pure CG liquids in the bulk state, indicating transferability of the current UCG force fields. For an even more complicated liquid/liquid interface where the MS-CG method exhibits clear limitations, the interactions obtained from the UCG model successfully recapitulate the interactions from

41



the bulk systems that correspond to each region. Surprisingly, transferability was identified in all systems; this is particularly important given that designing transferable CG interactions via bottom-up CG methods is still challenging. The current UCG-RLE approach suggests that one can design transferable CG interactions in different systems using proper order parameters that can be derived from the system of interest using the UCG framework. A future direction will be to study which variables should be used to construct internal states that are rapidly mixed in many different systems. On the other hand, in terms of more physical and practical applications, we believe that our UCG approach can treat more complex systems such as amphiphilic molecules at a liquid/liquid interface or even more collective systems with hydrophobic interactions. In spite of the complicated nature of amphiphilic and hydrophobic molecules in a given system, we expect that the local density modulates the interaction and this could be easily modeled by the UCG-RLE framework, including its extendibility and transferability.

Associated Content Supporting Information. The SI contains the speed-up analysis for the liquid/vapor interface and the implementation of homogeneous density UCG model to the liquid/liquid interface by comparing with the crossdensity UCG model in terms of M-statistics. Validation of the constructed bulk systems for liquid/liquid interfaces is also discussed. This material is available free of charge via the Internet at http://pubs.acs.org.

Acknowledgements This material is based upon work supported by the National Science Foundation (NSF Grant CHE-1465248). Simulations were performed using the resources provided by the University of Chicago Research Computing Center (RCC). J.J. acknowledges a graduate fellowship from the

42


Page 42 of 47

Page 43 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Kwanjeong Educational Foundation. We thank Dr. Alex Pak for his comments and a careful reading of the manuscript.

References 1. Robinson, G.; Singh, S.; Zhu, S.; Evans, M., Experimental overviews and computational methodologies. World Scientific: 1996. 2. Cox, M.; Hidalgo, M.; Valiente, M., Solvent extraction for the 21st century. SCI: 2001. 3. Kazarinov, V. E., The interface structure and electrochemical processes at the boundary between two immiscible liquids. Springer: Berlin: 1987. 4. Starks, C. M.; Liotta, C. L.; Halpern, M. E., Phase-transfer catalysis: Fundamentals I. In Phase-Transfer Catalysis, Springer: 1994; pp 23-47. 5. Volkov, A. G., Interfacial catalysis. CRC Press: 2002. 6. Welton, T., Room-temperature ionic liquids. Solvents for synthesis and catalysis. Chem. Rev. 1999, 99, 2071-2084. 7. Gennis, R. B., Biomembranes: molecular structure and function. Springer Science & Business Media: 2013. 8. Haymet, A.; Oxtoby, D. W., A molecular theory for the solid–liquid interface. J. Chem. Phys. 1981, 74, 2559-2565. 9. Wilson, M. A.; Pohorille, A.; Pratt, L. R., Molecular dynamics of the water liquid-vapor interface. J. Phys. Chem. 1987, 91, 4873-4878. 10. Dang, L. X.; Chang, T.-M., Molecular dynamics study of water clusters, liquid, and liquid–vapor interface of water with many-body potentials. J. Chem. Phys. 1997, 106, 81498159. 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. Sides, S. W.; Grest, G. S.; Lacasse, M.-D., Capillary waves at liquid-vapor interfaces: A molecular dynamics simulation. Phys. Rev. E 1999, 60, 6708. 13. Kuo, I.-F. W.; Mundy, C. J., An ab initio molecular dynamics study of the aqueous liquid-vapor interface. Science 2004, 303, 658-660. 14. Scheraga, H. A.; Khalili, M.; Liwo, A., Protein-folding dynamics: overview of molecular simulation techniques. Annu. Rev. Phys. Chem. 2007, 58, 57-83. 15. Liwo, A.; Czaplewski, C.; O"dziej, S.; Scheraga, H. A., Computational techniques for efficient conformational sampling of proteins. Curr. Opin. Struct. Biol. 2008, 18, 134-139. 16. Voth, G. A., Coarse-graining of condensed phase and biomolecular systems. CRC press: 2008. 17. Peter, C.; Kremer, K., Multiscale simulation of soft matter systems–from the atomistic to the coarse-grained level and back. Soft Matter 2009, 5, 4357-4366. 18. Murtola, T.; Bunker, A.; Vattulainen, I.; Deserno, M.; Karttunen, M., Multiscale modeling of emergent materials: biological and soft matter. Phys. Chem. Chem. Phys. 2009, 11, 1869-1892. 19. Noid, W., Perspective: Coarse-grained models for biomolecular systems. J. Chem. Phys. 2013, 139, 090901.

43



20. Saunders, M. G.; Voth, G. A., Coarse-graining methods for computational biology. Annu. Rev. Biophys. 2013, 42, 73-93. 21. Marrink, S. J.; Tieleman, D. P., Perspective on the Martini model. Chem. Soc. Rev. 2013, 42, 6801-6822. 22. Barker, J.; Fisher, R.; Watts, R., Liquid argon: Monte carlo and molecular dynamics calculations. Mol. Phys. 1971, 21, 657-673. 23. Duncan, J.; Kelly, R.; Nivellini, G.; Tullini, F., The emprical general harmonic force field of ethane. J. Mol. Spectrosc. 1983, 98, 87-110. 24. Lindorff-Larsen, K.; Maragakis, P.; Piana, S.; Eastwood, M. P.; Dror, R. O.; Shaw, D. E., Systematic validation of protein force fields against experimental data. PLoS ONE 2012, 7, e32131. 25. Wagner, J. W.; Dama, J. F.; Durumeric, A. E.; Voth, G. A., On the representability problem and the physical meaning of coarse-grained models. J. Chem. Phys. 2016, 145, 044108. 26. Reith, D.; Pütz, M.; Müller! Plathe, F., Deriving effective mesoscale potentials from atomistic simulations. J. Comput. Chem. 2003, 24, 1624-1636. 27. Lyubartsev, A. P.; Laaksonen, A., Calculation of effective interaction potentials from radial distribution functions: A reverse Monte Carlo approach. Phys. Rev. E 1995, 52, 3730. 28. Izvekov, S.; Voth, G. A., A multiscale coarse-graining method for biomolecular systems. J. Phys. Chem. B 2005, 109, 2469-2473. 29. Izvekov, S.; Voth, G. A., Multiscale coarse graining of liquid-state systems. J. Chem. Phys. 2005, 123, 134105. 30. Noid, W.; Chu, J.-W.; Ayton, G. S.; Krishna, V.; Izvekov, S.; Voth, G. A.; Das, A.; Andersen, H. C., The multiscale coarse-graining method. I. A rigorous bridge between atomistic and coarse-grained models. J. Chem. Phys. 2008, 128, 244114. 31. Noid, W.; Liu, P.; Wang, Y.; Chu, J.-W.; Ayton, G. S.; Izvekov, S.; Andersen, H. C.; Voth, G. A., The multiscale coarse-graining method. II. Numerical implementation for coarsegrained molecular models. J. Chem. Phys. 2008, 128, 244115. 32. Shell, M. S., The relative entropy is fundamental to multiscale and inverse thermodynamic problems. J. Chem. Phys. 2008, 129, 144108. 33. Chaimovich, A.; Shell, M. S., Anomalous waterlike behavior in spherically-symmetric water models optimized with the relative entropy. Phys. Chem. Chem. Phys. 2009, 11, 19011915. 34. Chaimovich, A.; Shell, M. S., Relative entropy as a universal metric for multiscale errors. Phys. Rev. E 2010, 81, 060104. 35. Chaimovich, A.; Shell, M. S., Coarse-graining errors and numerical optimization using a relative entropy framework. J. Chem. Phys. 2011, 134, 094112. 36. Dunn, N. J.; Foley, T. T.; Noid, W. G., Van der Waals perspective on coarse-graining: Progress toward solving representability and transferability problems. Acc. Chem. Res. 2016, 49, 2832-2840. 37. Louis, A., Beware of density dependent pair potentials. J. Phys.: Condens. Matter 2002, 14, 9187. 38. Espanol, P.; Revenga, M., Smoothed dissipative particle dynamics. Phys. Rev. E 2003, 67, 026705. 39. Whitelam, S.; Rogers, C.; Pasqua, A.; Paavola, C.; Trent, J.; Geissler, P. L., The impact of conformational fluctuations on self-assembly: Cooperative aggregation of archaeal chaperonin proteins. Nano Lett. 2008, 9, 292-297.

44


Page 44 of 47

Page 45 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


40. Allen, E. C.; Rutledge, G. C., A novel algorithm for creating coarse-grained, density dependent implicit solvent models. J. Chem. Phys. 2008, 128, 154115. 41. Wagner, J. W.; Dannenhoffer-Lafage, T.; Jin, J.; Voth, G. A., Extending the range and physical accuracy of coarse-grained models: Order parameter dependent interactions. J. Chem. Phys. 2017, 147, 044113. 42. Dama, J. F.; Sinitskiy, A. V.; McCullagh, M.; Weare, J.; Roux, B.; Dinner, A. R.; Voth, G. A., The theory of Ultra-Coarse-Graining. 1. General principles. J. Chem. Theory Comput. 2013, 9, 2466-2480. 43. Davtyan, A.; Dama, J. F.; Sinitskiy, A. V.; Voth, G. A., The theory of Ultra-CoarseGraining. 2. Numerical implementation. J. Chem. Theory Comput. 2014, 10, 5265-5275. 44. Dama, J. F.; Jin, J.; Voth, G. A., The theory of Ultra-Coarse-Graining. 3. Coarse-grained sites with rapid local equilibrium of internal states. J. Chem. Theory Comput. 2017, 13, 10101022. 45. Willard, A. P.; Chandler, D., Instantaneous liquid interfaces. J. Phys. Chem. B 2010, 114, 1954-1958. 46. Ghoufi, A.; Malfreyt, P., Mesoscale modeling of the water liquid-vapor interface: A surface tension calculation. Phys. Rev. E 2011, 83, 051601. 47. Jochum, M.; Andrienko, D.; Kremer, K.; Peter, C., Structure-based coarse-graining in liquid slabs. J. Chem. Phys. 2012, 137, 064102. 48. Wang, Y.; Jiang, W.; Yan, T.; Voth, G. A., Understanding ionic liquids through atomistic and coarse-grained molecular dynamics simulations. Acc. Chem. Res. 2007, 40, 1193-1199. 49. Jiang, W.; Wang, Y.; Yan, T.; Voth, G. A., A multiscale coarse-graining study of the liquid/vacuum interface of room-temperature ionic liquids with alkyl substituents of different lengths. J. Phys. Chem. C 2008, 112, 1132-1139. 50. Wang, Y.; Feng, S.; Voth, G. A., Transferable coarse-grained models for ionic liquids. J. Chem. Theory Comput. 2009, 5, 1091-1098. 51. Merlet, C.; Salanne, M.; Rotenberg, B., New coarse-grained models of imidazolium ionic liquids for bulk and interfacial molecular simulations. J. Phys. Chem. C 2012, 116, 7687-7693. 52. Sega, M.; Fábián, B.; Horvai, G.; Jedlovszky, P., How is the surface tension of various liquids distributed along the interface normal? J. Phys. Chem. C 2016, 120, 27468-27477. 53. Humphrey, W.; Dalke, A.; Schulten, K., VMD: visual molecular dynamics. J. Molec. Graphics 1996, 14, 33-38. 54. 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. 55. Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J., Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids. J. Am. Chem. Soc. 1996, 118, 11225-11236. 56. Kaminski, G. A.; Friesner, R. A.; Tirado-Rives, J.; Jorgensen, W. L., Evaluation and reparametrization of the OPLS-AA force field for proteins via comparison with accurate quantum chemical calculations on peptides. J. Phys. Chem. B 2001, 105, 6474-6487. 57. Hockney, R. W.; Eastwood, J. W., Computer simulation using particles. CRC Press: 1988. 58. Plimpton, S., Fast parallel algorithms for short-range molecular dynamics. J. Comput. Phys. 1995, 117, 1-19.

45



59. Brown, W. M.; Wang, P.; Plimpton, S. J.; Tharrington, A. N., Implementing molecular dynamics on hybrid high performance computers–short range forces. Comput. Phys. Commun. 2011, 182, 898-911. 60. Brown, W. M.; Kohlmeyer, A.; Plimpton, S. J.; Tharrington, A. N., Implementing molecular dynamics on hybrid high performance computers–Particle–particle particle-mesh. Comput. Phys. Commun. 2012, 183, 449-459. 61. Nosé, S., A unified formulation of the constant temperature molecular dynamics methods. J. Chem. Phys. 1984, 81, 511-519. 62. Hoover, W. G., Canonical dynamics: equilibrium phase-space distributions. Phys. Rev. A 1985, 31, 1695. 63. Andersen, H. C., Molecular dynamics simulations at constant pressure and/or temperature. J. Chem. Phys. 1980, 72, 2384-2393. 64. Zhang, Y.; Feller, S. E.; Brooks, B. R.; Pastor, R. W., Computer simulation of liquid/liquid interfaces. I. Theory and application to octane/water. J. Chem. Phys. 1995, 103, 10252-10266. 65. Van Buuren, A. R.; Marrink, S.-J.; Berendsen, H. J., A molecular dynamics study of the decane/water interface. J. Phys. Chem. 1993, 97, 9206-9206. 66. Carpenter, I. L.; Hehre, W. J., A molecular dynamics study of the hexane/water interface. J. Phys. Chem. 1990, 94. 67. Linse, P., Monte Carlo simulation of liquid–liquid benzene–water interface. J. Chem. Phys. 1987, 86, 4177-4187. 68. Lu, L.; Izvekov, S.; Das, A.; Andersen, H. C.; Voth, G. A., Efficient, regularized, and scalable algorithms for multiscale coarse-graining. J. Chem. Theory Comput. 2010, 6, 954-965. 69. Benjamin, I., Chemical reactions and solvation at liquid interfaces: A microscopic perspective. Chem. Rev. 1996, 96, 1449-1476. 70. Matsumoto, M.; Kataoka, Y., Molecular orientation near liquid–vapor interface of methanol: Simulational study. J. Chem. Phys. 1989, 90, 2398-2407. 71. Cipcigan, F. S.; Sokhan, V. P.; Jones, A. P.; Crain, J.; Martyna, G. J., Hydrogen bonding and molecular orientation at the liquid–vapour interface of water. Phys. Chem. Chem. Phys. 2015, 17, 8660-8669. 72. Yeh, Y. L.; Zhang, C.; Held, H.; Mebel, A.; Wei, X.; Lin, S.; Shen, Y., Structure of the acetone liquid/vapor interface. J. Chem. Phys. 2001, 114, 1837-1843. 73. Das, A.; Lu, L.; Andersen, H. C.; Voth, G. A., The multiscale coarse-graining method. X. Improved algorithms for constructing coarse-grained potentials for molecular systems. J. Chem. Phys. 2012, 136, 194115. 74. Chapela, G. A.; Saville, G.; Thompson, S. M.; Rowlinson, J. S., Computer simulation of a gas–liquid surface. Part 1. J. Chem. Soc. Faraday Trans. 2 1977, 73, 1133-1144. 75. Henderson, R., A uniqueness theorem for fluid pair correlation functions. Phys. Lett. A 1974, 49, 197-198. 76. Fletcher, A. N., Effect of carbon tetrachloride upon the self-association of 1-octanol. J. Phys. Chem. 1969, 73, 2217-2225. 77. Durov, V. A.; Shilov, I. Y., Supramolecular structure and physicochemical properties of the mixture tetrachloromethane-methanol. J. Mol. Liq. 2001, 92, 165-184. 78. Torii, H., Atomic quadrupolar effect in the methanol–CCl 4 and water–CCl 4 intermolecular interactions. Chem. Phys. Lett. 2004, 393, 153-158.

46


Page 46 of 47

Page 47 of 47 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



Ultra-Coarse-Grained Models Allow for an Accurate and Transferable

Recommend Documents