Mesoscale model of polymer melt structure: self-consistent mapping of molecular correlations to coarse-grained potentials

Development and application of coarse-graining methods to condensed phases of macromolecules is an active area of research. Multiscale modeling of polymeric systems using coarse-graining methods presents unique challenges. Here we apply a coarse-graining method that self-consistently maps structural correlations from detailed molecular dynamics (MD) simulations of alkane oligomers onto coarse-grained potentials using a combination of MD and inverse Monte Carlo methods. Once derived, the coarse-grained potentials allow computationally efficient sampling of ensemble of conformations of significantly longer polyethylene chains. Conformational properties derived from coarse-grained simulations are in excellent agreement with experiments. The level of coarse graining provides a control over the balance of computational efficiency and retention of chemical identity of the underlying polymeric system. Challenges to extension and application of this and similar structure-based coarse-graining methods to model dynamics and phase behavior in polymeric systems are briefly discussed.     

The dynamo library for molecular simulations using hybrid quantum mechanical and molecular mechanical potentials

The Dynamo module library has been developed for the simulation of molecular systems using hybrid quantum mechanical (QM) and molecular mechanical (MM) potentials. Dynamo is not a program package but is a library of Fortran 90 modules that can be employed by those interested in writing their own programs for performing molecular simulations. The library supports a range of different types of molecular calculation including geometry optimizations, reaction‐path determinations and molecular dynamics and Monte Carlo simulations. This article outlines the general structure and capabilities of the library and describes in detail Dynamo's semiempirical QM/MM hybrid potential. Results are presented to indicate three particular aspects of this implementation—the handling of long‐range nonbonding interactions, the nature of the boundary between the quantum mechanical and molecular mechanical atoms and how to perform path‐integral hybrid‐potential molecular dynamics simulations. © 2000 John Wiley & Sons, Inc. J Comput Chem 21: 1088–1100, 2000     

Rigid quantum Monte Carlo simulations of condensed molecular matter: water clusters in the n=2-->8 range

The numerical advantage of quantum Monte Carlo simulations of rigid bodies relative to the flexible simulations is investigated for some simple systems. The results show that if high frequency modes in molecular condensed matter are predominantly in the ground state, the convergence of path integral simulations becomes nonuniform. Rigid body quantum parallel tempering simulations are necessary to accurately capture thermodynamic phenomena in the temperature range where the dynamics are influenced by intermolecular degrees of freedom; the stereographic projection path integral adapted for quantum simulations of asymmetric tops is a significantly more efficient strategy compared with Cartesian coordinate simulations for molecular condensed matter under these conditions. The reweighted random series approach for stereographic path integral Monte Carlo is refined and implemented for the quantum simulation of water clusters treated as an assembly of rigid asymmetric tops.     

Monte Carlo simulations of polymer dynamics: Recent advances

A brief review is given of applications of Monte Carlo simulations to study the dynamical properties of coarse-grained models of polymer melts, emphasizing the crossover from the Rouse model toward reptation, and the glass transition. The extent to which Monte Carlo algorithms can mimic the actual chain dynamics is critically examined, and the need for the use of coarse-grained rather than fully atomistic models for such simulations is explained. It is shown that various lattice and continuum models yield qualitatively similar results, and the behavior agrees with the findings of corresponding molecular dynamics simulations and experiments, where available. It is argued that these simulations significantly enhance our understanding of the theoretical concepts on the dynamics of dense macromolecular systems. © 1997 John Wiley & Sons, Inc.     

New phase for one-component hard spheres

A completely new phase for one-component hard spheres is reported in an unexpected region of the phase diagram. The new phase is observed at compressibility factors intermediate between the solid and the metastable branches. It can be obtained from either Monte Carlo simulations alone or a combination of Monte Carlo and molecular dynamics calculations. An analysis of the intermediate scattering function data shows that the new phase is in a stable equilibrium. Radial distribution function data, configurational snapshots, bond order parameters, and translational order parameters obtained from molecular simulations indicate that the new phase is significantly different from the isotropic liquid, metastable, or crystalline phases traditionally observed in hard sphere systems. This result significantly changes our previous understanding of the behavior of hard spheres.     

Application of combined Monte Carlo and molecular dynamics method to simulation of dipalmitoyl phosphatidylcholine lipid bilayer

We describe a new equilibration procedure for the atomic level simulation of a hydrated lipid bilayer. The procedure consists of alternating molecular dynamics trajectory calculations in a constant surface tension and temperature ensemble with configurational bias Monte Carlo moves to different regions of the configuration space of the bilayer, in a constant volume and temperature ensemble. The procedure is described in detail and is applied to a bilayer of 100 molecules of dipalmitoyl phosphatidylcholine (DPPC) and 3205 water molecules. We find that the hybrid simulation procedure enhances the equilibration of the bilayer as measured by the convergence of the area per molecule and the segmental order parameters, as compared with a simulation using only molecular dynamics (MD). Progress toward equilibration is almost three times as fast in CPU time, compared with a purely MD simulation. Equilibration is complete, as judged by the lack of energy drift in three separate 200-ps runs of continuous MD started from different initial states. Results of the simulation are presented and compared with experimental data and with other recent simulations of DPPC. ©1999 John Wiley & Sons, Inc. J Comput Chem 20: 1153–1164, 1999     

Molecular simulations of aqueous electrolyte solubility: 1. The expanded-ensemble osmotic molecular dynamics method for the solution phase

We have developed a molecular-level simulation technique called the expanded-ensemble osmotic molecular dynamics (EEOMD) method, for studying electrolyte solution systems. The EEOMD method performs simulations at a fixed number of solvent molecules, pressure, temperature, and overall electrolyte chemical potential. The method combines elements of constant pressure-constant temperature molecular dynamics and expanded-ensemble grand canonical Monte Carlo. The simulated electrolyte solution systems contain, in addition to solvent molecules, full and fractional ions and undissociated electrolyte molecular units. The fractional particles are coupled to the system via a coupling parameter that varies between 0 (no interaction between the fractional particle and the other particles in the system) and 1 (full interaction between the fractional particle and the other particles in the system). The time evolution of the system is governed by the constant pressure-constant temperature equations of motion and accompanied by random changes in the coupling parameter. The coupling-parameter changes are accepted with a probability derived from the expanded-ensemble osmotic partition function corresponding to the prescribed electrolyte chemical potential. The coupling-parameter changes mimic insertion/deletion of particles as in a crude grand canonical Monte Carlo simulation; if the coupling parameter becomes 0, the fractional particles disappear from the system, and as the coupling parameter reaches unity, the fractional particles become full particles. The method is demonstrated for a model of NaCl in water at ambient conditions. To test our approach, we first determine the chemical potential of NaCl in water by the thermodynamic integration technique and by the expanded-ensemble method. Then, we carry out EEOMD simulations for different specified values of the overall NaCl chemical potential and measure the concentration of ions resulting from the simulations. Both computations give consistent results, validating the EEOMD methodology.     

Coarse-grained molecular dynamics modeling of strongly associating fluids: thermodynamics, liquid structure, and dynamics of symmetric binary mixture fluids

Symmetric binary mixtures capable of strong association via a highly directional and saturable specific interaction between unlike molecules are investigated by canonical molecular dynamics simulations. The specific interaction of the molecules is defined in a new coarse-grained pair potential that can be applied in continuous molecular dynamics as well as in Monte Carlo simulations. The thermodynamic, structural, and dynamic properties of the associating mixture fluids are investigated as a function of density, temperature, and association strength of the specific interaction. Detailed analysis of the simulation data confirms a two-stage mechanism in the formation of specific bonds with increasing interaction strength, including a fast dimerization process and a subsequent stage of perfecting the bonds. A large heat capacity peak is found during the formation or breaking of the bonds, reflecting the large energy fluctuation introduced by the strong association. The fractions of nonbonded molecules obtained from the simulations as a function of density, temperature, and interaction strength are in excellent agreement with the predictions of Wertheim's thermodynamic perturbation theory. The translational and rotational dynamics of the Tmer mixture are effectively retarded with increasing association strength and are analyzed in terms of autocorrelation functions and a non-Gaussian parameter for the translational dynamics. The lifetimes of molecules in bonded and nonbonded states provide detailed information about the transformation of molecules between the bonded state and the nonbonded state. Finally, simulation sampling problems inherent to strongly interacting systems are easily overcome using the parallel tempering simulation technique. This latter result confirms that with the new continuous coarse-grained simulation potential we have a versatile and flexible interaction potential that can be used with many available molecular dynamics and Monte Carlo algorithms under various ensembles.     

A rapidly convergent simulation method: Mixed Monte Carlo/stochastic dynamics

Although Monte Carlo and molecular dynamics are the primary methods used for free energy simulations of molecular systems, their application to molecules that have multiple conformations separated by energy barriers of ≥ 3 kcal/mol is problematic because of slow rates of convergence. In this article we introduce a hybrid simulation method termed MC-SD which mixes Monte Carlo (MC) and stochastic dynamics (SD). This new method generates a canonical ensemble via alternating MC and SD steps and combines the local exploration strengths of dynamics with the barrier-crossing ability of large-step Monte Carlo. Using calculations on double-well potentials and long simulations (10⁸ steps of MC and 1 μs of SD) of the simple, conformationally flexible molecule n-pentane, we find that MC-SD simulations converage faster than either MC or SD alone and generate ensembles which are equivalent to those created by classical MC or SD. Using pure SD at 300 K, the conformational populations of n-pentane are shown to be poorly converged even after a full microsecond of simulation. © 1994 by John Wiley & Sons, Inc.     

MCDP: an advanced tool to simulate comb‐like polymers

A strategy to study polymeric systems with ordered structures, and in particular comb‐like polymers, is presented. These are dense systems for which atomistic simulations with conventional methods are difficult or even impracticable. The strategy, which has been incorporated into a computer program named MCDP, is based on a Configuration Bias Monte Carlo algorithm and a method to investigate the structure of crystalline polymers using force‐field calculations. To obtain a maximum efficiency, the MCDP computer program has been optimized and parallelized. The ability of MCDP to investigate ordered polymers have been tested by simulating two complex systems: (1) the crystal structure of poly(4‐methyl‐1‐pentene), and (2) the biphasic structure of poly(α‐octyl‐β‐L‐aspartate), a comb‐like polyamide derived from β‐amino acids. The results obtained from MCDP simulations demonstrates the efficiency and reliability of this method to study both the NVT and NPT behavior of ordered dense polymers. © 2000 John Wiley & Sons, Inc. J Comput Chem 22: 162–171, 2001     

Microsecond simulations of DNA and ion transport in nanopores with novel ion–ion and ion–nucleotides effective potentials

We developed a novel scheme based on the grand‐canonical Monte Carlo/Brownian dynamics simulations and have extended it to studies of ion currents across three nanopores with the potential for single‐stranded DNA (ssDNA) sequencing: solid‐state nanopore Si₃N₄, α‐hemolysin, and E111N/M113Y/K147N mutant. To describe nucleotide‐specific ion dynamics compatible with ssDNA coarse‐grained model, we used the inverse Monte Carlo protocol, which maps the relevant ion–nucleotide distribution functions from all‐atom molecular dynamics (MD) simulations. Combined with the previously developed simulation platform for Brownian dynamics simulations of ion transport, it allows for microsecond‐ and millisecond‐long simulations of ssDNA dynamics in the nanopore with a conductance computation accuracy that equals or exceeds that of all‐atom MD simulations. In spite of the simplifications, the protocol produces results that agree with the results of previous studies on ion conductance across open channels and provide direct correlations with experimentally measured blockade currents and ion conductances that have been estimated from all‐atom MD simulations. © 2014 Wiley Periodicals, Inc.     

Hybrid Monte Carlo: An efficient algorithm for condensed matter simulation

A detailed comparison has been made of the performance of molecular dynamics and hybrid Monte Carlo simulation algorithms for calculating thermodynamic properties of 2D Lennard-Jonesium. The hybrid Monte Carlo simulation required an order of magnitude fewer steps than the molecular dynamics simulation to calculate reproducible values of the specific heat. The ergodicity of the two algorithms was compared via the use of intermediate scattering functions. For classical systems the intermediate scattering functions should be real; however, a simple analysis demonstrates that this function will have a significant imaginary component when ergodicity breaks down. For q vectors near the zone boundary, the scattering functions are real for both algorithms. However, for q vectors near the zone center (i.e., harmonic, weakly coupled modes), the scattering function calculated via molecular dynamics had a significantly larger imaginary component than that calculated using hybrid Monte Carlo. Therefore, the hybrid Monte Carlo algorithm is more ergodic and samples phase space more efficiently than molecular dynamics for simulations of 2D Lennard-Jonesium. © 1994 by John Wiley & Sons, Inc.     

Calculation of pressure using the virtual-volume-variation method and the virial method from chain conformations obtained by Monte Carlo simulations on the second nearest neighbor diamond lattice

For a model system of polyethylene of chain lengths 40 and 100 carbon atoms, we calculated the pressure at different densities and compared them with the experimental values. The simulation was conducted on the second nearest neighbor diamond lattice, and the pressure was calculated using the virtual-volume-variation method after the system was reverse mapped to its fully atomistic form in continuous space and energy minimized. In addition, the pressure was also calculated from the virial route by conducting a short molecular dynamics simulation starting from the energy minimized structure. We show that the pressure obtained from our simulations is quite reasonable in the length of simulation time (in Monte Carlo steps) normally employed in our group. These results provide additional evidence for the equilibration of our model systems, and methodology to calculate the pressure in our lattice models.     

A new form for intermolecular potential energy functions

A new intermolecular potential energy function is presented for Monte Carlo and molecular dynamics simulations of liquids, solutions and other molecular assemblies. The potential energy function is expressed in terms of intermolecular overlap integrals over localized molecular orbitais of isolated molecules and Coulomb potentials between fractional point charges placed on the nuclei. The potential function is easy to generate and is applicable to a wide range of molecules. As examples the potential functions are generated for the water, ammonia and hydrogen fluoride dimers.     

Hybrid Monte Carlo implementation of the Fourier path integral algorithm

This paper formulates a hybrid Monte Carlo implementation of the Fourier path integral (FPI-HMC) approach with partial averaging. Such a hybrid Monte Carlo approach allows one to generate collective moves through configuration space using molecular dynamics while retaining the computational advantages associated with the Fourier path integral Monte Carlo method. In comparison with the earlier Metropolis Monte Carlo implementations of the FPI algorithm, the present HMC method is shown to be significantly more efficient for quantum Lennard-Jones solids and suggests that such algorithms may prove useful for efficient simulations of a range of atomic and molecular systems.     

Molecular silverware. I. General solutions to excluded volume constrained problems

General mathematical solutions to excluded volume constrained problems in computational chemistry are reported. The solutions have been used to create a new family of molecular modeling algorithms to facilitate the study of molecular interactions in condensed phases. The new algorithms, collectively known as Molecular Silverware, are for the most part interactive and designed for packing, solvating, and sampling molecules embedded in simple or complex topological environments. Multifolded, disconnected, or porous molecular structures are permitted. Molecular Silverware assists the preparation of Monte Carlo and molecular dynamics simulations at a small fraction of the total simulation time. Primary targets for applications include the study of molecular recognition mechanisms and the selective binding of DNA, RNA, peptides, saccharides and other biopolymers in solution as well as the prediction of phase separation behavior and physical properties of non-crystalline condensed phases such as bulk polymers, polymer blends, organic liquids, membranes, micelles, gels, crosslinked networks, glasses, and amorphous heterogeneous catalysts. As a result of this new approach to excluded volume constraints, the computer simulation of noncrystalline condensed phases is no longer hampered by the lack of a general and efficient method for the creation and configurational sampling of small and large molecular assemblies at high densities.     

Using internal and collective variables in Monte Carlo simulations of nucleic acid structures: chain breakage/closure algorithm and associated Jacobians

This article describes a method for solving the geometric closure problem for simplified models of nucleic acid structures by using the constant bond lengths approximation. The resulting chain breakage/closure equations, formulated in the space of variable torsion and bond angles, are easy to solve, and have only two solutions. The analytical simplicity is in contrast with the high complexity of the closure problem in the torsion angle space with at most 16 solutions, which has been dealt with by several authors and was solved analytically by Wu and Deem (J. Chem. Phys. 1999, 111, 6625). The discussion on the choice of variables and associated Jacobians is focussed on the question of how conformational equilibration is affected in Monte Carlo simulations of molecular systems. In addition to the closure of the phosphate backbone, it is necessary to also solve the closure problem for the five-membered flexible furanose sugar ring. Explicit closure equations and the resulting Jacobians are given both for the complete four-variable model of the furanose ring and simulations in the phase-amplitude space of the five-membered ring, which are based on the approximate two-variable model of furanose introduced by Gabb et al. (J. Comput. Chem. 1995, 16, 667). The suggested closure algorithm can be combined with collective variables defined by translations and rotations of the monomeric nucleotide units. In comparison with simple internal coordinate moves, the resulting concerted moves describe local structural changes that have high acceptance rates and enable fast conformational equilibration. Appropriate molecular models are put forward for prospective Monte Carlo simulations of nucleic acids, but can be easily adapted to other biomolecular systems, such as proteins and lipid structures in biological membranes.     

MOLSIM: A modular molecular simulation software

The modular software MOLSIM for all‐atom molecular and coarse‐grained simulations is presented with focus on the underlying concepts used. The software possesses four unique features: (1) it is an integrated software for molecular dynamic, Monte Carlo, and Brownian dynamics simulations; (2) simulated objects are constructed in a hierarchical fashion representing atoms, rigid molecules and colloids, flexible chains, hierarchical polymers, and cross‐linked networks; (3) long‐range interactions involving charges, dipoles and/or anisotropic dipole polarizabilities are handled either with the standard Ewald sum, the smooth particle mesh Ewald sum, or the reaction‐field technique; (4) statistical uncertainties are provided for all calculated observables. In addition, MOLSIM supports various statistical ensembles, and several types of simulation cells and boundary conditions are available. Intermolecular interactions comprise tabulated pairwise potentials for speed and uniformity and many‐body interactions involve anisotropic polarizabilities. Intramolecular interactions include bond, angle, and crosslink potentials. A very large set of analyses of static and dynamic properties is provided. The capability of MOLSIM can be extended by user‐providing routines controlling, for example, start conditions, intermolecular potentials, and analyses. An extensive set of case studies in the field of soft matter is presented covering colloids, polymers, and crosslinked networks. © 2015 The Authors. Journal of Computational Chemistry Published by Wiley Periodicals, Inc.     

Quantum effects on adsorption and diffusion of hydrogen and deuterium in microporous materials

Monte Carlo and molecular dynamics simulations and neutron scattering experiments are used to study the adsorption and diffusion of hydrogen and deuterium in zeolite Rho in the temperature range of 30-150 K. In the molecular simulations, quantum effects are incorporated via the Feynman-Hibbs variational approach. We suggest a new set of potential parameters for hydrogen, which can be used when Feynman-Hibbs variational approach is used for quantum corrections. The dynamic properties obtained from molecular dynamics simulations are in excellent agreement with the experimental results and show significant quantum effects on the transport at very low temperature. The molecular dynamics simulation results show that the quantum effect is very sensitive to pore dimensions and under suitable conditions can lead to a reverse kinetic molecular sieving with deuterium diffusing faster than hydrogen.     

Free energies and phase equilibria of chain molecules

In this Article, a review is presented of recent developments in Monte Carlo simulations of chain molecules. The Rosenbluth chain insertion technique is used to calculate the free energy of the chain molecules. Furthermore, this insertion method is used to generate biased Monte Carlo moves. It is shown that this bias can be removed by adjusting the acceptance rules such that configurations are generated with their correct Boltzmann weight. This configurational-bias Monte Carlo method can be combined with the Gibbs-ensemble technique which results in an efficient method to simulate phase equilibria of chain molecules.