Isomolar-semigrand ensemble molecular dynamics: application to vapor-liquid equilibrium of the mixture methane/ethane

The isomolar-semigrand ensemble molecular dynamics (iSGMD) method is applied to the simulation of the binary system methane/ethane. The vapor-liquid equilibrium properties of this system at a temperature of 192.37 K are computed using the Gibbs-Duhem integration method. The iSGMD method, which resembles conventional hybrid Monte Carlo (MC) but is applicable to phase equilibrium calculations, is designed to overcome the difficulties associated with performing standard Monte Carlo-type particle transformations in liquid systems that are very dense and/or are comprised of complex molecules with many intramolecular degrees of freedom. This work shows that particle transformations using the iSGMD method for the simple system methane/ethane are at least 25 times more successful than standard MC-type transformations. The P-x-y curve for the system methane/ethane at 192.37 K computed using iSGMD simulations agrees very well with the experimental P-x-y curve as well as results of a previous MC study.     

Multibaric-multithermal ensemble molecular dynamics simulations

We present new generalized-ensemble molecular dynamics simulation algorithms, which we refer to as the multibaric-multithermal molecular dynamics. We describe three algorithms based on (1) the Nosé thermostat and the Andersen barostat, (2) the Nosé-Poincaré thermostat and the Andersen barostat, and (3) the Gaussian thermostat and the Andersen barostat. The multibaric-multithermal simulations perform random walks widely both in the potential-energy space and in the volume space. Therefore, one can calculate isobaric-isothermal ensemble averages in wide ranges of temperature and pressure from only one simulation run. We test the effectiveness of the multibaric-multithermal algorithm by applying it to a Lennard-Jones 12-6 potential system.     

A new method of semigrand canonical ensemble to calculate first-order phase transitions for binary mixtures

First-order phase transitions of binary mixtures at the given pressure (P) and temperature (T) are studied by taking into account the composition fluctuations. Isothermal-isobaric semigrand canonical ensemble is adopted to find the relations among the total number of molecules, the composition fluctuations and Gibbs free energy density. By combining two identical subsystems of mixtures successively, the free energy density is transformed until being stable and its linear segments represent phase transitions. A new method is developed to calculate the phase equilibriums of binary mixtures. The method handles multiple types and number of phase equilibriums at single time and its solutions are physically justified. One example is shown for calculating the phase diagram of binary Lennard-Jones mixture. It demonstrates that the fluctuations of the total number of molecules in mixtures are fundamental behind phase transitions and the van der Waals loops in Gibbs free energy are reasonable.     

Ab initio centroid path integral molecular dynamics: application to vibrational dynamics of diatomic molecular systems

An ab initio centroid molecular dynamics (CMD) method is developed by combining the CMD method with the ab initio molecular orbital method. The ab initio CMD method is applied to vibrational dynamics of diatomic molecules, H2 and HF. For the H2 molecule, the temperature dependence of the peak frequency of the vibrational spectral density is investigated. The results are compared with those obtained by the ab initio classical molecular dynamics method and exact quantum mechanical treatment. It is shown that the vibrational frequency obtained from the ab initio CMD approaches the exact first excitation frequency as the temperature lowers. For the HF molecule, the position autocorrelation function is also analyzed in detail. The present CMD method is shown to well reproduce the exact quantum result for the information on the vibrational properties of the system.     

Equilibrium distributions of dipalmitoyl phosphatidylcholine and dilauroyl phosphatidylcholine in a mixed lipid bilayer: atomistic semigrand canonical ensemble simulations

Conventional molecular dynamics (MD) simulations are seriously limited by the slow rate of diffusive mixing in their ability to predict lateral distributions of different lipid types within mixed-lipid bilayers using atomistic models. A method to overcome this limitation, using configuration-bias Monte Carlo (MC) "mutation" moves to transform lipids from one type to another in dynamic equilibrium, is demonstrated in binary fluid-phase mixtures of lipids whose tails differ in length by four carbons. The hybrid MC-MD method operates within a semigrand canonical ensemble, so that an equilibrium composition of the mixture is determined by a constant difference in chemical potential (Delta(mu)) chosen for the components. Within several nanoseconds, bilayer structures initiated as pure dipalmitoyl phosphatidylcholine (DPPC) or pure dilauroyl phosphatidylcholine (DLPC) converge to a common composition and structure in independent simulations conducted at the same Delta(mu). Trends in bilayer thickness, area per lipid, density distributions across the bilayer, and order parameters have been investigated at three mixture compositions and compared with results from the pure bilayers at 323 K. The mixtures showed a moderate increase in DPPC acyl tail sites crossing the bilayer midplane relative to pure DPPC. Correlations between lateral positions of the two lipid types within or across the bilayer were found to be weak or absent. While the lateral distribution is consistent with nearly ideal mixing, the dependence of composition on Delta(mu) indicates a positive excess free energy of mixing.     

Statistical temperature molecular dynamics: application to coarse-grained beta-barrel-forming protein models

Recently the authors proposed a novel sampling algorithm, "statistical temperature molecular dynamics" (STMD) [J. Kim et al., Phys. Rev. Lett. 97, 050601 (2006)], which combines ingredients of multicanonical molecular dynamics and Wang-Landau sampling. Exploiting the relation between the statistical temperature and the density of states, STMD generates a flat energy distribution and efficient sampling with a dynamic update of the statistical temperature, transforming an initial constant estimate to the true statistical temperature T(U), with U being the potential energy. Here, the performance of STMD is examined in the Lennard-Jones fluid with diverse simulation conditions, and in the coarse-grained, off-lattice BLN 46-mer and 69-mer protein models, exhibiting rugged potential energy landscapes with a high degree of frustration. STMD simulations combined with inherent structure (IS) analysis allow an accurate determination of protein thermodynamics down to very low temperatures, overcoming quasiergodicity, and illuminate the transitions occurring in folding in terms of the energy landscape. It is found that a thermodynamic signature of folding is significantly suppressed by accurate sampling, due to an incoherent contribution from low-lying non-native IS in multifunneled landscapes. It is also shown that preferred accessibility to such IS during the collapse transition is intimately related to misfolding or poor foldability.     

Thermodynamic properties in the molecular dynamics ensemble applied to the gaussian core model fluid

The thermodynamic properties of pressure, energy, isothermal pressure coefficient, thermal expansion coefficient, isothermal and adiabatic compressibilities, isobaric and isochoric heat capacities, Joule-Thomson coefficient, and speed of sound are considered in a classical molecular dynamics ensemble. These properties were obtained using the treatment of Lustig [J. Chem. Phys. 100, 3048 (1994)] and Meier and Kabelac [J. Chem. Phys. 124, 064104 (2006)], whereby thermodynamic state variables are expressible in terms of phase-space functions determined directly from molecular dynamics simulations. The complete thermodynamic information about an equilibrium system can be obtained from this general formalism. We apply this method to the gaussian core model fluid because the complex phase behavior of this simple model provides a severe test for this treatment. Waterlike and other anomalies are observed for some of the thermodynamic properties of the gaussian core model fluid.     

Constrained dynamics and extraction of normal modes from ab initio molecular dynamics: application to ammonia

We present a methodology for extracting phonon data from ab initio Born-Oppenheimer molecular dynamics calculations of molecular crystals. Conventional ab initio phonon methods based on perturbations are difficult to apply to lattice modes because the perturbation energy is dominated by intramolecular modes. We use constrained molecular dynamics to eliminate the effect of bond bends and stretches and then show how trajectories can be used to isolate and define in particular, the eigenvalues and eigenvectors of modes irrespective of their symmetry or wave vector. This is done by k-point and frequency filtering and projection onto plane wave states. The method is applied to crystalline ammonia: the constrained molecular dynamics allows a significant speed-up without affecting structural or vibrational modes. All Gamma point lattice modes are isolated: the frequencies are in agreement with previous studies; however, the mode assignments are different.     

Determination of solvation free energies by adaptive expanded ensemble molecular dynamics

A new method of calculating absolute free energies is presented. It was developed as an extension to the expanded ensemble molecular dynamics scheme and uses probability density estimation to continuously optimize the expanded ensemble parameters. The new method is much faster as it removes the time-consuming and expertise-requiring step of determining balancing factors. Its efficiency and accuracy are demonstrated for the dissolution of three qualitatively very different chemical species in water: methane, ionic salts, and benzylamine. A recently suggested optimization scheme by Wang and Landau [Phys. Rev. Lett. 86, 2050 (2001)] was also implemented and found to be computationally less efficient than the proposed adaptive expanded ensemble method.     

Quaternary liquid-liquid equilibria: Experimental determination and correlation of equilibrium data. Part III. New methods of representation and correlation of liquid-liquid equilibria in quaternary systems

Part III of this work describes two new methods of representation and correlation of liquid-liquid equilibria in quaternary systems. These methods are appropiate for the interpolation of new nonexperimental equilibrium data. The methods suggested are applied to the quaternary systems presented in Parts I and II of this series (Ruiz and Prats, 1982a,b).     

Self-consistent mean-field model based on molecular dynamics: application to lipid-cholesterol bilayers

We have developed a dynamic self-consistent mean-field model, based on molecular-dynamics simulations, to study lipid-cholesterol bilayers. In this model the lipid bilayer is represented as a two-dimensional lattice field in the lipid chain order parameters, while cholesterol molecules are represented by hard rods. The motion of rods in the system is continuous and is not confined to lattice cells. The statistical mechanics of chain ordering is described by a mean field derived from an extension of a model due to Marcelja. The time evolution of the system is governed by stochastic equations. The ensemble of chain configurations required in partition sums, and the energies of interaction, are taken from atomistic level molecular-dynamics simulations of lipid bilayers. The model allows us to simulate systems 500 nm in lateral size for 20 micros time scales, or greater. We have applied the model to dipalmitoyl-phosphatidylcholine-cholesterol (Chol) bilayers at 50 degrees C for Chol concentrations between 2% and 33%. At low concentrations of Chol (2%-4%), the model predicts the formation of isolated clusters of Chol surrounded by relatively ordered lipid chains, randomly dispersed in the disordered bilayer. With increasing Chol composition, regions of Chol-induced order begin to overlap. Starting from about 11% Chol this ordering effect becomes system wide and regions unaffected by Chol are no longer detectable. From the analysis of properties of the model we conclude that the change in lipid chain order with increasing Chol concentration is continuous over the 20-mus scale of the simulations. We also conclude that at 50 degrees C no large-scale Chol-rich and Chol-depleted coexisting phase-separated regions form at any concentration. At no point in any of the simulations do we observe a higher degree of lateral organization, such as Chol-based superlattice structures.     

Low-temperature vapor-liquid equilibria from parallelized molecular dynamics simulations. Application to 1- and 2-methylnaphthalene

A parallelized sampling version of the Gibbs Ensemble (Mol. Phys. 2000, 98, 1887) has been implemented to predict low-temperature vapor-liquid equilibria of 1- and 2-methylnaphthalene modeled by anisotropic united atom potentials. The simulation were performed at the low temperature of 364.2 K at which common direct simulation methods fail due to particle transfer problems. The simulation results are compared with published results obtained from the Gibbs-Duhem integration method and with experimental data. Both methods are compared and discussed in terms of computational efficiency and with respect to their future use at other thermodynamic conditions.     

Many-body interaction analysis: algorithm development and application to large molecular clusters

A completely automated algorithm for performing many-body interaction energy analysis of clusters (MBAC) [M. J. Elrodt and R. J. Saykally, Chem. Rev. 94, 1975 (1994); S. S. Xantheas, J. Chem. Phys. 104, 8821 (1996)] at restricted Hartree-Fock (RHF)/MA Plesset 2nd order perturbation theory (MP2)/density functional theory (DFT) level of theory is reported. Use of superior guess density matrices (DM's) for smaller fragments generated from DM of the parent system and elimination of energetically insignificant higher-body combinations, leads to a more efficient performance (speed-up up to 2) compared to the conventional procedure. MBAC approach has been tested out on several large-sized weakly bound molecular clusters such as (H(2)O)(n), n=8, 12, 16, 20 and hydrated clusters of amides and aldehydes. The MBAC results indicate that the amides interact more strongly with water than aldehydes in these clusters. It also reconfirms minimization of the basis set superposition error for large cluster on using superior quality basis set. In case of larger weakly bound clusters, the contributions higher than four body are found to be repulsive in nature and smaller in magnitude. The reason for this may be attributed to the increased random orientations of the interacting molecules separated from each other by large distances.     

Thermal Gaussian molecular dynamics for quantum dynamics simulations of many-body systems: application to liquid para-hydrogen

A new method, here called thermal Gaussian molecular dynamics (TGMD), for simulating the dynamics of quantum many-body systems has recently been introduced [I. Georgescu and V. A. Mandelshtam, Phys. Rev. B 82, 094305 (2010)]. As in the centroid molecular dynamics (CMD), in TGMD the N-body quantum system is mapped to an N-body classical system. The associated both effective Hamiltonian and effective force are computed within the variational Gaussian wave-packet approximation. The TGMD is exact for the high-temperature limit, accurate for short times, and preserves the quantum canonical distribution. For a harmonic potential and any form of operator A?, it provides exact time correlation functions C(AB)(t) at least for the case of B, a linear combination of the position, x, and momentum, p, operators. While conceptually similar to CMD and other quantum molecular dynamics approaches, the great advantage of TGMD is its computational efficiency. We introduce the many-body implementation and demonstrate it on the benchmark problem of calculating the velocity time auto-correlation function for liquid para-hydrogen, using a system of up to N = 2592 particles.     

Mapping long-range interactions in alpha-synuclein using spin-label NMR and ensemble molecular dynamics simulations

The intrinsically disordered protein alpha-synuclein plays a key role in the pathogenesis of Parkinson's disease (PD). We show here that the native state of alpha-synuclein consists of a broad distribution of conformers with an ensemble-averaged hydrodynamic radius significantly smaller than that expected for a random coil structure. This partial condensation is driven by interactions between the highly charged C-terminus and a large hydrophobic central region of the protein sequence. We suggest that this structure could inhibit the formation of alpha-synuclein aggregates, which are thought to be the cytotoxic species responsible for neurodegeneration in PD.     

Grand canonical ensemble molecular dynamics simulation of water solubility in polyamide-6,6

Grand canonical ensemble molecular dynamics simulation is employed to calculate the solubility of water in polyamide-6,6. It is shown that performing two separate simulations, one in the polymeric phase and one in the gaseous phase, is sufficient to find the phase coexistence point. In this method, the chemical potential of water in the polymer phase is expanded as a first-order Taylor series in terms of pressure. Knowing the chemical potential of water in the polymer phase in terms of pressure, another simulation for water in the gaseous phase, in the grand canonical ensemble, is done in which the target chemical potential is set in terms of pressure in the gas phase. The phase coexistence point can easily be calculated from the results of these two independent simulations. Our calculated sorption isotherms and solubility coefficients of water in polyamide-6,6, over a wide range of temperatures and pressures, agree with experimental data.     

Conformational sampling and ensemble generation by molecular dynamics simulations: 18-Crown-6 as a test case

Using the crown ether 18-crown-6 as a test system, molecular dynamics has been evaluated as a technique for conformational searching and thermodynamic ensemble generation. By running a series of 200 ps and 2 ns simulations, an “optimum” temperature range for conformational searching, i.e., the temperature at which one finds the largest number of low energy structures, was demonstrated to be dependent on the time interval at which one examines the structure. By considering conformational degeneracy and entropy with the rigid rotor harmonic oscillator approximation we have been able to demonstrate that the ensemble generated approaches thermodynamic equilibrium in about 6 ns of simulation. To our knowledge this is the first time this has been demonstrated for a complex organic molecule and it highlights the power and usefulness of molecular dynamics as a method for thermodynamic ensemble generation and conformational searching.     

Wavelet transform analysis of ab initio molecular dynamics simulation: application to core-excitation dynamics of BF3

We propose a novel analysis method of ab initio molecular dynamics (AIMD) simulation using a continuous wavelet transform (c-WT) technique. The c-WT technique, one of the time-frequency signal analysis methods, provides a clear view of the dynamical information in time developments. Combined with the auto-correlation function of velocity by AIMD simulation, c-WT analysis enables us to well understand dynamical distribution, such as the vibrational properties following a change of electronic structure in a molecular system. As a practical application, AIMD simulation of core-excited BF(3) (B1s --> 2a(2) (')) is illustrated. AIMD simulation leads to the change of vibrational motion as well as structural deformation by core-excitation. The c-WT analysis clarifies the relationship between structural deformation and the related significant vibrational modes in core-excitation within 50 fs.     

Synergistic effect of dicarbollide anions in liquid-liquid extraction: a molecular dynamics study at the octanol-water interface

We report a molecular dynamics study of chlorinated cobalt bis(dicarbollide) anions [(B(9)C(2)H(8)Cl(3))(2)Co](-)"CCD(-)" in octanol and at the octanol-water interface, with the main aim to understand why these hydrophobic species act as strong synergists in assisted liquid-liquid cation extraction. Neat octanol is quite heterogeneous and is found to display dual solvation properties, allowing to well solubilize CCD(-), Cs(+) salts in the form of diluted pairs or oligomers, without displaying aggregation. At the aqueous interface, octanol behaves as an amphiphile, forming either monolayers or bilayers, depending on the initial state and confinement conditions. In biphasic octanol-water systems, CCD(-) anions are found to mainly partition to the organic phase, thus attracting Cs(+) or even more hydrophilic counterions like Eu(3+) into that phase. The remaining CCD(-) anions adsorb at the interface, but are less surface active than at the chloroform interface. Finally, we compare the interfacial behavior of the Eu(BTP)(3)(3+) complex in the absence and in the presence of CCD(-) anions and extractant molecules. It is found that when the CCD(-)'s are concentrated enough, the complex is extracted to the octanol phase. Otherwise, it is trapped at the interface, attracted by water. These results are compared to those obtained with chloroform as organic phase and discussed in the context of synergistic effect of CCD(-) in liquid-liquid extraction, pointing to the importance of dual solvation properties of octanol and of the hydrophobic character of CCD(-) for synergistic extraction of cations.     

Vapor-liquid and liquid-liquid equilibria of water, 2-methoxyethanol and cyclohexanone: Experiment and correlation

Vapor-liquid equilibria and liquid-liquid equilibria of a ternary mixture consisting of water, 2-methoxyethanol and cyclohexanone and in addition of all binary subsystems were studied experimentally at several temperatures. A ternary corrective term in the expression for the Gibbs free energy based on the NRTL model improves simultaneous representation of binary and ternary phase equilibria.