Accuracy, Transferability, and Efficiency of Coarse-Grained Models of

Subscriber access provided by University of South Dakota

Feature Article

On the Accuracy, Transferability, and Efficiency of Coarse-Grained Models of Molecular Liquids Marina Giuseppina Guenza, Mohammadhasan Dinpajooh, James McCarty, and Ivan Lyubimov J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.8b06687 • Publication Date (Web): 28 Aug 2018 Downloaded from http://pubs.acs.org on September 1, 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.

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 68 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 Journal of Physical Chemistry

On the Accuracy, Transferability, and Efficiency of Coarse-Grained Models of Molecular Liquids M. G. Guenza,ú,† M. Dinpajooh,† J. McCarty,,‡ and I. Y. Lyubimov,¶ †Department of Chemistry and Biochemistry, and Institute of Theoretical Science, University of Oregon, Eugene, OR 97403, USA

‡Present Address: Department of Chemistry and Biochemistry, University of California Santa Barbara, Santa Barbara, CA 93106, USA E-mail: [email protected] Abstract Coarse-graining (CG) approaches are becoming essential tools in the study of complex systems because they can considerably speed up computer simulations, with the promise of determining properties in a range of lengthscales and timescales never before possible. While much progress in this field has been achieved in recent years, application of CG methods is still inhibited by the limited understanding of a number of conceptual points that need to be resolved to open up the field of CG to a wide range of applications in material science and biology. In this paper we present some of the key findings that emerged from the development of the Integral Equation theory of Coarse-Graining (IECG), which addresses some of the fundamental questions in coarse-graining. Although the IECG method pertains to the CG of polymer liquids, and specifically homopolymer melts illustrated here, many of the results that emerge from the study of the IECG approach are general and apply to the CG of any molecular liquid. Through this method, we developed a formal relation between the statistical mechanics of CG and a number of predicted physical properties. Based on the theory

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 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

of liquids, the IECG affords the analytical solution of the intermolecular potential for macromolecules represented by a Markov chain of CG sites, thus providing a transparent tool for analysis of the properties in coarse-graining. We identify three key requirements that render a CG model useful: accuracy, transferability, and computational efficiency. When these three requirements are fulfilled, the CG model becomes widely applicable and useful for studying regions in the phase space that are not covered by atomistic simulations. In the process, the IECG answers formally a number of relevant questions on how structural, thermodynamic, and dynamical properties are modified during coarse-graining. It sheds light upon how the level of CG affects the shape of the CG potential, and how, in turn, the shape of the potential affects the physical properties. It tests the validity of selecting the potential-of-mean force as the effective pairwise CG potential, and the role of higher-order many-body corrections to the pairwise potential to recover structural and thermodynamic consistency of the CG model. Because the IECG theory can be analytically formalized, it does not suffer from the problem of transferability, and in the canonical ensemble, leads to consistent pair distribution functions, pressure, isothermal compressibility, and excess free energy at variable levels of CG from the atomistic to the ultra-CG model, where macromolecules are represented as inter-penetrable soft spheres.

1

Introduction

Computer simulations are well established as a tool to study the microscopic effects of molecular structure, such as local order and interchain packing, on the macroscopic properties sampled experimentally. 1,2 However, while atomistic Molecular Dynamics (MD) simulations are parameterized to reproduce a number of experimental properties, they are still computationally too costly to address many properties in the range of time and lengthscales of experimental interest. Thus, despite the rapid advances in computer technology, atomisticscale MD simulations are still limited in their range of applicability, so that it is often

2


Page 2 of 68

Page 3 of 68 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


impossible to directly connect experimental results to the microscopic structures through atomistic simulations. 3 A powerful strategy to overcome these limitations is to adopt a reduced, CG, representation of the molecules to be simulated, where local degrees of freedom are averaged out. 4–13 Because of the simplified representation, which reduces the level of detail and the number of degrees of freedom, the computational time of a CG simulation is decreased, opening up the opportunity of studying a system at a much longer timescale and larger lengthscale than the corresponding atomistic simulation. Once the CG simulation is combined with a local atomistic simulation one can obtain the full-scale information that is needed. 14–20 While this avenue has been pursued at length in recent years, a number of less understood points still remain, which have hindered the practical application of CG methods in computer simulations. In this paper we address some of these questions such as the accuracy, transferability, and efficiency of CG models with the help of the Integral Equation theory of Coarse-Graining or IECG method. 7,21–23 The advantage of this approach is that it is analytically solved for liquids of long polymer chains providing a transparent formalism, useful to evaluate the effects of coarse-graining on the structural, thermodynamic, and dynamical properties of molecular liquids in general. In the body of this paper, properties that are specific of polymeric liquids will be mentioned accordingly. The consequences of coarse-graining that are unveiled by the IECG analysis emerge, instead, from the statistical mechanics of coarse-graining and as such are general and pertain to any CG formalism: those will be referred as pertaining to molecular liquids. For example, it is often the case that CG models can reproduce with accuracy atomistic structural quantities, such as the pair distribution function, g(r), but cannot reproduce with the same accuracy thermodynamic quantities, such as the pressure or the excess free energy, which, in most cases, are not even reported. This observation seems to be in contradiction with what is known about the statistical mechanics of liquids, namely that the pair distribution function is unique, and that any thermodynamic quantity can be calculated from 3



the knowledge of g(r). 24,25 As we describe in more detail in the body of this paper, these observed inconsistencies can have multiple origins. When a system is coarse-grained, the interaction potential becomes a free energy and state dependent. Unless a well-defined statistical mechanics formalism is in place, the derivation of the effective CG potential becomes a complex numerical optimization problem. 26–30 In general, CG models are either bottom-up approaches, which aim at reproducing quantities from atomistic simulations, or top-down approaches, optimized to agree with some large-scale experimental quantity. In both cases, the numerical optimization of the potential depend on the appropriate selection of the initial conditions, and the choice of the physical quantities used to optimize the potential, sometimes leading to different solutions depending on the physical quantity that is optimized or the path selected during the optimization scheme. 31 While this problem is complex, the statistical mechanics of CG liquids can give helpful insights on its solution. Thus another point, confirmed by the IECG approach, is the fact that the observed inconsistency of the pressure in structurally-consistent CG models can be explained by considering that the pair distribution function in dense liquids is sensitive to the short-range, mostly repulsive, part of the potential, while the pressure is sensitive, through the virial equation, to the tail of the potential. 32,33 When numerical CG methods focus on reproducing g(r), or the related force, only the short-ranged part of the potential is optimized, and consistency in the pressure is not guaranteed. 21,34 Furthermore, given that the pair distribution function is largely insensitive to the precise shape of the potential, while different potentials can give pair distribution functions that are similar, 35 it is not surprising that a potential that appears to be correctly optimized to reproduce g(r) still finds some difficulty reproducing other physical quantities like the equation of state. Other observed inconsistencies can be explained in the framework of the analytical IECG formalism. One question that has emerged in the field of coarse-graining is if the potential of mean force (PMF) is a viable choice as the trial function for the optimized effective CG potential in a molecular liquid. The answer to this question is directly related to the issue 4


Page 4 of 68

Page 5 of 68 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


of how to include many-body effects in the training CG potential, and to the question if a pairwise CG potential can be, in any circumstance, a good representation of the real CG potential. In the IECG formalism, which is based on liquid theory, many-body contributions of the CG potential are easily identified in the closure equation and can be investigated. This IECG analysis shows that it is possible to derive pairwise CG potentials that correctly account for the many-body effects, thus reproducing with accuracy both structural and thermodynamic properties. In the IECG formalism, however, not all the thermodynamic properties are conserved. 7,21–23,36 The question of which properties are conserved and which properties are not conserved in CG models is important when numerical optimization procedures are used to optimize the CG potential. 37,38 In general it is difficult to sort out which property should be conserved, and which property should not be conserved (if any) by a given CG model, unless a formal analytical approach is available to distinguish these cases. In the IECG approach, physical quantities that depend on the number of degrees of freedom are modified as a consequence of coarsegraining, and the IECG shows that quantities such as the entropy and the internal energy are not conserved during CG, while the excess free energy is. McCarty et al. have shown that pair distribution functions, pressure, excess free energy, and isothermal compressibility are conserved by the IECG potential as the granularity of the CG model is varied from the atomistic resolution to the soft sphere representation. 22 They also observed that other thermodynamical quantities, like internal energy and entropy, are not conserved. In short, in the IECG approach, structural and thermodynamical consistencies are observed for the properties for which the change in the number of degrees of freedom results in the compensation of intramolecular and intermolecular contributions. As an example, Figure 1 shows the consistency of pressure and radial distribution functions for polymer melts with a number of monomers N = 300, from MD simulations that use the IECG potential. The comparison is made with the atomistic MD simulations of the same system. Also included in the figure are the lines representing the analytical solution of the IECG equation. As can be seen, the 5



differences between the atomistic and CG properties, such as pair correlation functions and pressure, are within the atomistic statistical uncertainties. 23,36 Furthermore, many studies have shown that diffusivity is also not conserved during coarse-graining. 39–42 In a series of papers Lyubimov et al. derived from the IECG approach, in the framework of Mori-Zwanzig projection operators, 43 the needed analytical corrections that reconstruct the atomistic diffusion from the mean-square displacements measured in CG simulations. Those corrections include both the rescaling of the friction coefficient and an entropy correction that derives from the smoothing of the free energy landscape due to coarse-graining. 39,40 With these corrections implemented, the data of diffusivity that are collected from the IECG simulations accurately reproduce the atomistic values of diffusivity. This example demonstrates the usefulness of an analytically formalized CG approach, like the IECG, which allows one to calculate the proper corrections that need to be implemented to recover the atomistic properties that are not conserved during coarse-graining. Finally, from the study of which properties are conserved and which properties are not conserved, the IECG suggests that when developing a CG model, care needs to be taken in selecting the physical quantities used in the numerical optimization of the CG potential among the properties that are conserved during coarse-graining. Quantities that are not conserved during coarse-graining should not be included in the optimization procedure as they could interfere with the proper convergence of the optimization method. As a consequence of coarse-graining, computer simulations speed up, giving the opportunity of studying a system at a much longer timescale and larger lengthscale than the corresponding atomistic simulation. Because of the simplified representation, which reduces the level of detail and the number of degrees of freedom, the computational time of a CG simulation is decreased. In general, the larger the degree of CG, the more substantial the computational gain. Considering the consequences of coarse-graining on the computational speed-up, one notices that they are complex and sometimes have opposite effects. 44 For example, Clark et al. have shown that as the level of CG increases, the range of the poten6


Page 6 of 68

Page 7 of 68

0.06

0.95 0.9

0.05

sbb

g (r)

0.04

Probability(P)

1

0.05

Probability(P)

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


342

344

346

0.8

atomistic nb = 4

0.75

348

P, atm

0.03

0.85

theory

0.7 0.65

atomistic simulations nb=1

0.02

0.6

nb=4

0

15

nb=6

45

60

nb=10

0.01

0

30 r, Å

-200

0

200

P, atm

400

600

800

Figure 1: Structural and thermodynamical consistencies for the IECG model, illustrated for a polymer melt with N = 300 at 503 K and a monomer density of 0.03296Å≠3 . Pressure distribution for the atomistic simulation (black curve). The blue line shows the average pressure for atomistic simulations with the simulation error bars obtained from block averages (see left panel). The red, orange, cyan, and magenta lines show the average pressure for coarse-grained (CG) simulations with 1, 4, 6, and 10 blobs per chain (nb ), respectively. Theoretical and simulated radial distribution functions (right panel) when the polymer is represented by four blobs. Atomistic simulations (blue triangle) are compared with CG simulations (orange square) and with the analytical solution (black line). The theoretical data and the data from the CG simulations are both within the error of the atomistic simulation. tial increases, 21,34 with the consequence that the IECG simulation requires larger simulation boxes, which contain a higher number of molecules than the atomistic simulations. 44 To improve computational efficiency one would be tempted to truncate the CG potential at short range. However, Dinpajooh at al. have shown that the truncation of the CG potential at relatively small distances has effects on both the structure and the thermodynamics of the liquid. 44 Additionally, they showed that at distances close to the first extremum of the energy, the structure is consistent with its atomistic counterpart, but the thermodynamics can be completely off unless the CG tail corrections are applied. 21,34,36,44 When structural and/or thermodynamic properties of the system are incorrectly predicted by a CG model, for example because of an incorrect truncation of the potential or, more generally, because of the difficulty in deriving the correct CG potential, having a large computational gain in the CG model is of limited advantage: a CG simulation that predicts properties in disagreement 7



with the corresponding atomistic description, and thus in disagreement with the experiments used to parameterize the atomistic potential, is useful only qualitatively. For example, if the CG model predicts incorrectly the excess free energy, it would be impossible for that CG model to identify the physical phase that is energetically stable at the given state point. This has important consequences in the studies of molecular liquids, where one would like to use simulations to predict the correct phase morphology at a given set of thermodynamic conditions. Summarizing, a CG model is most useful if it fulfills the following requirements: 1) the CG method is accurate in predicting physical quantities that are in agreement with the description at the atomistic level. The problem is more complex than this obvious statement indicates: the IECG approach shows that while some quantities are conserved during CG, other quantities that depend on the number of degrees of freedom are bound to change as a result of coarse-graining. A CG model that has a rigorous statistical-mechanical foundation allows one to evaluate formally the corrections that need to be applied to recover the atomistic properties from the ones that are measured in the CG-MD simulations. One example of this process for the IECG method is the reconstruction of the atomistic diffusion coefficient from the measured CG mean-square-displacement of the center-of-mass in the work of Lyubimov at al.’s on the reconstruction of the diffusion from soft-sphere IECG-MD simulations. 39,40 2) the CG method is transferable and predictive. Because CG models are useful when applied to study properties in conditions different from the ones used when deriving the potential, it is important for the CG potential to be transferable and applicable to other systems and/or in other state points. This requirement is not easy to fulfill because the CG potential is a free energy, and so it is parameter dependent. The CG method needs to be predictive in the sense that a CG-MD simulation has to be able to evaluate new physical properties that cannot be measured directly in the atomistic simulation, for example by allowing simulations of samples with higher molecular weight, or by extending the timescales 8


Page 8 of 68

Page 9 of 68 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


and/or the lengthscales sampled in the simulation. 3) the CG method is computationally efficient, with a net computational gain with respect to the atomistic simulations of a given system. CG simulations contain less information on the local scale than atomistic simulations. If CG models are not computationally advantageous with respect to the atomistic simulations, then there is little to gain in performing CG simulations. The IECG approach agrees with the three requirements listed above, with the advantage that it provides both a numerical and the analytical solution of the structural, thermodynamic and dynamical properties of the IECG description, as well as the solution of the pairwise IECG potential. In the IECG method the potential is derived from a generalized Ornstein-Zernike (OZ) equation for a liquid of macromolecules that are represented by atomistic and CG sites. The equation is based on the Polymer Reference Interaction Site Model (PRISM) by Schweizer and Curro, 45 which in turn is an extension of the Reference Interaction Site Model (RISM) integral equation theory to polymeric liquids. 46 The potential depends on one non trivial parameter, c0 , which is the direct correlation function at wave vector k = 0. The IECG potential is transferable to other state points, as well as applicable to other systems of the same type, because the dependence of the one parameter c0 is explicitly expressed by an equation of state and straightforward functional forms. The IECG method is also computationally efficient because coarse-graining smooths the free energies, thus reducing the number of degrees of freedom, and increases the diffusivity of the model. Because of these characteristics, the IECG approach is a useful method to study the effects of CG on the structural, thermodynamic, and dynamic properties of the system. Because of its structural and thermodynamic consistency, and because of its considerable computational speed-up, the IECG model may be a natural choice for Multiscale Modeling (MM) simulations where a system is represented, often in the same simulation box, by CG models at different resolutions. The advantage of MM simulations is to collect information on the atomistic scale in some spatial region of a system where atomistic information 9



is necessary, while taking advantage of the computational speed-up in the other parts of the simulation box, which are described at lower resolution, and where information at the atomistic scale is not needed. An example of this type of simulation is the Adaptive Resolution Simulation (AdResS) method, which has been proposed and developed by Delle Site, Paprotnik, Kremer and others. 47,48 The IECG method can also be used to perform large-scale, large-box, simulations of the CG system, which are combined a posteriori with small-scale, small-box, atomistic simulations of the local structure using Eq.2. 18,19 These types of simulations are particularly useful in the study of systems approaching a second-order phase transition, where fluctuations diverge in legthscale and timescale, and which are computationally prohibitive for atomistic resolution simulations. 49 Furthermore, in the soft-sphere representation, the IECG representation can be used to bridge into a continuum description, in analogy to recent work based on the Dissipative Particle Dynamics (DPD) model. 50,51 Contrasting from the DPD model, however, the IECG description has the advantage of conserving the pressure, and thus the related thermodynamic quantities, while DPD does not: this is a consequence of the fact that the DPD method adopts an unrealistic short-ranged interaction potential. In the mean time, IECG differs from the DPD model because while the DPD theory was constructed to reproduce the dynamical properties of the atomistic liquid, such as hydrodynamics and dynamical properties, the IECG displays a distinct speed up in the dynamics. 52 Diffusion is highly enhanced in the IECG formalism, while in DPD it is realistically reproduced. The IECG speed-up has the advantage of enhancing the computational efficiency of the CG simulation, while the realistic dynamics can be recovered after formal rescaling of the dynamical properties as measured in the IECG simulation. Obviously, the IECG potential can be used directly in other types of simulations, for example in the Widom method, where the softness of the potential facilitates the acceptance of configurations on the basis of their statistical weight. 1 Finally, it is worth mentioning that the IECG approach presented here applies to any 10


Page 10 of 68

Page 11 of 68 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


liquid of homopolymer systems, with side-chains and/or different monomeric structure, as far as the size of the CG units is taken to be larger than the polymer persistence length. 21,53 While the IECG theory is summarized here for liquids of compositionally homogeneous polymers, a version of this approach that models a polymer as an infinitesimally thin, infinitely long, chain (the so-called thread model), instead of the more realistic Markov chain treated in this paper, for block-copolymer liquids and for mixtures of polymers at different temperatures and compositions have been the subject of previous literature. 54,55 Therefore, in what follows the ”IECG-thread” term is used when the IECG approach is based on the thread model. The paper is structured as follows: after a brief overview of the IECG formalism in Section 2, Section 3 summarizes the analytical solution of the IECG potential, and provides some useful insights on the general properties that characterize CG potentials. Section 4 analyzes the relation between partition function, equation of state, and thermodynamic consistency, with a close look at the effect of optimization procedures in CG models that rely on finding an effective, pairwise, CG potential by numerical procedures. Given that the many-body nature of the CG potential is an important point of consideration when analyzing the accuracy of a CG model, Section 5 specifically focuses on the characteristics of the two-body potential of mean force, and on the effects that selecting this potential as the initial function in an optimization procedure, have on the accuracy of the CG model. In that Section, the paper also analyzes the consequence of treating the density as an active variable in the free-energy CG potential, and shows that in this approach, quite conveniently, the density dependence is included in the one parameter that is optimized against simulations, giving consistent results for the IECG potential with atomistic simulations. Because the IECG approach conserves the Equation Of State (EOS) but not the partition function, some thermodynamic and dynamical quantities are not conserved in CG simulations: however Section 6 shows how it is possible to use the IECG formalism to calculate the corrections that need to be implemented to recover the atomistic quantities from the ones predicted by 11



Page 12 of 68

the simulation that uses the IECG potential. A CG potential is of limited use if it is not transferable to other state points: in Section 7 this issue is discussed in the framework of the IECG theory with examples of structural and thermodynamic consistency for samples at varying densities and varying molecular sizes. Section 8 presents an analysis of the factors that need to be considered to fulfill the third requirement for a useful CG model, namely computational efficiency: how the choice in building a CG model affects not only structural and thermodynamic consistency, but also the shape of the CG potential, and how this in turn defines the speed up in the computational time are briefly analyzed in the context of the IECG approach. A summary with a brief discussion concludes the paper.

2

The Integral Equation theory of Coarse-Graining (IECG)

In the present and following sections we summarize the Integral Equation theory of CoarseGraining (IECG). This presentation is not intended to be an exhaustive description of the method and its applications, but rather a way to introduce the physical quantities that are used in the following sections. For a more complete description of the IECG method we refer to the papers previously published on this subject. 7,21–23,34,56 For a liquid of n molecules of chain length N in a volume V , the density of chains, flch = n/V , is related to the liquid monomer density by flm = flch N . Every chain in the liquid is partitioned in a variable number of CG units or blobs, nb , each containing a number Nb = N/nb of monomers, with the blob density flb = flm /Nb . When the number of blobs is one, nb = 1, the sample is represented as a liquid of soft spheres, where each polymer is a point particle centered in its center-of-mass and interacting with a bound potential. 34,53–58 Ô Ô The effective blob segment length is defined as ‡b = R/ N where Rb = R/ nb is the blob size and R is the polymer end-to-end distance. The IECG formalism for a CG liquid of polymers was derived starting from the generalized Ornstein-Zernike (OZ) equation, 25 as a formal approach to coarse-grain macromolecular

12


Page 13 of 68 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


liquids. In the OZ equation CG sites are considered as fictitious, while monomer sites are real sites. Coarse-grained or fictitious interacting sites are located at the center of mass of the polymer chain for the soft sphere model or at the centers of mass of several monomers along the same chain for the connected blob model (see Figure 2). 34

Figure 2: Illustration of polymer chains represented by four blobs (nb = 4), where each polymer chain consists of 300 monomers (N = 300) and the number of monomers in each blob is large enough to allow one to use the Gaussian statistics to obtain the intramolecular correlation (see Section 2). In the integral equation coarse-graining (IECG) method the blobblob correlation function is obtained in terms of monomer-monomer correlation functions (see Eq.1). The OZ equation is solved for a system that contains monomers and CG sites, 56,59 to give the intermolecular blob-blob (bb) pair distribution function of the polymer melts, which is related to the atomistic, monomer-monomer, pair distribution by ˆ bm (k)]2 ˆhbb (k) = [ ˆ mm (k) . h mm 2 ˆ [ (k)]

(1)

Here ˆ bm (k) and ˆ mm (k) are the normalized intramolecular blob-monomer (bm) and monomermonomer (mm) pair distribution functions. 23,56 Because the blobs represent the CG sites and the monomers correspond to the atomistic sites, these are the distributions of atomistic units around a CG site, and the relative distribution of atomistic sites in the molecule. 13



Page 14 of 68

The OZ equation thus reduces to 2

ˆ bb (k) = ≠ Nb nb h flm where

b

b

[ ˆ bm (k)]2 , 1 + nb b ˆ mm (k)

(2)

= ≠Nb flm c0 , and c0 is the monomer-monomer direct correlation function at k = 0,

which captures the large lengthscale behavior of monomeric units.

The OZ equation can then be used for a system consisting of the CG units (blobs) to calculate the blob-blob direct correlation function as cˆbb (k) =

ˆ bb (k) h i . ˆ bb (k) nb ˆ bb (k) nb ˆ bb (k) + flb h h

(3)

Spatial and wavevector representations of the structural pair correlation functions are connected through the Fourier transform.

2.1

Effective CG potentials and the potential of mean force

In liquid state theory, the correlation functions can be used to determine the properties of homogeneous classical liquids. The OZ equation (see Eq.3) relates the total, direct, and intramolecular correlation functions. However, in order to solve the OZ equation, one needs to apply an additional relation between the total and direct correlation functions, which is the so-called ‘closure equation’ bb bb cbb (r) = hbb (r) ≠ —Ueff (r) + —UPMF (r) + B(r) ,

(4)

with — the inverse temperature. Here B(r) is the bridge function, originated from the bb (r) is the effective connected structure of the clusters in the diagrammatic expansion, Ueff bb (r) is the PMF, CG potential, and UPMF

bb (r) = ≠kB T ln [hbb (r) + 1] . UPMF

14


(5)

Page 15 of 68 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


Following the traditional definition in statistical mechanics, the PMF is a function of the two-body radial distribution function. 25 The closure equation accounts for the many-body contributions to the pairwise potential, and is so-called because it projects the higher-than-second order contributions to the expansion in the density of the pair correlation function into the second order term, leading to the pairwise interaction potential. The choice of the most appropriate closure for the system under consideration depends to a large extent on the shape of the potential. For example, the Percus-Yevick (PY) closure is known to work well for hard repulsive potentials, such as the Lennard-Jones potential or the Hard Sphere potential. 25 Soft repulsive potentials, which are the typical potentials of highly CG descriptions, are well represented by the HyperNetted-Chain (HNC) closure, which ignores the bridge function. 60 Note that while the HNC closure is the approximation of choice of the bridge-function when the potential is soft and long-ranged, other closures are known to work better when the potential is short-ranged and more repulsive. This has to be accounted for when the granularity level of coarse-graining is small, just a few atoms are grouped together, and a fine-graining model is developed. Thus the CG potential, using the HNC closure is evaluated as bb Ueff (r)/kB T = ≠ ln [hbb (r) + 1] + hbb (r) ≠ cbb (r) .

(6)

Direct inspection of Eqs. 5 and 6 shows that the potential contains corrections to the PMF, which represent three-body and higher-order interactions. In this way the effective CG potential of Eq.6 correctly accounts for the many-body interactions due to the propagation of the forces through the liquid, because the high-order interactions are properly projected onto the effective pairwise potential of Eq.6, which can be conveniently used as an input in MD simulations.

15



2.2

Page 16 of 68

Partition function and pressure

The effective Hamiltonian of the IECG description is defined as the sum of kinetic energy (K), intramolecular (U intra ), and intermolecular (U inter ) potentials as H = K + U intra + U inter , where U

inter

=

nb X nb n≠1 X n X X i

j>i

“

(7)

(“–)

bb Ueff (rij , flch , T, N ), and the intramolecular potential contains

–

contributions from bond, angle, dihedral angle, and long-range interactions. The configuration integral of the CG system may be written as 61,62 Z=

where dr =

nb n Q Q

i=1 “=1

Z

···

Z

e≠— (U

) dr,

(8)

dr“i , and the corresponding canonical partition function is given by Z

Q= where

intra +U inter

= 2fi—¯ h2 /mb

1/2

3n nb n!

(9)

,

, the mean thermal deBroglie wavelength. The pressure is calcu-

lated directly from the configurational integral as P = kB T

✓

ˆ ln Z ˆV

◆

,

(10)

n,T

where P and V are the pressure and volume of the system, where the density can be assumed to be either an implicit (passive) or an explicit (active) variable. 36,63 In either cases, the theory depends on only one non-trivial parameter, namely the direct correlation function at zero wavevector (k æ 0, c0 ), in contrast to more traditional CG approaches that rely on

optimizing several adjustable parameters (see Section 5.2). 17,64–66

16


Page 17 of 68 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


3

Analytical solution of the CG potential for macromolecular liquids

One advantage of the IECG method is that it allows for the formal analytical solutions of the CG potential. When the granularity of the CG model is large enough, the statistical position of CG sites along the chain becomes pair correlated, and the IECG model can be further simplified by adopting Markovian statistics. Clark and coworkers derived in this limit the analytical forms for both the intramolecular and the intermolecular distribution functions, which have been shown to be in remarkable agreement with atomistic simulations. 21,34,56 From the direct correlation function of Eq.3, the interaction potential is easily calculated by solving the HNC closure in Eq.6. At melt densities, the HNC equation, (Eq.6), is accurately simplified into the Mean Spherical Approximation (MSA), defined below, which is obtained from the HNC by linearization of the logarithmic fuction. 67 The MSA equation was originally proposed by Lebowitz and Percus and then analytically solved for a number of specific systems. 25 Both the HNC and the MSA closure equations have been shown to give consistent results for the EOS calculated following different thermodynamic routes for a number of test systems. When a polymer is coarse-grained, for example at the soft sphere level, the potential becomes a smoothly-varying function of the distance, which slowly decays within a range ‡. At liquid density, the function flch ‡ 3 >> 1, and the average interparticle distance is a = flch , so that ‡ >> a. Thus, each CG site is interacting ≠1/3

in this regime with a very large number of surrounding CG sites, supporting a mean-field approximation to the HNC closure. 60,67 It is worth mentioning that the word mean-field here indicates a weak, long-range CG potential, which is in some aspects analogous to a Kac potential. 67 Thus, the excess free energy of the system can be well approximated by the RR drdrÕ U cc (|r ≠ rÕ |)flch (r)flch (rÕ ), where U cc (r) is the mean-field equation Fex [flch (r)] ¥ 1/2

CG potential and flch (r) is the position-dependent density profile. Given that the potential

is a slowly decaying function and the liquid is homogeneous, in the range of the potential, 17



Page 18 of 68

‡, the density can be considered constant or ccc (|r ≠ rÕ |; flch ) = ≠

lim

flch (r)æflch

ˆ 2 —Fex [flch (r)] , ˆflch (r)ˆflch (rÕ )

(11)

thus recovering the MSA, 68 where U cc (r) ¥ ≠kB T ccc (r) .

(12)

The MSA is well justified in the regime where most of our calculations apply, namely at the density of a polymer melt where the system is almost incompressible, and when the potential is a long-ranged, slowly-varying function of the interparticle distance, which is the shape typical of highly CG models. In those conditions the IECG potential becomes analytically solvable. Considering the multiblob models, for large separations in real space, where r/Rb >> 1, the inverse transform integral of the direct correlation function is sufficiently dominated by the small wave vector limit that the large wave vector contribution can be entirely neglected. Since the expansion for small wave vectors is bounded at large k, the error incurred in using the small k form with the integral bounds extended to infinity is small. This approximation along with the assumption of Gaussian statistics for the intramolecular correlations leads to 21,34,56 cbb (r) ¥ =

⇣

h

RŒ ≠Nb b k sin(kr) 45+45b k4 3 r 0 2fi 2 flm Rgb ✓ Ô ◆ 1/4 ⇥ Õ r) Õ 45 2Nb b Ô Ô e≠Q r ≠ 8fi 3 4 5fl R3 sin(Q QÕ r m gb +

✓

where QÕ = 51/4

◆

r sin(QÕ r) +

≠1/4 b

and Q © QÕ

945+13Q4

p

1/4 b

3/2

+ +

i⌘

5k2 13 b k4 ≠3780 dk 4 2 b k +45) ◆ ✓28 ( Ô ⇥ 5Nb (13Q3 (QÕ r 1/4 3 672fiflm b Rgb

≠ 4)) cos(QÕ r)

1/4 b

⇤ ≠QÕ r ⇤ cos(QÕ r) e QÕ r ,

1/4 b

. The key quantity of interest is

945r

(13) b

= Nb flm |c0 |,

which is defined once the molecular and thermodynamic parameters are known, but also depends on the direct correlation function at k = 0, which has to be determined. 18


Page 19 of 68 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 function

b,

which defines the decay of the repulsive part of the potential, increases

with decreasing CG granularity. It is also a function of the liquid density, explicitly through the parameter flm and indirectly through the density dependence of |c0 |. As will be discussed

later, the degree of polymerization also enters the function

b

indirectly through the param-

eter |c0 | (see Figure 8), however it is important to notice that |c0 | does not depend on the granularity of the CG model adopted.

Interestingly, a simple thread model representation of the chain statistics, 54,55,57 which is the conventional model used in polymer field theory, gives an incorrect density dependence of |c0 |, as noted in McCarty et al. 22 when compared to United Atom simulations of

polyethylene melts. Of course the thread limit is a highly idealized model, and the correct density dependence for the thermodynamic properties is predicted for realistic chain models. This case highlights an important point that using an inappropriate monomer-level theory, or incompatible closure relations, as the starting point for a CG model cannot be expected to yield accurate thermodynamic predictions. Interestingly, the thread model performs well in the description of polymers in solution, correctly reproducing deGennes scaling laws for the mesh length, osmotic pressure, and more. 69

3.1

Some considerations on the properties of CG potentials

From the direct inspection of the analytical equation for the IECG potential a number of observations emerge, which are general enough to apply to the potential of any possible CG model. Here is a summary of the observations: • The range of the potential increases with the number of atomistic sites that are averaged

into the CG unit, so that the more extended the level of CG, the longer ranged the potential. Interestingly, the potential doesn’t vanish even for infinitely long chains. This observation

is important because the range of the potential is the quantity that ultimately defines the computational speed-up of the CG model with respect to the atomistic MD simulations. 44 • The potential has a long-ranged, slowly decaying repulsive component and a second 19



attractive part that is smaller in absolute value than the repulsive part. This attractive contribution is important when one evaluates the thermodynamic properties of the system and cannot be discarded. Truncating the potential to the short-ranged repulsive contribution does not lead to consistent thermodynamic observables, unless a correction for the tail of the potential is added a posteriori. Higher order terms, which are present in the equation of the potential, tend to give increasingly more negligible contributions, so that including the first attractive contribution, while neglecting higher order terms, gives quantitative agreement of both pressure and structural distribution functions between the IECG and the atomistic representations. 22,23,36 • The attractive contribution is present in the effective potential even when the inter-

molecular atomistic potential, from which the CG potential is derived, is purely repulsive (hard-sphere potential). This indicates that the attractive contribution to the intermolecular potential is, at least in part, a consequence of coarse-graining and of the many-body propagation of the interactions through the liquid. Being the resultant of the projection of many-body interactions onto the pair of CG units, the attractive component of the potential is, at least in part, entropic in nature. Interestingly, the PMF, which is missing the many-body originated attractive component of the potential, does not correctly reproduce the thermodynamics of the system. Note that while the IECG potential described here is specific to CG polymer liquids, the overall shape and scaling behavior of the CG potential are general: one should expect that similar features emerge as a consequence of coarse-graining in any CG model of molecular liquids. Because the IECG model in its present version uses a high level CG representation, with low granularity, the characteristic features of the CG potential are enhanced while a less pronounced behavior should be expected in high-granularity CG models that represent only a few atoms in the CG site. 13,70–75

20


Page 20 of 68

Page 21 of 68 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

Accuracy in CG models: the problem of predicting consistent structural and thermodynamic properties

The first requirement for building a useful CG model is that the simulation that uses the CG description needs to be able to predict with accuracy atomistic physical quantities on the CG lengthscale and larger. Given that atomistic simulations use force-fields that have been parameterized to reproduce experiments, it is safe to assume that CG descriptions that accurately reproduce properties from atomistic simulations, are also in agreement with the experiments. Thus, the consistency between atomistic and CG descriptions also implies reproducibility of experimental data in the CG description. The difficulty in predicting with accuracy atomistic properties from CG simulations is rooted in the process of coarse-graining. During the process of coarse-graining, a number of units are defined as CG sites. These can be selected as a subensemble of the total position coordinates, or can be just the resultant of a direct combination of atomistic units into an effective CG unit, for example a unit centered in the center of mass of a group of atoms like in the IECG model (see Eq.1). In both cases, the consequence of averaging over the coordinates of some atomistic particles leads to a mean-field description of the CG units in the field of the atomistic ones that are now eliminated from the description. By coarse-graining the atomistic degrees of freedom, the local scale information is lost, and the potential between CG units becomes a free energy of the coordinates that are eliminated. Thus, the resulting effective potential becomes a function of the molecular and thermodynamic parameters that define the system. For example, in a molecular liquid that contains n molecules each having N atoms, with each atom i positioned at the coordinates ri , at the monomer/atomistic resolution the partition function is defined as QAT =

1 3n N n!

Z

dr1

Z

dr2 · · ·

21

Z

drnN e≠—UAT ,


(14)


Page 22 of 68

where the atomistic potential contains both intra and intermolecular contributions as UAT = intra inter +UAT . The partition function for the same system where each molecule coarse-grained UAT

into a number nb of blobs is defined as 1

QCG =

3n nb n!

Z

dR1

Z

dR2 · · ·

Z

(15)

dRnnb e≠—UCG .

The two partition functions are, in general, not equivalent. However, if one search for a CG model that has a partition functions identical to the underlying atomistic model, so that all the structural and thermodynamic properties are the same at any level of CG, one can formalize the problem as follows. One can identify in an ensemble of N atomistic units a subensemble of n atomistic units that are identified as effective CG units, so that the partition function of the CG ensemble, in the field of the averaged-out N ≠ n atomistic particles, is given by QCG =

1 3n n!

with UCG,eff = ≠kB T ln

✓

Z

dr1

Z

1 3n(N ≠1) n!

dr2 · · ·

Z

Z

drn+1

drn e≠—UCG,eff = QAT ,

Z

drn+2 · · ·

Z

(16)

≠—UAT

drnN e

◆

.

(17)

Eq.17 shows that coarse-graining implies an effective CG potential that becomes a function of the properties of the atomistic particles that are integrated out. Given this definition of the effective CG potential, Eq.17, the CG description for this model has the same partition function as the atomistic description with the consequence that all the thermodynamic and structural properties at different resolution are conserved. In fact those properties could be directly calculated starting from the CG partition function. In most cases, however, a CG unit is defined as a particle centered at the center-of-mass of a group of atoms. In that case to recover the same partition function one usually defines the so-called mapping function (M ), which transforms the atomistic description of n molecules

22


Page 23 of 68 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


including N sites, into the CG representation of n centers-of-mass as Z

dRn ”(M (rnN ) ≠ Rn ) = 1 .

(18)

From the definition of the mapping function, it follows that the atomistic partition function can be rewritten as

1

QAT =

3n n!

Z

dRn e≠— (UCG,eff ) = QCG , Õ

(19)

if the effective CG potential is defined as Õ UCG,eff

= ≠kB T ln

1 3n(N ≠1)

Z

dr1

Z

dr2 · · ·

Z

drnN ”(M (rnN ) ≠ Rn )e≠—UAT .

(20)

Eq.20 is formally correct, but cannot be directly implemented in a simulation because the effective CG potential is inherently a many-body potential, which is hard to reduce to the pairwise form needed to perform conventional MD simulations. The importance of higherorder contributions in the CG potential has been proven in a number of studies, which have shown that already the inclusion of three-body corrections to the two-body contributions improves the consistency of the pair distribution function. 76 Furthermore, in this description the effective potential is formalized by enforcing the property that the partition function has to be independent of the level of CG. In practice, when one performs a MD simulation of the reduced representation, one finds that several properties such as the internal energy of the system or its entropy are not conserved. In general Eqs.18-20 are not easily translated into formal functions that can be calculated and used as an input in simulations. Ideally, if these equations were solved and the transformation found, then the properties that depend on the partition function, i.e. all the structural and thermodynamic properties of a liquid, would be equivalent in the atomistic and CG representations, and could be directly calculated through Eq.18. The current method to solve this problem is, for many CG models, to select a trial potential, which

23



Page 24 of 68

is optimized numerically by trying to reproduce all the possible physical quantities, structural, thermodynamic, and dynamic, of the atomistic model. This procedure is repeated until consistency of the reference properties is achieved (see Section 4.2 for a more complete discussion). Different from other approaches, the IECG method is based on the formal relation of Eq.1, which connects the atomistic and the CG radial distribution functions through the distribution of intramolecular atomistic sites contained in one CG unit. In the IECG approach the condition of equivalency between the atomistic and the CG properties is based on the EOS, which is shown analytically to be independent of the degree of CG.

4.1

IECG equation of state

The pressure (EOS) for the multiblobs representation is obtained from Eq.10, which gives a decomposition of the pressure into kinetic and intramolecular and intermolecular contributions as 36

P = Pkin + Pinter + Pintra .

(21)

While the kinetic and intermolecular contributions are defined as Pkin = flb kB T , and Pinter

2fifl2b =≠ 3

Z

Œ 0

r3 g bb (r)

Û bb (r, flc , T, N ) dr , ˆr

(22)

(23)

the intramolecular contribution to the pressure contains both a harmonic (Pbharm ) and a repulsive (Pbrep ) part for the bond interactions, plus a non-bonded interaction contribution (Pnb ) as Pintra = Pbharm + Pbrep + Pnb . 24


(24)

Page 25 of 68 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 contribution from the harmonic bonds can be integrated analytically and combined with the kinetic contributions to give Pkin + Pbharm = flc kB T ,

(25)

suggesting that the kinetic and harmonic bond contributions in the multiblob models can compensate each other. The contributions can also be solved numerically as reported in Table 1 (from ref 23). The table summarizes how different contributions to the EOS compensate each other, giving a constant total pressure independent of the CG granularity of the representation, as it was also illustrated in Figure 1. Table 1: The calculated pressure for a polyethylene melts with N = 300 when the polymer is represented by various number of blobs. In the IECG theory, Pkin , Pinter , Pbharm , Pbrep , and Pnb are the contributions from the kinetic, harmonic bond, repulsive part of bond, intramolecular non-bonded interactions, and internal virial, respectively. The atomistic pressure for this polymer melt is 343 ± 4 atm. Also see ref 23. IECG theory nb ÈPkin Í ÈPinter Í ÈPbharm Í ÈPbrep Í ÈPnb Í ÈPvir Í ÈP Í 1 7.5 335.8 0 0 0 335.8 343.3 2 15.1 335.0 ≠7.5 0.7 0 328.2 343.3 4 30.1 332.2 ≠22.6 2.0 0.3 311.9 342.0 6 45.2 329.6 ≠37.7 2.9 1.6 296.4 341.6 10 75.3 326.6 ≠67.8 4.0 3.6 266.4 341.7 IECG simulations nb ÈPkin Í ÈPinter Í ÈPb Í ÈPnb Í ÈPvir Í ÈP Í 1 7.5 336.0 0 0 336.0 343.5 2 15.1 335.0 ≠6.6 0 328.4 343.5 4 30.1 331.9 ≠19.2 0.4 313.1 343.2 6 45.2 328.9 ≠32.3 1.6 298.2 343.4 10 75.3 325.3 ≠59.1 4.5 270.7 346.1 The total pressure for the polymer partitioned in 2, 4, 6, and 10 blobs is consistent, within the numerical precision of the atomistic simulations, with the soft sphere representation (nb = 1), which is given by the intermolecular contribution to the pressure as P 2fiflc =1≠ flc k B T 3kB T

Z

Œ 0

25

r3 g cc (r)

Û cc (r) dr . ˆr


(26)


Page 26 of 68

Eq. 26 can be solved in the MSA of Eq. 12, which implies g cc (r) ¥ 1 in the melt density

(see ref 36 for details about this approximation), yielding the solution P flm N c 0 nb =1≠ =1+ flch KB T 2 2 where the EOS is written in terms of

b

b

(27)

to stress that for a relatively high level of CG the

intramolecular contributions almost balance the intermolecular ones and results in an overall pressure that is independent of the CG description. Table 1 shows how different parts of the virial contribute to the final pressure, which almost remains constant in a range of CG that bridges from the soft sphere to the multiblob descriptions. Note that on a very local scale, which is beyond the scope of the current IECG formalism, the contributions of the intramolecular terms become relatively complex and the MSA of the intermolecular potential becomes incorrect. The equivalence in the EOS, shown here, does not imply equivalence in the partition function, and, in fact, the IECG approach shows that the internal energy is not conserved, but the excess free energies along an isotherm is kept constant. 22 The consistency of Eq.27 is fulfilled at any level of CG because, at any given state point, the parameter c0 is independent of the granularity of the CG representation. This is reasonable considering that the parameter c0 is defined in the atomistic representation of the polymer liquid, and the latter is the reference system of any IECG representation, independent of the degree of CG adopted in the given reduced description. Because the EOS is consistent across different levels of CG, this implies that the excess free energy is also conserved. In this way, because in the IECG approach the EOS is conserved during coarse-graining, not all the thermodynamic quantities that are calculated are consistent with the atomistic ones. The granularity dependence of the partition function can provide the corrections to the properties measured in CG simulations, yielding accurate predictions of all the atomistic properties measured in the CG-MD simulations (see Section 6 for a more

26


Page 27 of 68 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


detailed discussion).

4.2

Challenges in the numerical search for the optimal CG potential

While the IECG approach has a formal solution of the CG potential for molecular liquids, Eq.6, for other CG approaches only a numerical solution of the potential is possible. In those CG models, the functional form of the potential is not easy to formalize, and one has to resort to a trial potential that is then refined self-consistently by matching atomistic data. The trial potential is a free energy and as such is parameter dependent. In general, the trial potential contains a considerable number of variables, rendering the numerical derivation of the effective CG potential computationally non trivial. Thus care needs to be taken in selecting the proper starting potential and in selecting the proper minimization procedure to be used in the parameter optimization. If the CG model allows for a non-uniform level of CG, where different parts of the same molecule can be described with different degrees of resolution, selecting the proper multiresolution CG model requires testing all the possible choices of how to group the atoms, which grows exponentially with the number of units one wants to group. The goal of this nonuniform coarse-graining is, for example, to describe the dynamics of proteins by recognizing that their secondary structure defines regions that move in a coordinate way and that could be represented by one bead, with the purpose of most accurately reproducing atomistic data or experiments. The number of possible ways atoms in a structure can be grouped together becomes rapidly prohibitive, rendering the numerical problem of finding the optimal multiresolution, parameter-dependent, CG description even more computationally challenging. 77 Finally, while the optimization of the trial function can be performed using different methods, for example with the relative entropy method 28 or using machine learning techniques, 26 the success of the optimization method depends on the appropriate selection of the initial trial function representing the potential. If the trial function is unphysical, it 27



is unlikely for the optimization method to be able to converge to the proper function. An example of this situation was shown starting from a version of the IECG approach where the incorrect selection of the initial potential was found to lead, through the application of the relative entropy method, to an optimized potential that doesn’t match the correct IECG potential, thus compromising the accuracy of the IECG predictions regarding consistency in structure and thermodynamics. 78 The traditional choice of the trial potential is the PMF: we argue that unless the system is a liquid at low density or a dilute solution of a large molecule, the PMF is a poor choice as a starting function (see Section 5). Summarizing, while effective numerical optimization procedures lead to potentials that reproduce with accuracy the physical properties used as reference in the optimization scheme, remarkably there is no numerically-opimized CG model that is able to reproduce with accuracy all of the structural, thermodynamic, and dynamical properties of the system. When the model has a very low granularity the error is small a a better agreement is observed (see Section 3.1). 70,72–75,79 Even in the case of the analytical solution of the IECG model, not all the physical properties are consistent with the atomistic values, specifically the properties that are a function of the degrees of freedom are not conserved during coarse-graining. 7,22,23 Examples of a successful low-granularity CG model are the United Atom model, where two hydrogens and one carbon are modeled as one CG unit, or the MARTINI model: in these cases the analysis of the IECG potential (see Section 3.1) suggests that the inconsistency between predicted CG properties and atomistic properties that depend on the degrees of freedom is smaller the higher the resolution in the CG model, indicating that high granularity CG models are viable and useful options for coarse-graining. 70,73–75 Those models, however, allow a relatively small improvement in the simulation timestep, and yield a smaller computational gain than more coarse-grained models (see Section 8).

28


Page 28 of 68

Page 29 of 68

5

CG potential, potential of mean force, and manybody effects

In the process of minimizing a trial function, care is required to understand the role of the effective CG pairwise potential and how this potential differs from the PMF. Effective CG pair potentials are in principle designed to reproduce key underlying atomistic simulation results and thus they need to implicitly capture many-body contributions to the two-body potentials. In the IECG approach this is done through the closure equation, Eq.6. On the other hand, the PMF is defined for pair interactions in the mean field of the surrounding molecules (see Eq.5), and thus ignores all the higher-order contributions to the two-body effects. 1 nb = 2 nb = 4 nb = 6

0.002

bb

U (r), kcal mol

0.6

0

bb

-1

-1

0.8

U (r), kcal mol

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


0.4 -0.002

0.2

0

60

120

180

240

300

r, Å

0

0

100

50

150

r, Å

Figure 3: Comparison of the CG potentials obtained from the HyperNetted-Chain (HNC) closure (solid lines) to the PMF (dashed lines) for polyethylene melts with degrees of polymerization of 192 at a monomer density of 0.03296 Å≠3 and a temperature of 503 K: the polymer melts are represented by two, four, and six blobs. The inset shows the long-range behavior of the potentials. As an example, Figure 3 compares the HNC and PMF potentials for polymer melts at a monomer density of 0.03296 Å≠3 with degrees of polymerization of N = 192 when they are represented as two, four, and six blobs. The ranges of both potentials increase as the level of CG increases, while the PMF misses the entropy-generated attractive well at large distances depicted in the inset. The attractive contribution is generated by the many-body 29



interactions that are not included in the PMF but are included in the effective CG potential. Additionally, note that at this monomer density the CG potentials are in general more repulsive and more long-ranged than the corresponding PMFs at any given level of CG. The IECG theory is based on the statistical mechanics of molecular liquids, thus one expects that the CG potential becomes more similar to the PMF at low densities because many-body contributions become less important. Figure 4 shows how the CG potentials obtained from the HNC and MSA closures as well as the PMF change as the density changes. The left panel of the Figure compares the PMF and the CG potentials obtained from the HNC and MSA closures at relatively low densities. Indeed, the figure shows that there is less difference between the CG potentials and the PMF at low densities, while the deviations become increasingly more pronounced for higher densities as shown in the right panel of Figure 4. U (r), kcal mol × 10

0.6

-1

HNC MSA PMF

0.8

-3

ρm = 0.03296 Å

0.4

bb

0.2 0

10

20

30 r, Å

bb

40

120

150

50

60

1.5

0.15

1.2

-3

ρm = 0.03125 Å

-3

ρm = 0.03439 Å

0.9

0.1

0.6

0.05 0

90 r, Å

-1

-1

0.2

60

bb

0

2 1 0 -1 -2 30

U (r), kcal mol

bb

U (r), kcal mol

-1

-3

1

U (r), kcal mol

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 30 of 68

0.3 0

10

20

30 r, Å

40

50

60

0

0

10

20

30 r, Å

40

50

60

Figure 4: Top: the CG potential obtained from HyperNetted-Chain (HNC) closure (black line) is compared with the one obtained from the mean spherical approximation (MSA) closure (red line) as well as the potential of mean force (PMF) (blue line) at a polyethylene monomer density of 0.03296 Å≠3 and a temperature of 503 K. The inset shows the long-range behavior of the potentials where HNC and MSA superimpose. Left bottom: comparison of the potentials at a monomer density of 0.03125 Å≠3 . As the density decreases the PMF deviates less from the CG potentials obtained from the HNC and MSA closures and the attractive contribution disappears. Right bottom: comparison of the potentials at a monomer density of 0.03439 Å≠3 . As the density increases the PMF deviates more from the CG potentials obtained from the HNC and MSA closures. As discussed previously, the IECG potentials allow one to get consistent results between 30


Page 31 of 68

the CG simulations and the atomistic ones. However, the PMF that ignores many-body effects is not able to reproduce the aforementioned consistencies. Figure 5 illustrates this point, where the RDF obtained from the atomistic simulations are compared with the ones obtained from the IECG potential using the HNC closure and the PMF (see the Appendix for simulation details). As can be seen, excellent agreement is observed between the atomistic simulations and the IECG simulations; however, the correlation hole is significantly overestimated by the PMF because it is less repulsive, likely due to ignoring higher-order forces. In addition, the pressure values that are obtained from the IECG simulations are within the statistical uncertainties of the atomistic simulations, while the PMF results are significantly underestimated. 1

atomistic simulations PMF-CG simulations IECG simulations

bb

g (r)

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


0.8

0.6

0

10

20

30

40

50

r, Å

Figure 5: The radial distribution functions for polyethylene melts with a degree of polymerization of 300 at a monomer density of 0.03296 Å≠3 at 503 K. The polymer melts are represented as four blobs. The radial distribution function obtained from the IECG potential with the HyperNetted-Chain (HNC) closure (squares) is compared with the ones obtained from the atomistic simulations (triangles) and the potential of mean force (PMF) (diamonds). Summarizing, when the pair potential is calculated starting from the PMF, which does not include many-body interactions, the agreement with the physical quantities is not good, but could be improved by including three-body corrections. 17,28,76 However, if the pairwise potential is properly calculated by including the appropriate closure correction to the PMF, three and higher-order body corrections are effectively included in the pairwise poten31



tial, and consistency for structural and thermodynamics properties is preserved (see Figure 1). 7,21–23,36,56,58

5.1

Consistency in the radial distribution function does not imply consistency in the thermodynamics: challenges in using a selfconsistent PMF

As shown in Figure 5, selecting the PMF as the CG potential in a simulation of a molecular liquid does not lead to either the correct pair distribution function or the correct thermodynamics, and gives a description of the system that is different from the related atomistic simulation. The reason is because the PMF does not account properly for the many-body effects that at liquid density are important. This shortcoming of the PMF has been noticed earlier, and has lead to the development of CG models that start from the PMF, and optimize a function of the PMF by iterative methods until self-consistency with the pair distribution function, or a function of the pair distribution function like the force, is obtained. Because the resulting potential is self-consistently optimized to reproduce a physical quantity, usually the radial distribution function, these methods display a good agreement with the functions used to optimize the potential. 9,10,80–82 Still, physical quantities that are not included in the optimization scheme, such as thermodynamic properties like the pressure, are not correctly reproduced. This observation has been illustrated also by Akkermans et al. 83 by means of an inverse Monte Carlo method, where the potential shows that it can consistently reproduce the structure of a liquid of polymers, while it does not reproduce the pressure. This result seems somewhat surprising because the statistical mechanics of liquids directly relates the thermodynamic properties to the structural properties. Specifically, Henderson’s theorem states that the radial distribution function is unique up to a constant, 24 while one of the properties of the radial distribution function is that all the thermodynamic quantities can

32


Page 32 of 68

Page 33 of 68

be directly calculated from this function following the rules of statistical mechanics. 25 One would expect that once g(r) is precisely calculated, all the statistical and thermodynamic properties of the liquid are defined. 600 CG, full cc CG, h (r) cut at 3.5Rg

0

atomistic 500

cc

N100-full N100-3.5Rg N192-full N192-3.5Rg

-0.2

Pressure, atm

N78-full N78-3.5Rg

-0.1

h (r)

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


400

atomistic 300

-0.3

-0.4

0

10

20

30

40 50 r, Å

60

70

80

200

80

120

N

160

200

Figure 6: Left panel: The center of mass total correlation function (tcf), hcc (r), for mesoscale soft sphere simulations of polyethylene melts with N =78 (partially-filled circles), N =100 (partially-filled squares), and N =192 (partially-filled diamonds) is compared with the results from theory (solid lines) and mesoscale soft sphere simulations with the potential derived by cutting hcc (r) at r=3.5Rg , N =78 (unfilled up-pointing triangles), N =100 (unfilled downpointing triangles), and N =192 (unfilled right-pointing triangles). The United Atom results are shown by plus symbols. Right panel: the pressure values for various approaches for polymer melts with different degree of polymerization. Figure 6 shows an example of calculations of the radial distribution functions and pressures for a number of polymer liquids at fixed thermodynamic conditions and increasing degree of polymerization, N . The left panel illustrates that by truncating the radial distribution function at 3.5 times the size of a polymer, the IECG truncated potential is still capable of predicting with high accuracy the structure of the liquid through the radial distribution function, for all the samples considered. This result proves that the radial distribution function is largely insensitive to the tail of the CG potential, and that good predictions for the structure of the liquid can be obtained even when the CG potential is truncated. The panel on the right shows that truncating the potential has dramatic effects on the value of the pressure and thus on the thermodynamic properties of the simulated CG systems. This result suggests that the optimization of the potential by matching solely the radial 33



distribution function, or the related force, cannot ensure a proper account of the many-body contributions to the CG potential, and prediction of the correct equation-of-state. Considering a one-component non-polar liquid, it is known that its structure is defined by the packing of the molecules and thus by the short-range, repulsive part of the potential, while thermodynamic quantities, like the pressure, are sensitive to the long-range part of the potential, particularly to the presence of an attractive contribution. 32,33 Thus the selection of the PMF as the trial function in the numerical optimization procedure, tends to underestimate the contribution due to the long-range part of the potential, and can be, in fact, a poor starting choice in the optimization scheme. The relative importance of the attractive and repulsive components of the potential on the structure and thermodynamics of dense liquids, at the atomistic level, is well-known. Barker and Henderson mapped the repulsive part into a short-ranged potential that accurately reproduces the structure, while the attractive component was treated as a perturbation to the short-range repulsion. 84 Widom has shown that structural properties close to the triple point of the phase diagram depend on the hard-core repulsive part of the potential, while in the proximity of the critical point the liquid displays long-range fluctuations, which are dominated by the attractive tail of the potential. The compressibility, which diverges at the critical point, is dominated by the long-ranged part of the potential and so is the liquid pressure. 32,33 Weeks, Chandler, and Andersen’s approach separates repulsive from attractive forces, leading to a perturbation series that is rapidly converging and can be truncated after the first-order perturbation. They showed that neglecting the attractive part of the potential has a small effect on the pair distribution function of medium density liquids and negligible effect at melt-like densities. 85 The need for higher order corrections to the radial distribution function can be seen by direct comparison between the equation for the CG potential, Eq.6, and the PMF, Eq.5. In the framework of liquid state theory, and of the IECG approach, the PMF represents interactions before many-body effects are included through the closure equation: the closure 34


Page 34 of 68

Page 35 of 68 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


effectively projects the many-body interactions into a CG pairwise potential. It is not surprising that if a proper closure is applied, the resulting pairwise potential correctly reproduces both structural and thermodynamic properties. 21–23,36 In agreement with this statement is the work by Akkermans and Briels 83 who argued that “many-body interactions cannot be excluded” in a proper description of a liquid. Later, Zwicker and Lovett 86 showed that the typical radial distribution function calculation fails to fix a unique thermodynamic state, but that if the appropriate closure relation to the integral equation is selected, that correct closure relation “fixes not only the distribution functions, but also the complete thermodynamic state of the system.” In the same paper, the correctly-closed integral equation is shown to possess a “unique solution.” The quest is then for a CG approach of molecular liquids to find the appropriate closure that brings consistency to both structural and thermodynamic properties. The selection of the correct closure makes the need of higher order corrections to the potential unnecessary, as it leads to the proper effective two-body potential. Finally, another typical source of inconsistency in numerically optimized potentials that use the radial distribution function as the test quantity, and thus the associated PMF as the trial potential, is due to the truncation of g(r) in the data, whose effects are illustrated in Figures 5 and 6. Those can be atomistic simulation data that are truncated due to the finite size of the box, or experimental data that are truncated at large distances (small k). Missing the large-scale contributions to g(r) can have large effects on the quality and ability of correctly predicting thermodynamic properties of the resulting optimized effective potential.

5.2

Density dependence of the CG potential

When potentials are density dependent, as is the case for CG potentials, the density dependence leads to two different “routes” that can be used to calculate some thermodynamic quantities, such as the pressure. 87,88 If the density is assumed to be a “passive” variable, the 35



equations are consistent with the ones used for density independent potentials. If the density is considered as an “active” variable, then thermodynamic quantities that are calculated as a function of the density, contain extra terms in their definition. This issue was formally discussed by Stillinger, Sakai, and Torquato. 63 The two routes, passive and active, lead to different values of a number of thermodynamic quantities even in the canonical ensemble. For example, even in the canonical ensemble, where n, V , and T are constant, the virial equation for the pressure, which in the “passive" description is defined as in Eq.26, when the density is assumed to be an “active" variable is modified and follows Eq.28. Thus, in the case that the density is considered as an active variable, the calculation of the virial pressure conventionally used in MD codes, would need to be accordingly mofidied from the traditional equation (Eq.26). Most conveniently, in the framework of the IECG approach, considering the density as an active variable amounts to just a trivial rescaling of the one free parameter, c0 , which is accordingly redefined (see Eq.29). In this way, consistency of the pressure can be ensured without need of redefining the virial calculation of the pressure in the MD simulation. Once the parameter c0 is defined to reproduce the atomistic pressure, MD simulations of IECG models correctly agree with the equation of state of the underlying atomistic system, even when the granularity of the model is varied as illustrated in Section 4.1. Note that c0 in both passive and active (ceff 0 ) routes is density dependent and care is required to use them in the IECG method to obtain the thermodynamical and structural consistency. 36 Specifically, for coarse-graining up to the soft-sphere level the resulting pressure is calculated following the method described in Section 4.1, while maintaining the density-dependence of the formalism as 23

36


Page 36 of 68

Page 37 of 68 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


Z Û cc (r, flch , T, N ) P 2fiflch Œ 3 cc =1≠ dr r g (r, flch , T, N ) eff flch kB T 3kB T 0 ˆr Z Û cc (r, flch , T, N ) 2fifl2ch Œ 2 cc r g (r, flch , T, N ) eff dr, + kB T 0 ˆflch

(28)

cc where g cc and Ueff are the radial distribution function and the effective density-dependent

potential between the soft spheres. The analytical solution of Eq.28 with the MSA closure yields the EOS for the soft sphere CG model P flm N ceff 0 (flm , T, N ) =1≠ , flch kB T 2

(29)

where ceff 0 (flm , T, N ) = (c0 (flm , T, N ) + flm ˆc0 (flm , T, N )/ˆflm ). The inclusion of the density dependence in the IECG formalism effectively corresponds to a renormalization of the direct correlation function at k = 0 as ccc eff (k) = ≠ where

eff

=

+

fl

and

N eff flm 1 +

= ≠N flm c0 and

(k)2 eff ( mm (k) ≠ cm

fl

cm (k)2 )

,

(30)

= ≠N fl2m ˆc0 /ˆflm .

The EOS in the density dependent description results in eff P =1+ , flch kB T 2

(31)

which has the same form as the density independent formalism, where the density is a “passive” variable, but with a renormalized

parameter. Similarly, the EOS for multiblobs

may be obtained from Eq.10 and one can define

eff b

following the discussions after Eq.27.

Once the effective ceff 0 (flm , T, N ) parameter is optimized against the atomistic data, the IECG approach recovers both structural and thermodynamic properties of the molecular liquid. Note that in the calculation of the pressure, a correction is applied a posteriori to recover the correct pressure, accounting for truncation of the effective CG potential in simu37



lations. 1 Therefore, in the IECG approach at melt density, the long-range forces originated from the coarse-graining procedure are necessary for thermodynamical properties. 23

6

IECG corrections for physical quantities that depends on the local degrees of freedom: the case of diffusion

While physical quantities that depend on the number of degrees of freedom, like entropy and internal energy, are bound to change as a result of IECG coarse-graining, if the CG method is based on a rigorous statistical mechanic approach, like ours, it is possible to evaluate the corrections that need to be applied to recover the correct atomistic quantities from the modified properties measured in CG simulations. A point that needs to be considered with care is the fact that the process of coarsegraining eliminates internal degrees of freedom, inevitably modifying some of the physical properties of the system, which become different from the atomistic description. 89,90 Thus, a numerical optimization procedure that aims at reproducing in the CG simulation all the properties of the underlying atomistic description is likely bound to be unsuccessful. For example, by reducing the degrees of freedom, CG changes the entropy of the system. When the CG model conserves the free energy at variable granularity, the internal energy becomes a function of the CG resolution. In this way, by deriving a reduced description of the molecule, one should expect to modify, as a consequence, both the internal energy and the entropy of the system. 91 The latter is the mapping entropy, which needs to be distinguished from the information entropy (or relative entropy), which monitors the loss of information due to a non-exact mapping by coarse-graining. 28,92 In fact, in the IECG model, the loss of information entropy is zero, while the mapping entropy is different and changes with the level of CG. 22

38


Page 38 of 68

Page 39 of 68 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


Another physical quantity that is predicted to be different in CG simulations from the atomistic description is the diffusion coefficient. The dynamics of the CG system is accelerated because both the friction coefficient and the free energy landscape are modified by coarse-graining. The free energy landscape is smoothed out because local energy minima are eliminated by averaging over the local degrees of freedom, thus reducing the number of microstates that are available and modifying the entropy of the system. 41,42,89,90,93 Furthermore, the site friction coefficient is a function of the particle hydrodynamic radius, which is different depending on the level of CG given that the surface of the molecule exposed to the solvent changes with the granularity of the model (a soft sphere representation has a smaller surface exposed to the “solvent” than in the atomistic description). Both factors together modify the diffusivity of the molecule, in a way that depends on the granularity of the model. Thus, CG models that aim at reproducing structure, thermodynamics, and dynamics of the underlying atomistic description are bound to fail at reproducing with accuracy all these properties by one CG potential. Because the inaccuracy of the model increases as the granularity of the CG model decreases, models that only group together a few atoms are less prone to error than larger granularity models, as it was discussed in Section 3. While much progress has been made in this direction, still CG theories that only account for either the friction correction or the entropy correction struggle to make accurate, quantitative predictions of the dynamics in agreement with the atomistic simulations. 89,90,93 The IECG approach based on the thread model has been used to reconstruct the atomistic diffusion from the dynamics measured in IECG simulations of soft spheres. 39 The dynamics in this soft-sphere simulation is many orders of magnitude faster than its atomistic counterpart; however, starting from a Mori-Zwanzig projection operator formalism, 43 Lyubimov et al. were able to recover with unprecedented precision the diffusion coefficient of the atomistic simulations by applying the calculated IECG correction coefficients to the dynamics. 39,40 It is important to notice in this regard, that atomistic simulations are not used to derive the numerical correction factor to rescale the diffusion. Instead these factors are calculated from 39



the solution of the memory function, which gives the friction, and from the calculation of the change in entropy of the intramolecular potential. The solution of the memory function proceeds through a number of conventional steps, which include i) the substitution of the projected dynamics with the real dynamics, which is correct in the long-time limit of diffusion, ii) (Kirkwood-like) factorization of the four-point correlation functions into products of pair distribution functions, which are analytically solved in the IECG approach, and iii) assumption that the center-of-mass diffusion is the leading contribution in the dynamical relaxation. The friction coefficient is calculated from the memory function for both the atomistic and the CG descriptions: their ratio defines the correction factor that needs to be applied to recover from the friction measured in the CG simulation the friction of the atomistic description. Note that the atomistic simulations of liquids of short chains are used only to define the effective diameter in the hard-sphere model that is used in the memory function calculation of the monomer friction coefficient. 39,40 Reconstruction of the atomistic diffusion coefficient with this method yields results that are in quantitative agreement with atomistic simulations and with experiments 94 (see Figure 7). Both terms are solved analytically using the time-dependent structure factors and pair distribution functions from integral equation theory in the atomistic (PRISM) 95–97 and CG (IECG) 7,22,23 resolutions. Thus they provide the analytical expressions for the correction factors that need to be applied to reconstruct from the mesoscale CG simulation the atomistic quantities. To prove that the method is fully predictive, Lyubimov et al. also calculated atomistic diffusion coefficients for new systems for which experimental data or atomistic simulations were not available. 40 Thus, the analytical solution of the IECG theory has the advantage of identifying, without ambiguity, whether the lack of consistency is a consequence of the application of the proper procedure of coarse-graining or is a consequence of numerical error generated by a poor optimization scheme, for example due to an unfortunate selection of the initial trial function in the numerical solution of a CG potential as discussed earlier. Furthermore, the IECG 40


Page 40 of 68

Page 41 of 68

experiment IECG-MD rescaled

UA MD IECG-MD rescaled 2

2

10

-1

2

N

2

-1

N

D, Å /ns

10

D, Å /ns

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

1

10

10

-2.1

cis-polybutadiene

N

-2.5

N

polyethylene

0

0

10

10 2

10

3

N

10

1

10

2

10 N

3

10

Figure 7: Left panel: center-of-mass diffusion coefficient as a function of degree of polymerization, N , for cis-1,4-polybutadiene melts in the unentangled and the entangled regime. Diffusion coefficients reconstructred from the IECG-thread molecular dynamics simulations (open symbols) are directly compared with data from United Atom simulations (filled squares), and show remarkable agreement (See ref 40). Right Panel: diffusion coefficients for polyethylene melts as a function of degree of polymerization, N . Comparison between the theoretically predicted values (triangle) and experiments (circle) from ref 97 and references therein. Also shown is the scaling for unentangled N ≠1 (dot-dashed line) and entangled N ≠2.5 (dashed line) systems. method formally describes how properties are modified when the granularity of the CG model is changed across the different levels of resolution, and provides accurate information on when the consistency across the different levels of CG breaks down and how to correct for this loss of information.

7

Transferability of CG potentials

It is often observed that CG potentials that are calibrated to reproduce properties at a given state point are not able to reproduce the same property at different state points. Furthermore, the accuracy of the CG parameters optimized for a set of thermodynamic states is only limited to state points close to the calibration state and do not apply to other thermodynamic conditions. 27,98 This is known as the transferability problem. In the IECG method, the CG potential is transferable due to its dependence on the accurate EOS and 41



principles of polymer physics. 7,21–23,36,56 The top panel of Figure 8 reports values of the direct correlation function at k æ 0, c0 , and

of the mean-square end-to-end distance ÈR2 Í for a number of polyethylene samples: the values

of the parameters fall on functional forms that are easily extrapolated, making it simple to evaluate the parameters for new systems that are not simulated at the atomistic level. The bottom panel of Figure 8 shows that the density dependence c0 , follows a Carnahan-Starling equation of state as a function of density, making the IECG method transferable to systems

-10

20000

-15

15000

-20

10000

-25

5000

2

25000

2

, Å

-5

eff

c0 , Å

3

at different densities. 99

-30

0

200

400

N

600

800

0 1000

0 N=192

3

-10 eff

c0 , Å

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 42 of 68

-20 -30 -40 0.0315

0.033

0.0345 -3 ρm, Å

0.036

2 Figure 8: Top panel: the values of ceff 0 (red circles) and ÈR Í (blue squares) computed for polyethylene melts of various degrees of polymerization, N , at 503 K and a monomer density of 0.03296 Å≠3 . The red and blue solid lines are fits of the form a + b/N and A + B N , respectively. Bottom panel: the values of ceff 0 computed for a polyethylene melt with N = 192 at 503 K and at various monomer densities. The red line is a fit to a Carnahan Starling-type equation. 99 The value of ÈR2 Í does not change significantly (see the main text).

The IECG parameter c0 is, in general, optimized from atomistic simulations as shown in Figure 8, but it could also be directly calculated from the numerical solution of the PRISM equation 95,96 after adjusting the hard-sphere diameter in the formalism, 22 thus completely avoiding the need to perform atomistic simulations for the states of interest. Alternatively, it can be directly taken from the experimental data of pressure or isothermal compressibility when IECG is used as a top-down approach. 36 Note that c0 in the IECG is related to the 42


Page 43 of 68 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


thermodynamic conditions, for example density and temperature in the canonical ensemble, at which the simulation is performed. The only molecular parameter that enters the IECG potential is the mean-square end-toend distance ÈR2 Í, which is related to the chemical structure of the molecule and its degree of

internal flexibility. This parameter can be taken from experiments or it can be extrapolated to the long chains from the data obtained from a few atomistic simulations of shorter chains, following the statistical properties of polymers. 100,101 Considering the uncertainties in the end-to-end distance, it is worth mentioning that in the range of degree of polymerization studied by us the scaling of ÈR2 Í with N follows the theoretical predictions, while for shorter

chains it will have to be adjusted. Furthermore, at a given degree of polymerization, the change in the values of ÈR2 Í with density at the densities studied so far is small and can

be ignored. For instance, the values of ÈR2 Í for a polymer melt with N = 192 at monomer

densities of 0.03201, 03296, and 0.03343 Å≠3 are directly obtained from atomistic MD simulations as 3625 ± 116, 3653 ± 93, and 3668 ± 135 Å2 , respectively. 36 Alternatively, the values

of ÈR2 Í can be predicted from the principles of polymer physics and the related extensive

work. 102–104

In general, it is not necessary to perform atomistic simulations of the system to coarsegrain, in order to extract the parameters that enter the CG potential at the state points of interest, which is a clear advantage of the IECG method when compared to other methods. To illustrate this point, the transferability of the IECG method is demonstrated for the polymer melts at 503 K at various N in Table 2 and Figure 9. Independent atomistic simulations were performed to show consistencies of the RDFs and the pressure values with the IECG theory and/or IECG simulation results, which did not use any information from atomistic simulations directly. The values of pressure obtained from the IECG theory and IECG simulations in soft sphere and multiblob representations are listed in Table 2 and compared with the results from independent atomistic simulations. Excellent agreement is observed between the IECG theory, IECG simulations, and atomistic simulations demonstrating the 43



transferability of the IECG method. Note that all IECG simulation results in Table 2 are within the statistical uncertainties of atomistic simulation results showing the consistency for various levels of CG. In addition, Figure 9 shows excellent agreement between the RDFs obtained from the IECG approaches and the independent atomistic simulations for N = 44 and N = 124. Table 2: Comparison of the pressure values from independent atomistic simulations and the IECG theory and IECG simulations for polyethylene melts at 503 K and various N . The IECG parameters are computed from the fits shown in Figure 8. The units of the pressure are in atm and the uncertainties in the IECG simulations are less than one atm. N nb IECG theory IECG simulations atomistic simulations 44 1 846 846 849 ± 11 124 1 467 469 466 ± 5 124 2 467 469 466 ± 5 400 1 322 322 325 ± 4 400 2 322 322 325 ± 4 400 4 321 322 325 ± 4 400 8 320 322 325 ± 4

1

1

N124, nb=2

N44 0.9

g (r)

0.8 cc

bb

g (r)

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 44 of 68

0.8 0.6 0.7

0.4

0

0.5

2 1/2

1

1.5

0.6

r/

0

0.2

0.4 0.6 2 1/2 r/

0.8

1

Figure 9: Comparison of RDFs obtained from independent atomistic simulations (filled triangles) and the IECG theory (solid lines) and IECG simulations for polyethylene melts with N = 44 represented by soft spheres and and N = 124 represented by di-blobs (filled squares) at 503.17 K and a monomer density of 0.03296 Å≠3 . The IECG theory/simulation results are within the uncertainties of independent atomistic simulation. Therefore, in the IECG approach, the state-dependence is explicitly introduced through the input parameters of the theory: these have a clear definition, such as density, tempera44


Page 45 of 68 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


ture, end-to-end chain distance, and c0 . As stated above, the compressibility and pressure are closely related to the c0 parameter, which at the monomer-level serves to parameterize the strength of the monomer excluded volume. This can be explicitly seen by considering the second virial coefficient B2 ≥ ≠c0 /2. However, in practice, c0 has a nontrivial temperature

and density dependence. 7,21–23,36,56

8

Computational efficiency

The third important requirement for a CG model is that by reducing the description of the molecule the simulation gains in computational efficiency. Simply put, if a CG simulation affords a very modest improvement in computational efficiency with respect to the related atomistic simulation, the latter is in comparison more useful because it contains the localscale information that is averaged out in the CG description. It is typical in CG models to trade off the local-scale information for an increase in the computational efficiency of the MD simulation. How the speed-up in MD simulations is affected by the degree of CG is an important issue. The consequences of coarse-graining on the computational time are, actually, quite complex because reducing the degrees of freedom has multiple and sometimes competing effects on the efficiency of a MD-CG simulation. In a recent paper, Dinpajooh et al. performed a careful and exhaustive analysis of the factors that contribute to increase efficiency of the CG simulations in comparison with the atomistic ones. 44 In this respect, the IECG formalism is useful in analyzing the effects of CG on the computational efficiency because it provides a rigorous method to monitor the changes in the range of potential as the level of CG changes. Therefore, given that the IECG shows structural and thermodynamic consistency at variable levels of CG, it is possible to analyze how the granularity of the model affects the computational efficiency, through the changes in the shape and range of the potential. Finally, given that the IECG potential has been

45



Page 46 of 68

solved for the highest-levels of CG, where one molecule is reduced to one CG site, all the effects due to coarse-graining are amplified and made clearly visible. In general, the averaging over the local degrees of freedom modifies the pairwise potential such that it becomes progressively more long-ranged and less repulsive at contact, with respect to the initial atomistic Lennard-Jones potential (see Section 3.1). The increase in the range of the potential is the resultant of the contributions due to many-body interactions that propagate through the liquid and are projected onto the pairwise potential (see Section 5). Furthermore, excluded volume repulsion become less relevant with decreasing resolution, as the CG site is located on the center-of-mass of the chain, which is a fictitious site and allows for superimposing of the CG units, while the positions of the atomistic sites are averaged out. Thus, low-resolution CG potentials are bound at zero separation and are soft, long-ranged potentials: the lower the granularity of the CG unit, the softer the potential. The computational efficiency of a CG approach, at a given density and temperature, can be quantified as Á(N, ns ) = Npr

tr –,

(32)

where Npr is the ratio of the total number of pairwise interactions in the atomistic simulations to the ones in the CG simulations,

tr is the ratio of the timesteps used in the CG simulations

to the ones used in the atomistic simulations, and – is the dynamical scaling factor, which corresponds to the speed-up of the dynamics due to coarse-graining. (see Figure 10.) Coarse-graining can have competing effects on the efficiency of simulations. Decreasing the granularity of the CG model leads to an increase in the range of the potential, which demands for a larger cut-off distance and skin distance than in the atomistic simulation, thus requiring an increase of the box size and of the number of pairs of particles that need to be accounted for in the calculation of the intermolecular forces. The computational time is thus increased with respect to the atomistic simulation, i.e. Npr decreases, because at each timestep the simulation needs to calculate the force acting between any pair of particles comprised in the volume defined by the cut-off distance, which scales as the square of the 46


Page 47 of 68 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


Figure 10: Schematic overview of the physical quantities that determine the computational efficiency of a CG method: As the range of CG potential increases, the number of polymers interacting increases while the number of CG interaction sites in a given polymer decreases and one can use larger timesteps in molecular dynamics simulations. This is accompanied by the enhancement of the CG dynamics (increased smoothness of the free energy and reduced friction). Therefore, competing effects exist on the computational efficiency (see the main text). The computational efficiency of polymer melts with N = 300 represented by four blobs is shown in the bottom left considering the aforementioned effects. 44 number of CG units. 44 To overcome this problem, one would be tempted to adopt a short-ranged CG potential, which would decrease the number of pairwise interactions that need to be calculated, thus speeding up the simulation by avoiding a large simulation box. This convenient first approximation needs to be treated carefully, because truncation of the potential can decrease the accuracy of the predictions of the CG simulation. In fact, the tail of the potential is important in the calculations of the thermodynamics (pressure) of the system through the virial equation, while the structure of the liquid is dominated by the short-range part of the potential, as we discussed earlier (see Sections 4.2 and 5.1). If the potential is properly truncated, the tail correction can be added a posteriori, and the pressure can be fixed, after the simulation is completed, through the use of the virial correction. However, if the potential is truncated at too short a distance, not only the pressure but also the structure of the liquid (pair distribution function) is affected by the approximation, and simulation predictions are no longer accurate (see Sections 5.1 and 5.2). As discussed by Dinpajooh et al., a reasonable truncation of the potential would be at distances r that are slightly larger than the distance 47



at which the intermolecular force is zero. 44 This truncation accounts for the repulsive contribution to the potential, plus a small fraction of the attractive contribution. Both repulsive and attractive contributions need to be accounted for to maintain consistency of structure and thermodynamics. Because of the increase in the range of the CG potential, increasing the simulation box requires the calculation of an increasing number of pairwise interactions, thus decreasing the computational efficiency. In some instances, however, this issue becomes less relevant when the type of physical problem that one wants to simulate requires large simulation boxes already at the atomistic level. This is the case when, for example, one wants to study physical quantities that need to be measured at large scale, such as the compressibility of the system or the building up of concentration fluctuations as a system approaches the second-order phase transition leading to demixing. 105,106 These types of problems require as large a box size in the simulation as possible also at the atomistic resolution. Thus the use of CG descriptions not only has little effect on the computational time, but even more, it makes it possible to study processes that cannot be approached by atomistic simulations, because the number of particles involved would be so high that the simulations would become slow to the point that they are computationally unfeasible. The IECG method has the advantage, for example, that it makes predictions for the behavior of polymer mixtures in conditions that cannot be approached by atomistic simulations, while ensuring consistency of the thermodynamic properties at a given state point in the phase diagram. 19,54,55 Competing effects: speeding up MD simulations as a consequence of coarse-graining. One consequence of coarse-graining is that excluded volume interactions are less prominent in CG simulations, and the CG units can superimpose and cross each other, leading to a rapid equilibration of the system and a speed-up of the simulation. 107 This is particularly useful when considering liquids of macromolecules, which can be nearly impossible to equilibrate when the chains have a high degree of polymerization due to entanglements. Once equilibrated at the CG resolution, the atomistic chains can be rebuilt using a backmapping procedure. 64,108 48


Page 48 of 68

Page 49 of 68 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


When a molecule is coarse-grained the CG sites are on average at an intramolecular distance larger than the atomistic bond length, and the bond potential in the CG representation is fluctuating with a lower spring constant than in the atomistic representation. 109 Still the CG bond vibration is the fastest motion in the CG-MD. Thus it defines the time step of the MD simulation,

tr , which has to be captured to correctly simulate the CG dynamics.

When a larger timestep is used in the CG-MD, the simulation allows for the access to longer timescales. If the potential is accurately designed to reproduce the properties of the atomistic description, then the CG-MD is useful in sampling longer timescales in the process of interest. A final computational advantage of CG simulations is given by the acceleration of the dynamics, –, which can reach many order of magnitude speed up. Lyubimov et al. have shown that by coarse-graining the molecular representations two important effects come into play. 39 The first effect of coarse-graining is that the system samples a smoothed free energy landscape, where local energy minima have been averaged. Because of the reduced free energy barriers of the CG representation, sampling of the phase space is accelerated. The decrease in the number of microstates available for sampling is related to the change in entropy due to coarse-graining that was discussed earlier. Another important effect in speeding up the dynamics is due to the change in shape of the molecule when represented by CG sites instead of atomistic sites. In that case, the surface available to the solvent and the hydrodynamic radius of each unit is reduced by coarse-graining, thus modifying the friction coefficient and increasing the diffusivity of the molecule. Those two effects speed-up the dynamics of the CG system and increase the efficiency of the diffusive sampling of the configurational space. When the two contributions are properly accounted for, the CG dynamics can be rescaled to recover the atomistic diffusion. Interestingly, in the IECG approach the structural and dynamical quantities entering the calculation of the memory function, which gives the friction coefficient, and the entropic 49



correction, are known and analytically expressed. Thus the IECG approach affords the direct calculation of these rescaling contributions from the theoretical formalism, without the need of performing first the atomistic simulations to derive the numerical rescaling parameter to the dynamics, as is often done. In this sense, the IECG approach is predictive of the diffusive properties of the system, starting from the atomistic simulations. 39,40 The accelerated dynamics of the CG model is important in speeding up the MD simulations because the dynamical sampling of the configurational space becomes more efficient when the granularity is lower in the CG model. This is an important contribution to the computational efficiency of the CG model: for example, for the soft sphere representation the alpha factor increases with polymer length, and reaches up to eight or nine order of magnitude for the system that was studied (see Table 3). 44 In this respect, it is interesting to note, as mentioned earlier, that in some instances the CG models are developed and parameterized to reproduce the structural and thermodynamic quantities as well as the molecular diffusion of the atomistic description. It seems clear from our discussion, that matching the dynamics of the atomistic model implies forcing the CG model to slow down, losing part of the computational efficiency gained by the adoption of the CG model, which seems a counter-intuitive decision. It is more convenient to harvest the computational gain of the IECG model by taking advantage of the many orders of magnitude of computational speed-up, and to rescale the dynamics a posteriori after the IECG-MD simulation is completed. As an example, the efficiency of a number of IECG-MD simulations are reported in Table 3. 44 This table shows how the IECG method allows one to treat the same polymeric liquid at different resolutions while offering millions of magnitude computational efficiency.

9

Summary and conclusions

One of the advantages of having an analytical formalism of CG, albeit for the specific case of a macromolecular liquid, is that it is possible to address from a fundamental perspective

50


Page 50 of 68

Page 51 of 68 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


Table 3: The efficiency of CG simulations compared to the atomistic ones, when Eq. 32 is used. The values of tr are the ratios of CG timestep to the atomistic timestep. The values of Npr are the ratios of the total number of pairwise forces in the atomistic ones to the ones in the IECG simulations (for more details see ref 44). The Nosé-Hoover thermostat is used for all molecular dynamics simulations. N nb – tr Npr Á 2 44 1 6.5 ◊ 10 900 1.9 1.0 ◊ 106 192 1 2 ◊ 104 3500 0.5 3.3 ◊ 107 192 2 7 ◊ 103 800 2.7 1.4 ◊ 107 192 4 2 ◊ 103 400 11.2 1.1 ◊ 106 3 192 6 1.5 ◊ 10 250 22.4 7.8 ◊ 106 300 1 7 ◊ 104 6000 0.3 1.1 ◊ 108 300 2 2 ◊ 104 1500 1.3 4.0 ◊ 107 300 4 8 ◊ 103 750 4.6 2.5 ◊ 107 300 6 3.5 ◊ 103 500 4.5 9.7 ◊ 106 300 10 2 ◊ 103 300 13.2 9.1 ◊ 106 the problem of how CG affects the shape of the potential, and how this in turn defines the structural, dynamical, and thermodynamical properties of a system once it is coarse-grained. The IECG approach is analytically solved for macromolecular systems at liquid density, but its numerical solution applies to a wider range of possible systems, and its predictions have been extensively numerically tested. By this method we investigated the effects of CG on the shape of the potential, in the form of range, slope, and sign of the tail, and sharpness and range of the short-range repulsive component. Those characteristics affect the structural, thermodynamic, and dynamical properties of the system, as well as the computational gain associated with adopting the CG description. While the IECG method is specific to macromolecular liquids, the properties of CG that emerge from the study of this method are general and apply to any CG model. Thus the discussion and conclusions presented in this paper pertain to the theory of CG in general, and can be used as guiding considerations when developing a formal CG model. For example, one should be aware that a consequence of coarse-graining is that not all the physical properties of the atomistic system can be conserved. Some characteristics of the IECG approach make this method useful to study the theory of CG in general. First of all, the IECG is formally expressed by a Hamiltonian and is defined by 51



a specific partition function, making a direct connection to the theory of statistical mechanics. Because the IECG approach has been shown to conserve both the structural distribution functions, and selected thermodynamic properties (EOS, pressure, and excess free energy) of the atomistic description, while the level of resolution of the CG model is changed, we know that those properties are conserved in the range of granularity that goes from the atomistic description to the soft sphere resolution. In practice, for the properties that are conserved, performing MD simulations of the liquid where molecules are represented either as point particles or at the atomistic level, gives consistent numerical values within the error of the atomistic simulation. The IECG-MD description loses the local-scale information, which can be included a posteriori through a multiscale modeling approach, but speeds up the computational time many orders of magnitude, depending on the size of the molecule coarsegrained, opening up the possibility of studying systems by MD simulations in a window of time and lengthscales impossible to approach by atomistic simulations. Because the method has a strong foundation in the statistical mechanics of liquids, many observations that naturally emerge from solving the CG formalism find their roots in the well-known literature on the theory of simple liquids. This observed parallelism between atomistic liquids and CG molecular liquids, which are expressed, for example, as liquids of soft spheres, further strengthens the impact of the IECG observations. The analytical formalism, furthermore, provides results that are immune from the possible numerical errors of optimization procedures, which are routinely used to derive CG potentials, and leads to the derivation of scaling exponents, which defines the general behavior of CG quantities, with respect to atomistic ones. In this paper we have reviewed and summarized a number of observations that emerged during the course of our studies on the effects of CG on the potential and on the physical properties of the system when simulated in the CG resolution. As a summary of those observations, we identified a number of selected requirements that, in our opinion, need to be fulfilled for a CG approach to be useful, when MD simulations of the CG description 52


Page 52 of 68

Page 53 of 68 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


are performed. These requirements are: accuracy of the predictions, transferability of the potential, and computational speed-up with respect to the atomistic simulations. For all three of them, we were able to identify how the CG potential needs to be constructed to obtain the desired properties. Furthermore, while doing so, we addressed a number of questions that have emerged in the literature on CG models, including for example, if a pairwise potential could ever be capable of predicting structural and thermodynamic properties consistent with the atomistic description; if a CG model is able to reproduce quantitatively all the thermodynamic and dynamic properties of the system, or more succinctly, if there are physical quantities that cannot be reproduced by a CG model; how it is possible to calculate the correction to the properties that are not quantitatively reproduced by the CG model; what are the effects of the truncation of the CG potential on the structural, thermodynamic and dynamical properties of the system, and so on. We hope that this manuscript will be useful to the practitioners of this field.

53



**

A

Simulation details

Here we present the simulation details for the IECG simulations and the CG simulations with the PMF (CG-PMF) of the polymer melts with N = 300 that are represented by four blobs. All simulations were performed with the LAMMPS software program 110 at a monomer density of 0.03296 Å≠3 at 503 K in the canonical ensemble using the NoséHoover thermostat and standard velocity-Verlet integrator. Periodic boundary conditions were applied in all three dimensions. Intra and inter molecular effective IECG potentials, as described in ref 36 were adopted, noting that, in the multi-site CG models, the nonbonded intrachain effective CG potential involves intramolecular CG sites that are separated more than two apart. Cutoff distances of 79.0 Å and 120 Å were used for IECG and CG-PMF simulations, respectively. Timesteps of 30 and 15 fs were used for the IECG and CG-PMF simulations, respectively. The timesteps were chosen noting that they should be smaller than the relaxation times of the fastest modes of vibration, · , which theoretically scales as for these sets of simulations. The combined Verlet neighbor list and the link-cell · Ã n≠1 b

binning algorithm were used to build the Verlet neighbor list, which was updated every 10 steps. 1 The equilibration period ranged from 30-75 ns while the production runs for the IECG and CG-PMF simulations consisted of 30 and 75 ns, respectively.

Acknowledgement This work was supported by the National Science Foundation (NSF) Grant No. CHE1665466. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), 111 which is supported by the National Science Foundation Grant No. ACI1548562. This work used the XSEDE COMET at the San Diego Supercomputer Center through allocation TG-CHE100082. We are grateful to Paula J. Seeger for reading/editing 54


Page 54 of 68

Page 55 of 68 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 manuscript.

References (1) Allen, M. P.; Tildesley, D. J. Computer Simulations of Liquids; Oxford University Press: Oxford, 1987. (2) Frenkel, D.; Smit, B. Understanding Molecular Simulations, 2nd ed.; Academic Press: San Diego, 2002. (3) Ramos, J.; Vega, J. F.; Martínez-Salazar, J. Molecular Dynamics Simulations for the Description of Experimental Molecular Conformation, Melt Dynamics, and Phase Transitions in Polyethylene. Macromolecules 2015, 48, 5016–5027. (4) Shinoda, W.; DeVane, R.; Klein, M. L. Multi-property Fitting and Parameterization of a Coarse Grained Model for Aqueous Surfactants. Mol. Simul. 2007, 33, 27–36. (5) Müller-Plathe, F. Coarse-graining in Polymer Simulation: From the Atomistic to the Mesoscopic Scale and Back. ChemPhysChem 2002, 3, 754–769. (6) Brini, E.; Algaer, E. A.; Ganguly, P.; Li, C.; Rodríguez-Ropero, F.; van der Vegt, N. F. A. Systematic Coarse-Graining Methods for Soft Matter Simulations - a Review. Soft Matter 2013, 9, 2108–2119. (7) Guenza, M. G. Coarse-Grained Modeling of Biomolecules; Taylor & Francis Group, LLC; CRC Press: Boca Raton, 2018; Chapter 2, p 27. (8) Lu, L.; Voth, G. A. In Adv. Chem. Phys., first edit ed.; A., R. S., Dinner Aaron, R., Eds.; John Wiley & Sons, Inc., 2012; Vol. 149; pp 21–40. (9) Faller, R. Automatic Coarse Graining of Polymers. Polymer 2004, 45, 3869–3876.

55



(10) Molinero, V.; Moore, E. B. Water Modeled As an Intermediate Element between Carbon and Silicon â

. J. Phys. Chem. B 2009, 113, 4008–4016.

(11) Takahashi, K. Z.; Nishimura, R.; Yamato, N.; Yasuoka, K.; Masubuchi, Y. Onset of Static and Dynamic Universality among Molecular Models of Polymers. Sci. Rep. 2017, 7, 12379. (12) Salerno, K. M.; Agrawal, A.; Perahia, D.; Grest, G. S. Resolving Dynamic Properties of Polymers through Coarse-Grained Computational Studies. Phys. Rev. Lett. 2016, 116, 3–7. (13) Wang, H.; Junghans, C.; Kremer, K. Comparative Atomistic and Coarse-grained Study of Water: What Do We Lose by Coarse-graining? Eur. Phys. J. E 2009, 28, 221–229. (14) Delle Site, L.; Praprotnik, M. Molecular Systems with Open Boundaries: Theory and Simulation. Phys. Rep. 2017, 693, 1–56. (15) Fritz, D.; Herbers, C. R.; Kremer, K.; van der Vegt, N. F. A. Hierarchical Modeling of Polymer Permeation. Soft Matter 2009, 5, 4556. (16) Fritz, D.; Koschke, K.; Harmandaris, V. A.; van der Vegt, N. F. A.; Kremer, K. Multiscale Modeling of Soft Matter: Scaling of Dynamics. Phys. Chem. Chem. Phys. 2011, 13, 10412. (17) 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. (18) McCarty, J.; Lyubimov, I. Y.; Guenza, M. G. Multiscale Modeling of Coarse-Grained Macromolecular Liquids. J. Phys. Chem. B 2009, 113, 11876–11886.

56


Page 56 of 68

Page 57 of 68 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


(19) McCarty, J.; Guenza, M. G. Multiscale Modeling of Binary Polymer Mixtures: Scale Bridging in the Athermal and Thermal Regime. J. Chem. Phys. 2010, 133 . (20) Guenza, M. G. Advancements in Multi Scale Modeling: Adaptive Resolution Simulations and Related Issues. Eur. Phys. J. ST 2015, 224, 2491–2495. (21) Clark, A. J.; McCarty, J.; Lyubimov, I. Y.; Guenza, M. G. Thermodynamic Consistency in Variable-Level Coarse Graining of Polymeric Liquids. Phys. Rev. Lett. 2012, 109, 168301. (22) McCarty, J.; Clark, A. J.; Copperman, J.; Guenza, M. G. An Analytical Coarsegraining Method which Preserves the Free energy, Structural correlations, and Thermodynamic State of Polymer Melts from the Atomistic to the Mesoscale. J. Chem. Phys. 2014, 140, 204913. (23) Dinpajooh, M.; Guenza, M. G. On the Density Dependence of the Integral Equation Coarse-Graining Effective Potential. J. Phys. Chem. B 2018, 122, 3426–3440. (24) Henderson, R. L. A Uniqueness Theorem for Fluid Pair Correlation Functions. Phys. Lett. A 1974, 49, 197–198. (25) McQuarrie, D. A. Statistical Mechanics; University Science: Sausalito, CA, 2000. (26) Weinane, E. Principles of Multiscale Modeling; Cambridge University Press: Cambridge, 2011. (27) Mashayak, S. Y.; Miao, L.; Aluru, N. R. Integral Equation Theory Based Direct and Accelerated Systematic Coarse-graining Approaches. J. Chem. Phys. 2018, 148, 214105. (28) Shell, M. S. The Relative Entropy is Fundamental to Multiscale and Inverse Thermodynamic Problems. J. Chem. Phys. 2008, 129, 144108.

57



(29) Reith, D.; Pütz, M.; Müller-Plathe, F. Deriving Effective Mesoscale Potentials from Atomistic Simulations. J. Comp. Chem. 2003, 24, 1624–1636. (30) Siu, S. W. I.; Pluhackova, K.; Böckmann, R. A. Optimization of the OPLS-AA Force Field for Long Hydrocarbons. J. Chem. Theory Comput. 2012, 8, 1459–1470. (31) Sun, Q.; Faller, R. Systematic Coarse-graining of a Polymer Blend: Polyisoprene and Polystyrene. J. Chem. Theory Comput. 2006, 2, 607–615. (32) Widom, B. Intermolecular Forces and the Nature of the Liquid State. Science 1967, 157, 375–382. (33) Widom, B. Potential-distribution Theory and the Statistical Mechanics of Fluids. J. Phys. Chem. 1982, 86, 869–872. (34) Clark, A. J.; Guenza, M. G. Mapping of Polymer Melts onto Liquids of Soft-colloidal Chains. J. Chem. Phys. 2010, 132, 044902. (35) Sakai, H.; Stillinger, F. H.; Torquato, S. Equi-g(r) Sequence of Systems Derived from the Square-well Potential. J. Chem. Phys. 2002, 117, 297–307. (36) Dinpajooh, M.; Guenza, M. Thermodynamic Consistency in the Structure-based Integral Equation Coarse-grained Method. Polymer 2017, 117, 282–286. (37) Zhang, L.; Han, J.; Wang, H.; Car, R.; E, W. Deep Potential Molecular Dynamics: A Scalable Model with the Accuracy of Quantum Mechanics. Phys. Rev. Lett. 2018, 120, 143001. (38) Zhang, L.; Han, J.; Wang, H.; Car, R.; E, W. DeePCG: Constructing coarse-grained models via deep neural networks. J. Chem. Phys. 2018, 149, 034101. (39) Lyubimov, I.; Guenza, M. G. First-principle Approach to Rescale the Dynamics of Simulated Coarse-grained Macromolecular Liquids. Phys. Rev. E 2011, 84, 031801. 58


Page 58 of 68

Page 59 of 68 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) Lyubimov, I. Y.; Guenza, M. G. Theoretical Reconstruction of Realistic Dynamics of Highly Coarse-grained cis -1,4-polybutadiene Melts. J. Chem. Phys. 2013, 138, 12A546. (41) Di Pasquale, N.; Hudson, T.; Icardi, M. Systematic Derivation of Hybrid Coarsegrained Models. 2018, preprint. (42) Yoshimoto, Y.; Li, Z.; Kinefuchi, I.; Karniadakis, G. E. Construction of NonMarkovian Coarse-grained Models Employing the Moriâ

Zwanzig Formalism and

Iterative Boltzmann Inversion. J. Chem. Phys. 2017, 147, 244110. (43) Zwanzig, R. Nonequilibrium Statistical Mechanics; Oxford University Press: New York, 2001. (44) Dinpajooh, M.; Guenza, M. G. Soft Matter doi:10.1039/C8SM00868J. (45) Schweizer, K.; Curro, J. PRISM Theory of the Structure, Thermodynamics, and Phase Transitions of Polymer Liquids and Alloys. Adv. Polym. Sci. 1994, 116, 319 – 377. (46) Chandler, D.; Andersen, H. C. Optimized Cluster Expansions for Classical Fluids. II. Theory of Molecular Liquids. J. Chem. Phys. 1972, 57, 1930–1937. (47) Praprotnik, M.; Site, L. D.; Kremer, K. Multiscale Simulation of Soft Matter: From Scale Bridging to Adaptive Resolution. Annu. Rev. Phys. Chem. 2008, 59, 545–571. (48) Krekeler, C.; Agarwal, A.; Junghans, C.; Praprotnik, M.; Delle Site, L. Adaptive Resolution Molecular Dynamics Technique: Down to the Essential. J. Chem. Phys. 2018, 149 . (49) McCarty, J.; Lyubimov, I. Y.; Guenza, M. G. Effective Soft-Core Potentials and Mesoscopic Simulations of Binary Polymer Mixtures. Macromolecules 2010, 43, 3964–3979. (50) Flekkøy, E.; Coveney, P. From Molecular Dynamics to Dissipative Particle Dynamics. Phys. Rev. Lett. 1999, 83, 1775–1778. 59



(51) Flekkoy, E.; Coveney, P. V. P.; De Fabritiis G,; Flekko, E. G.; Fabritiis, G. D. Foundations of Dissipative Particle Dynamics. Physical review. E, Statistical physics, plasmas, fluids, and related interdisciplinary topics 2000, 62, 2140–57. (52) Espanol, P. Statistical Mechanics of Dissipative Particle Dynamics . EPL 1995, 30, 191–196. (53) Sambriski, E. J.; Yatsenko, G.; Nemirovskaya, M. A.; Guenza, M. G. Bridging Length Scales in Polymer Melt Relaxation for Macromolecules with Specific Local Structures. J. Phys. Condens. Matter 2007, 19, 205115. (54) Yatsenko, G.; Sambriski, E. J.; Nemirovskaya, M. A.; Guenza, M. Analytical SoftCore Potentials for Macromolecular Fluids and Mixtures. Phys. Rev. Lett. 2004, 93, 257803. (55) Yatsenko, G.; Sambriski, E. J.; Guenza, M. G. Coarse-grained Description of Polymer Blends as Interacting Soft-colloidal Particles. J. Chem. Phys. 2005, 122, 054907. (56) Clark, A. J.; McCarty, J.; Guenza, M. G. Effective Potentials for Representing Polymers in Melts as Chains of Interacting Soft Particles. J. Chem. Phys. 2013, 139, 124906. (57) Sambriski, E. J.; Yatsenko, G.; Nemirovskaya, M. A.; Guenza, M. G. Analytical Coarse-grained Description for Polymer Melts. J. Chem. Phys. 2006, 125, 234902. (58) McCarty, J.; Clark, A. J.; Lyubimov, I. Y.; Guenza, M. G. Thermodynamic Consistency between Analytic Integral Equation Theory and Coarse-Grained Molecular Dynamics Simulations of Homopolymer Melts. Macromolecules 2012, 45, 8482–8493. (59) Krakoviack, V.; Hansen, J. P.; Louis, A. A. Relating Monomer to Centre-of-Mass Distribution Functions in Polymer Solutions. EPL 2002, 58, 53–59. (60) Likos, C. N. Phys. Rep.; 2001; Vol. 348; pp 267–439. 60


Page 60 of 68

Page 61 of 68 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


(61) Hansen, J. P.; McDonald, I. R. Theory of Simple Liquids; Academic Press: Amsterdam, 2003. (62) Honnell, K. G.; Hall, C. K.; Dickman, R. On the Pressure Equation for Chain Molecules. J. Chem. Phys. 1987, 87, 664–674. (63) Stillinger, F. H.; Sakai, H.; Torquato, S. Statistical Mechanical Models with Effective Potentials: Definitions, Applications, and Thermodynamic Consequences. J. Chem. Phys. 2002, 117, 288–296. (64) Milano, G.; Müller-Plathe, F. Mapping Atomistic Simulations to Mesoscopic Models: A Systematic Coarse-Graining Procedure for Vinyl Polymer Chains. J. Phys. Chem. B 2005, 109, 18609–18619. (65) Ganguly, P.; Van Der Vegt, N. F. A. Representability and Transferability of KirkwoodBuff Iterative Boltzmann Inversion Models for Multicomponent Aqueous Systems. J. Chem. Theory Comput. 2013, 9, 5247–5256. (66) Rudzinski, J. F.; Noid, W. G. Coarse-graining Entropy, Forces, and Structures. J. Chem. Phys. 2011, 135, 214101. (67) Likos, C. N.; Lang, A.; Watzlawek, M.; Löwen, H. Criterion for Determining Clustering Versus Reentrant Melting Behavior for Bounded Interaction Potentials. Phys. Rev. E 2001, 63, 1–9. (68) Evans, R. The Nature of Liquid-Vapour Interface and Other Topics in the Statistical Mechanics of Non-Uniform, Classical Fluids. Adv. Phys. 1979, 28, 143–200. (69) Chatterjee, A. P.; Schweizer, K. S. Liquid-state theory of semidilute and concentrated polymer solutions. Macromolecules 1998, 31, 2353–2367. (70) Martin, M. G.; Siepmann, J. I. Transferable Potentials for Phase Equilibria. 1. UnitedAtom Description of n -Alkanes. J. Phys. Chem. B 1998, 102, 2569–2577. 61



(71) AvendanÌ o, C.; Lafitte, T.; Galindo, A.; Adjiman, C. S.; Jackson, G.; MuÌ ller, E. A. SAFT-“ Force Field for the Simulation of Molecular Fluids. 1. A Single-Site Coarse Grained Model of Carbon Dioxide. J. Phys. Chem. B 2011, 115, 11154–11169. (72) Orsi, M. Comparative Assessment of the ELBA Coarse-grained Model for Water. Mol. Phys. 2014, 112, 1566–1576. (73) Marrink, S. J.; de Vries, A. H.; Mark, A. E. Coarse Grained Model for Semiquantitative Lipid Simulations. J. Phys. Chem. B 2004, 108, 750–760. (74) Marrink, S. J.; Risselada, H. J.; Yefimov, S.; Tieleman, D. P.; De Vries, A. H. The MARTINI Force Field: Coarse-Grained Model for Biomolecular Simulations. J. Phys. Chem. B 2007, 111, 7812–7824. (75) Marrink, S. J.; Tieleman, D. P. Perspective on the Martini Model. Chem. Soc. Rev. 2013, 42, 6801–6822. (76) Das, A.; Andersen, H. C. The Multiscale Coarse-graining Method. IX. A General Method for Construction of Three Body Coarse-grained Force Fields. J. Chem. Phys. 2012, 136, 194114. (77) Zhang, Z.; Voth, G. A. Coarse-grained Representations of Large Biomolecular Complexes from Low-Resolution Structural Data. J. Chem. Theory Comput. 2010, 6, 2990– 3002. (78) Yang, D.; Wang, Q. Systematic and Simulation-Free Coarse Graining of Homopolymer Melts: a Relative-Entropy-Based Study. Soft Matter 2015, 11, 7109–7118. (79) Dunn, N. J. H.; Lebold, K. M.; DeLyser, M. R.; Rudzinski, J. F.; Noid, W. BOCS: Bottom-up Open-source Coarse-graining Software. J. Phys. Chem. B 2018, 122, 3363– 3377.

62


Page 62 of 68

Page 63 of 68 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


(80) Agrawal, V.; Arya, G.; Oswald, J. Simultaneous Iterative Boltzmann Inversion for Coarse-Graining of Polyurea. Macromolecules 2014, 47, 3378–3389. (81) Schaettle, K.; Ruiz Pestana, L.; Head-Gordon, T.; Lammers, L. N. A Structural Coarse-Grained Model for Clays Using Simple Iterative Boltzmann Inversion. J. Chem. Phys. 2018, 148, 222809. (82) Dunn, N. J. H.; Noid, W. G. Bottom-up Coarse-Grained Models that Accurately Describe the Structure, Pressure, and Compressibility of Molecular Liquids. J. Chem. Phys. 2015, 143, 243148. (83) Akkermans, R. L. C.; Briels, W. J. A Structure-Based Coarse-Grained Model for Polymer Melts. J. Chem. Phys. 2001, 114, 1020–1031. (84) Barker, J. A.; Henderson, D. Perturbation Theory and Equation of State for Fluids: The Square-Well Potential. J. Chem. Phys. 1967, 47, 2856–2861. (85) Weeks, J. D.; Chandler, D.; Andersen, H. C. Role of Repulsive Forces in Determining the Equilibrium Structure of Simple Liquids. J. Chem. Phys. 1971, 54, 5237–5247. (86) Zwicker, J.; Lovett, R. When Does a Pair Correlation Function Fix the State of an Equilibrium System? J. Chem. Phys. 1990, 93, 6752–6755. (87) Ascarelli, P.; Harrison, R. J. Density-Dependent Potentials and the Hard-Sphere Model for Liquid Metals. Phys. Rev. Lett. 1969, 22, 385–388. (88) Montes-Saralegui, M.; Kahl, G.; Nikoubashman, A. On the Applicability of Density Dependent Effective Interactions in Cluster-forming Systems. J. Chem. Phys. 2017, 146, 54904. (89) Song, J.; Hsu, D. D.; Shull, K. R.; Phelan, F. R.; Douglas, J. F.; Xia, W.; Keten, S. Energy Renormalization Method for the Coarse-Graining of Polymer Viscoelasticity. Macromolecules 2018, 51, 3818–3827. 63



(90) Xia, W.; Song, J.; Hansoge, N. K.; Phelan, F. R.; Keten, S.; Douglas, J. F. Energy Renormalization for Coarse-Graining the Dynamics of a Model Glass-Forming Liquid. J. Phys. Chem. B 2018, 122, 2040–2045. (91) Karplus, M.; Kushick, J. N. Method for Estimating the Configurational Entropy of Macromolecules. Macromolecules 1981, 14, 325–332. (92) Gao, M. C.; Widom, M. Information Entropy of Liquid Metals. 2017, V, 5–10. (93) Armstrong, J. A.; Chakravarty, C.; Ballone, P. Statistical Mechanics of Coarse Graining: Estimating Dynamical Speedups from Excess Entropies. J. Chem. Phys. 2012, 136 . (94) Zamponi, M.; Wischnewski, A.; Monkenbusch, M.; Willner, L.; Richter, D.; Falus, P.; Farago, B.; Guenza, M. G. Cooperative Dynamics in Homopolymer Melts: A Comparison of Theoretical Predictions with Neutron Spin Echo Experiments. J. Phys. Chem. B 2008, 112, 16220–16229. (95) Schweizer, K. S.; Curro, J. G. Integral-Equation Theory of the Structure of Polymer Melts. Phys. Rev. Lett. 1987, 58, 246–249. (96) Schweizer, K. S.; Curro, J. G. RISM Theory of Polymer Liquids: Analytical Results for Continuum Models of Melts and Alloys. Chem. Phys. 1990, 149, 105–127. (97) Schweizer, K. S.; Gurro, J. G. Integral Equation Theories of the Structure, Thermodynamics, and Phase Transitions of Polymer Fluids. Adv. Chem. Phys. 1997, 98, 1–142. (98) Dunn, N. J. H.; Foley, T. T.; Noid, W. G. Van der Waals Perspective on CoarseGraining: Progress toward Solving Representability and Transferability Problems. Acc. Chem. Res. 2016, 49, 2832–2840.

64


Page 64 of 68

Page 65 of 68 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


(99) Carnahan, N. F.; Starling, K. E. Thermodynamic Properties of a Rigid-Sphere Fluid. J. Chem. Phys. 1970, 53, 600–603. (100) Strobl, G. R. The Physics of Polymers: Concepts for Understanding their Structure and Behavior, 2nd ed.; Springer: Berlin Heidelberg, 1997. (101) Flory, P. J. Principles of Polymer Chemistry, sixteenth ed.; Cornell University Press: Ithaca, 1953. (102) Debye, P.; Bueche, F. Distribution of Segments in a Coiling Polymer Molecule. J. Chem. Phys 1952, 20, 1337–1338. (103) Hsu, H.-P.; Kremer, K. Static and Dynamic Properties of Large Polymer Melts in Equilibrium. J. Chem. Phys. 2016, 144, 154907. (104) Grayce, C. J.; Yethiraj, A.; Schweizer, K. S. Liquid-state Theory of the Density Dependent Conformation of Nonpolar Linear Polymers. J. Chem. Phys. 1994, 100, 6857– 6872. (105) Stanley, H. E. Introduction to Phase Transitions and Critical Phenomena; Oxford University Press: New York, 1987. (106) Biben, T.; Hansen, J. P. Phase Separation of Asymmetric Binary Hard-Sphere Fluids. Phys. Rev. Lett. 1991, 66, 2215–2218. (107) Ozog, D.; McCarty, J.; Gossett, G.; Malony, A. D.; Guenza, M. G. Fast Equilibration of Coarse-grained Polymeric Liquids. J. Comput. Sci. 2015, (108) Zhang, G.; Moreira, L. A.; Stuehn, T.; Daoulas, K. C.; Kremer, K. Equilibration of High Molecular Weight Polymer Melts: A Hierarchical Strategy. ACS Macro Lett. 2014, 3, 198–203. (109) Laso, M.; Öttinger, H. C.; Suter, U. W. Bond-length and Bond-angle Distributions in Coarse-grained Polymer Chains. J. Chem. Phys. 1991, 95, 2178–2182. 65



(110) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comp. Phys. 1995, 117, 1–19. (111) Towns, J.; Cockerill, T.; Dahan, M.; Foster, I.; Gaither, K.; Grimshaw, A.; Hazlewood, V.; Lathrop, S.; Lifka, D.; Peterson, G. D.; Roskies, R.; Scott, J. R.; WilkinsDiehr, N. XSEDE: Accelerating Scientific Discovery. Comput. Sci. Eng. 2014, 16, 62–74.

66


Page 66 of 68

Page 67 of 68 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


Marina G. Guenza is a Professor of Physical Chemistry at the University of Oregon: her research involves the development of theoretical and computational approaches to predict the properties of complex macromolecular systems in their liquid state from their atomic-level structure. She received her Ph.D. in Physical Chemistry from the consortium of the Universities of Genoa, Pavia, and Turin, Italy in 1989, and she joined the Italian National Laboratory (CNR) in 1989. In 1994, as a visiting scientist in the Chemistry Department at the University of Chicago, she developed models for the dynamics of proteins in collaboration with Prof Karl Freed. From 1995 to 1997 she was a visiting scientist in the Material Science and Engineering Department at the University of Illinois at Urbana Champaign, working on integral equation theory and mode-mode coupling theory with Prof. Kenneth Schweizer. In 2002, Marina joined the Chemistry and Biochemistry Department at the University of Oregon, where she remains today. In 2011 she was elected Fellow of the American Physical Society in the Division of Polymer Physics. Mohammadhasan (Hadi) Dinpajooh obtained his Ph.D. under the supervision of Prof. D. V. Matyushov from Arizona State University in 2016, where he investigated the solvent electrostatic response of simple solutes and proteins with a focus on the role of polarizability. Subsequently, he joined the group of Marina Guenza in December 2016 as a postdoctoral researcher, where he is now involved in projects ranging from developing the Integral Equation theory of Coarse-Graining to modeling protein dynamics. James McCarty received his B.S. with honors in biochemistry from the California Polytechnic State University in San Luis Obispo, CA. He subsequently moved to Eugene, OR to pursue doctoral work in theoretical physical chemistry at the University of Oregon in the group of Prof. Marina Guenza. After receiving his Ph.D. in chemistry in 2013, he joined the group of Prof. Michele Parrinello as a postdoctoral associate at ETH Zürich, Switzerland. Since 2017, he has been working as a postdoctoral scholar in the department of Chemistry and Biochemistry at the University of California, Santa Barbara in the group of Prof. Joan-Emma Shea. His research interests include molecular dynamics simulation, enhanced sampling of rare events, and polymer physics. Ivan Lyubimov received his MS degree in Physics under the supervision of Prof. Nail F. Fatkullin from the Kazan Federal University, Russia. He subsequently moved to the University of Oregon, where he obtained his Ph.D. in Physical Chemistry under the supervision of Prof. Marina G. Guenza. Ivan did postdoctoral research at the University of Chicago, studying stable glasses by mean of computer simulations in Prof. Juan J. de Pablo group. Since 2016 he is a postdoctoral researcher at the University of Delaware in Prof. Arthi Jayaraman group. His current main research interests include development of approaches that combine liquid state theory and computer simulations to investigate polymeric materials.




Page 68 of 68

Accuracy, Transferability, and Efficiency of Coarse-Grained Models of

Recommend Documents