A study on dissociation of sII krypton hydrate and the effect of hydrocarbon guest molecules as stabilizer by molecular dynamics simulation

In this paper, the characteristics of structure II krypton hydrate are studied by molecular dynamics simulation under isobaric-isothermal (NPT) ensemble condition. The dissociation process of the hydrate is simulated and the effect of krypton (Kr) and various types of hydrocarbon guest molecules (HGMs) on the stability of the hydrate structure is investigated during the simulation time of 1 ns. The studied HGMs are propane, isobutane, neopentane, cyclopropane, cyclobutane, cyclopentane and cyclopentene. The structural change of the Kr-hydrate is analyzed with the radial distribution function, mean square displacement and diffusion coefficient. As temperature increases, the obtained results indicate a gradual increase in the Kr-hydrate cell size, which leads to distortion of the hydrate lattice and escaping of the encapsulated Kr molecules from the hydrate structure to form small bubbles of Kr aggregated in the aqueous solution.     

Molecular simulation of chemical reaction equilibria by Kinetic Monte Carlo

ABSTRACT

We describe a new algorithm for the molecular simulation of chemical reaction equilibria, which we call the Reactive Kinetic Monte Carlo (ReKMC) algorithm. It is based on the use of the equilibrium Kinetic Monte Carlo (eKMC) method (Ustinov et al., J. Colloid Interface Sci., 2012, 366, 216–223) to generate configurations in the underlying nonreacting system and to calculate the species chemical potentials at essentially zero marginal computational cost. We consider in detail the typical case of specified temperature, T and pressure, P, but extensions to other thermodynamic constraints are straightforward in principle. In the course of this work, we also demonstrate an alternative method for calculating simulation box volume changes in NPT ensemble simulations to achieve the specified P. We consider two sets of example reacting systems previously considered in the literature, and compare the ReKMC results and computational efficiencies with those of different implementations of the REMC algorithm (Turner et al., Molec. Simulation, 2008, 34, 119–146).     

Vapour—liquid equilibrium of the square-well fluid of variable range via a hybrid simulation approach

The equilibrium between vapour and liquid in a square-well system has been determined by a hybrid simulation approach combining chemical potentials calculated via the Gibbs ensemble Monte Carlo technique with pressures calculated by the standard NVT Monte Carlo method. The phase equilibrium was determined from the thermodynamic conditions of equality of pressure and chemical potential between the two phases. The results of this hybrid approach were tested by independent NPT and μPT calculations and are shown to be of much higher accuracy than those of conventional GEMC simulations. The coexistence curves, vapour pressures and critical points were determined for SW systems of interaction ranges λ = 1.25, 1.5, 1.75 and 2. The new results show a systematic dependence on the range λ, in agreement with results from perturbation theory where previous work had shown more erratic behaviour.     

Transitioning model potentials to real systems

The parameters of two pair potentials that describe argon over its entire liquid phase at a fixed pressure were optimized through a novel application of constant temperature and pressure molecular dynamics (NPT-MD) and Monte Carlo (NPT-MC) computer simulations. The forms of these potentials were those of a modified Lennard-Jones potential and a Lennard-Jones potential. The optimized potential determined using NPT-MD simulations reproduces experimental densities, internal energies and enthalpies with an error less than 1% over most of the liquid range and yields self-diffusion coefficients that are in excellent agreement with experiment. The results using the potential determined by NPT-MC simulations are in almost as good agreement with deviations from experiment of no more than 5.89% for temperatures up to vaporization. Additionally, molar volumes predicted using this potential at pressures in the range 100–600 atm and over temperatures in the range 100–140K were within 0.83% of experimental values. These results show that, when properly parametrized, Lennard-Jones-like potentials can describe a system well over a large temperature range. Further, the method introduced is easy to implement and is independent of the form of the interaction potential used.     

Direct determination of phase coexistence properties of fluids by Monte Carlo simulation in a new ensemble

A methodology is presented for Monte Carlo simulation of fluids in a new ensemble that can be used to obtain phase coexistence properties of multicomponent systems from a single computer experiment. The method is based on performing a simulation simultaneously in two distinct physical regions of generally different densities and compositions. Three types of perturbations are performed, a random displacement of molecules that ensures equilibrium within each region, an equal and opposite change in the volume of the two regions that results in equality of pressures, and random transfers of molecules that equalize the chemical potentials of each component in the two regions. The method is applied to the calculation of the liquid-gas coexistence envelope for the pure Lennard-Jones (6, 12) fluid for several reduced temperatures from the vicinity of the triple point to close to the critical point (T? = 0.75 to T? = 1.30). Good overall agreement with previously available literature results is obtained, with some deviations at the extremes of this temperature range.     

Atomistic simulations of liquid water using Lekner electrostatics

Equilibrium molecular dynamics simulations have been performed for liquid water using three different potential models in the NVT and NPT ensembles. The flexible SPC model, the rigid TIP4P model and the rigid/polarizable TIP4P-FQ potential were studied. The Lekner method was used to handle long range electrostatic interactions, and an efficient trivariate cubic spline interpolation method was devised for this purpose. A partitioning of the electrostatic interactions into medium and long range parts was performed, and the concomitant use of multiple timestep techniques led to substantially enhanced computation speeds. The simulations were carried out using 256 molecules in the NVT ensemble at 25°C and 997 kg m^?3 and in the NPT ensemble at 25°C and 1 bar. Various dynamic, structural, dielectric, rotational and thermodynamic properties were calculated, and it was found that the simulation methodologies performed satisfactorily vis-à-vis previous simulation results and experimental observations.     

Accurate Simulation Estimates of Phase Behavior in Ternary Mixtures with Prescribed Composition

This paper describes an isobaric semi-grand canonical ensemble Monte Carlo scheme for the accurate study of phase behavior in ternary fluid mixtures under the experimentally relevant conditions of prescribed pressure, temperature and overall composition. It is shown how to tune the relative chemical potentials of the individual components to target some requisite overall composition and how, in regions of phase coexistence, to extract accurate estimates for the compositions and phase fractions of individual coexisting phases. These estimates have finite-size errors that are exponentially small in the system size even when the phase fraction of one phase is vanishingly small, as occurs at the coexistence boundary. The method is illustrated by tracking a path through the composition space of a model ternary Lennard-Jones mixture.     

A new simulation method for the determination of phase equilibria in mixtures in the grand canonical ensemble

A grand canonical Monte Carlo (GCMC) simulation method is presented for the determination of the phase equilibria of mixtures. The coexistence is derived by expanding the pressure into a Taylor series as a function of the temperature and the chemical potentials that are the independent intensive variables of the grand canonical ensemble. The coefficients of the Taylor series can be calculated from ensemble averages and fluctuation formulae that are obtained from GCMC simulations in both phases. The method is able to produce the equilibrium data in a certain domain of the (T, p) plane from two GCMC simulations. The vapour-liquid equilibrium results obtained for a Lennard-Jones mixture agree well with the corresponding Gibbs ensemble Monte Carlo data.     

Site factors in N.M.R. solvent effects

An algorithm was developed enabling implementation of a Nosé—Hoover thermostat within the framework of grand canonical molecular dynamics [Lynch, C. G. and Pettitt, B. M., 1997, J. chem. Phys., 107, 8594]. The proposed algorithm could readily be extended to mixtures of molecular species with different chemical potentials as shown in the paper. This algorithm was first applied to simulate a μVT ensemble of TIP4P water molecules at 298 K by means of a system comprising a number of full particles and a single scaled (fractional) particle, with the scaling factor considered as a dynamic variable in its own right and chemical potential a pre-set parameter. Our finding showed that the scheme with a single fractional particle tended to freeze in metastable states as well as failed to reproduce either the real-life (?24.05 kJmol^?1) or the model-specific chemical potential of water (?23.0kJ mol^?1). In order to overcome the inadequacy of a single fractional particle as a chemical potential ‘probe’ the treatment of Pettitt and co-workers was extended to introduce multiple fractional particles. The extended scheme (with 4 fractional particles) was able to reproduce the actual density of water for the driving chemical potential of -24.0k mo^?1. The actual behaviour of the density as a function of the chemical potential also agreed quite well with both the results of thermo-dynamic integration and the findings of Pettitt and co-workers.     

The extrapolation of phase equilibrium curves of mixtures in the isobaric—isothermal Gibbs ensemble

A powerful extrapolation scheme is proposed to determine the vapour—liquid and liquid—liquid equilibrium curves of mixtures by performing a single isothermal-isobaric Gibbs ensemble Monte Carlo (GEMC) simulation. The coexistence curves for the mole fraction and the density are extrapolated as functions of the temperature and the pressure by second-order Taylor series. The coefficients of the Taylor series, which are the temperature and pressure derivatives of these quantities along the coexistence curves, can be calculated from the data produced by a single GEMC simulation on the basis of fluctuation formulas. We show that the application of a Padé approximant considerably widens the temperature and pressure range where the extrapolation is accurate. Using Lennard-Jones mixtures as test systems, we show that the technique is able to produce quite accurate equilibrium curves at fixed temperature in the function of the pressure and vice versa. The procedure yields good results not only for vapour—liquid but also for liquid—liquid coexistence curves. The calculation of the vapour pressure curves at a fixed composition of the liquid side is straightforward with the method.     

Grand Equilibrium: vapour-liquid equilibria by a new molecular simulation method

A new molecular simulation method for the calculation of vapour-liquid equilibria of mixtures is presented. In this method, the independent thermodynamic variables are temperature and liquid composition. In the first step, one isobaric isothermal simulation for the liquid phase is performed, in which the chemical potentials of all components and their derivatives with respect to the pressure, i.e. the partial molar volumes, are calculated. From these results, first-order Taylor series expansions for the chemical potentials as functions of the pressure μ _i (p) at constant liquid composition are determined. This information is needed, as the specified pressure in the liquid will generally not be equal to the equilibrium pressure, which has to be found in the course of a vapour simulation. In the second step, one pseudo grand canonical simulation for the vapour phase is performed, where the chemical potentials are set according to the instantaneous pressure p ^v using the previously determined function μ _i (p ^v). In this way, results for the vapour pressure and vapour composition are achieved which are consistent for the given temperature and liquid composition. The new method is applied to the pure Lennard-Jones fluid, a binary and a ternary mixture of Lennard-Jones spheres, and shows very good agreement with corresponding data from the literature.     

Studies of a primitive model of mixtures of nonpolar molecules with water

The paper presents calculations of the properties of binary mixtures of hard spheres and directionally associating hard spheres, a simple model for mixtures of nonpolar molecules with water that was developed by Nezbeda and his coworkers. Extensive results from Monte Carlo simulations in the isobaric, isothermal ensemble are presented for the density, configurational energy and chemical potentials in the mixtures for fluid states over a range of temperatures, pressures and compositions. A species exchange technique is used to compute the chemical potential difference between components in the mixtures. The results obtained are compared with the predictions of first-order thermodynamic perturbation theory (TPT). It is found that this theory provides an accurate picture of the system over most of the conditions considered. Calculations are also made of vapour–liquid coexistence for the model using TPT and calculations of solid–fluid coexistence for the model using TPT and existing results for the free energy of the pure component solids. It is found that the vapour–liquid coexistence for the model is pre-empted by the solid–fluid coexistence, as had previously been found for the pure component directionally associating hard sphere system.     

Measurement of chemical potentials of systems with strong excluded volume interactions by computing the density of states

At high densities and low temperatures, the conventional Widom test particle method to compute the chemical potential of a system of particles with excluded volume interactions fails owing to bad statistics. A way to circumvent this problem is the use of expanded ensemble simulation techniques or thermodynamic integration. In this article, we will describe an alternative method to compute the chemical potential which is conceptually much easier, by computing the density of states of systems of N and N + 1 particles directly; and by performing a test particle simulation at very high temperature. The advantage of our technique is that the densities of states of the N and N + 1 particle system are computed in an ensemble in which particles can pass each other, resulting in a more efficient sampling. We will demonstrate our method not only for single particles but also for chain molecules with intramolecular interactions. By using an infinite temperature expansion and an extension of the density of states to very high energies, we will show that it is also possible to compute the chemical potential without having to compute the density of states for the N + 1 particle system.     

Monte Carlo molecular simulation of the Na-, Mg-, and mixtures of Na/Mg-montmorillonites systems,in function of the pressure

In this paper, Monte Carlo (MC) simulation has been used to study the swelling pattern of Na-montmorillonite (Na-Mnt), Mg-montmorillonite (Mg-Mnt), and Na/Mg-mixture montmorillonite (4NaMg-Mnt; 2Na2Mg-Mnt). The molecular simulation was performed in the NVT (number of molecules, volume and temperature are constant) ensemble at normal temperature (300 K) and 225, 300, and 340 bar over an H₂O content 147, 196, and 294 mg g^?1 of clay. The simulations reproduce the swelling pattern of Na-Mnt and Mg-Mnt. The predicted spacing of the Na/Mg-Mnt mixtures is closely related to that of Mg-Mnt and confirms the results reported in the literature for Na-rich/Mg-poor Mnt. The results of the water adsorption and the swelling properties on the system Na-Mnt, Mg-Mnt, and the Na/Mg-Mnt mixtures are reflected with a transformation to two-hydrate stages. The probability of the coordination number of Na⁺, Mg²⁺, and mixtures tends to increase with an increasing amount of H₂O molecules, but decreases with increasing pressure. The cation–oxygen distances (Na–O or Mg–O) show two signals, corresponding to the first and second coordination shells, which indicates that the ions behave as in bulk water.     

Cavity-bias sampling in reaction ensemble Monte Carlo simulations

A methodology is presented to sample efficiently configurations of reacting mixtures in the reaction ensemble Monte Carlo simulation technique. A cavity-biasing scheme is used, which more effectively samples configurations than conventional random sampling. Akin to other biasing schemes that are implemented into insertion-based Monte Carlo methods such as Gibbs ensemble Monte Carlo, the method presented here searches for space in the reacting mixture whereby the insertion of a product molecule is energetically favoured. This sampling bias is then corrected in the acceptance criteria. The approach allows for the study of reacting mixtures at high density as well as for polyatomic molecular species. For some cases, the method is shown to increase the efficiency of the reaction steps by a factor greater than 20. The approach is shown to be readily generalized to other biasing schemes such as orientational-biasing of polar molecules and configurational-biasing of chain-like molecules.     

Mixtures of hard spherocylinders and hard spheres

Monte Carlo simulations have been performed for equimolar mixtures of hard prolate spherocylinders of length: breadth ratio 2:1 and hard spheres, in the fluid region. Two systems have been studied. In the first the breadth of the spherocylinder was equal to the hard sphere diameter, and in the second system both components were of equal molecular volume.

The compressibility factor, PV/NkT, has been obtained for both mixtures at four reduced densities (packing fractions) from 0·20 to 0·45. The results have been compared with the predictions of several analytical equations appropriate to mixtures of hard convex molecules, and an equation due to Pavlicek et al. was found to be very accurate. The results have been used to calculate the excess volumes of mixing at constant pressure, in an attempt to establish the relative importance of the effects of differences in molecular volume and shape on the thermodynamic properties.

The structural properties of the mixtures have also been investigated by calculating pair distribution functions for the three types of pair interactions present in these mixtures.     

Monte Carlo simulations of neon and argon using ab initio potentials

Gibbs ensemble Monte Carlo simulations of neon and argon have been performed with pair potentials taken from literature as well as with new ab initio potentials from just above the triple point to close to the critical point. The densities of the coexisting phases, their pair correlation functions, the vapour pressure and the enthalpy and entropy of vaporization have been calculated. The influence of the potential choice and of the addition of the Axilrod-Teller (AT) three-body potential on the above mentioned properties have been investigated. It turns out that an accurate ab initio two-body potential in connection with the AT potential yields very good results for thermodynamic properties of phase equilibria.     

The O-H distance in ice

The thermodynamic properties of normal and para-hydrogen are computed from multiple time-step path integral hybrid Monte Carlo (PIHMC) simulations. Four different isotropic pair potentials are evaluated by comparing simulation results with experimental data. The Silvera–Goldman potential is found to be the most accurate of the potentials tested for computing the density and internal energy of fluid hydrogen. Using the Silvera–Goldman potential, simulation and experimental data are compared on isobars ranging from 0.1 to 100 MPa and for temperatures from 18 to 300 K. The Gibbs free energy is calculated from the PIHMC simulations by an adaptation of Widom's particle insertion technique to a path integral fluid. A new method is developed for computing phase equilibria for quantum fluids directly by combining PIHMC with the Gibbs ensemble technique. This Gibbs–PIHMC method is used to calculate the vapour–liquid phase diagram of hydrogen from simulations. Agreement with experimental data is good.     

Thermodynamic and transport properties of simple fluids using lattice sums: bulk phases and liquid-vapour interface

Molecular dynamics simulations in the canonical ensemble have been performed to obtain the thermodynamic and transport properties of the Lennard-Jones fluid. The dispersion interactions were calculated using lattice sums. This method makes it possible to simulate the full potential avoiding the inclusion of the long range corrections (LRC) during or at the end of simulations. In the calculation of dynamic properties in bulk phases and thermodynamic quantities of inhomogeneous systems where the interface is physically present, in general the LRC cannot easily be included. By using the lattice sums method, the results are independent of the truncation of the potential. In the liquid-vapour interface simulations it is not necessary to make any pre-judgments about the form of the LRC formula to calculate coexisting properties such as the surface tension. The lattice sums method has been applied to evaluate how well the full interaction can be calculated in the liquid phase and in the liquid-vapour interface. In the liquid phase the pressure, configurational energy, diffusion coefficient and shear viscosity were obtained. The results of the thermodynamic properties are compared with those obtained using the spherically truncated and shifted (STS) potential with the LRC added at the end of simulations, and excellent agreement is found. The transport properties are calculated on different system sizes for a state near the triple point. The diffusion coefficient using the lattice sums method increases with the number of molecules, and the results are higher than those of the STS model truncated at 2.5σ (STS2.5). The shear viscosity does not show any system size dependence for systems with more than 256 molecules, and the lattice sums results are essentially the same as those for the STS2.5. In the liquid-vapour equilibria the coexisting densities and vapour pressures for the full potential agree well with those obtained using the Gibbs ensemble and the NPT + test particle methods. The surface tension using lattice sums and truncation of forces at 2.5σ agrees well with STS results using large system sizes and cutoff distances.