Adsorption of carbon dioxide of 1-site and 3-site models in pillared clays: A Gibbs ensemble Monte Carlo simulation

We report Canonical (NVTMC), Grand Canonical (GOMC) and Gibbs Ensemble Monte Carlo (GEMC) simulations for adsorption of methane on graphite. Lennard-Jones (LJ) potentials are used for the intermolecular interactions, and both structured and structureless (10–4–3) solid-fluid potentials are considered. Several sets of methane-methane L.I parameters have been used in the literature, and we compare results obtained with these sets throughout our simulations. The adsorption isotherm and isosteric heat curve are obtained at 77.5 K and found in good agreement with experiments. The commensurateincommensurate transition (CIT) of methane on a graphite substrate with periodically varying adsorbate-adsorbent potential at 40.0 K is studied and is in qualitative agreement with experiment. The effect of varying the corrugation of the fluid-wall potential on the commensurate and incommensurate phases is explored. The GEMC simulations have been carried out to study the vapor-liquid equilibrium (VLE) of a two-dimensional (2D) LJ fluid with system sizes up to 3000 particles. The effect of system sizes on the critical behavior is investigated. The GEMC method has also been successfully applied to study the VLE in 2D adsorbed films for the first time.

Since the seminal paper by Panagiotopoulos [Mol. Phys. 61, 813 (1997)], the Gibbs ensemble Monte Carlo (GEMC) method has been the most popular particle-based simulation approach for the computation of vapor-liquid phase equilibria. However, the validity of GEMC simulations in the near-critical region has been questioned because rigorous finite-size scaling approaches cannot be applied to simulations with fluctuating volume. Valleau [Mol. Simul. 29, 627 (2003)] has argued that GEMC simulations would lead to a spurious overestimation of the critical temperature. More recently, Patel et al. [J. Chem. Phys. 134, 024101 (2011)] opined that the use of analytical tail corrections would be problematic in the near-critical region. To address these issues, we perform extensive GEMC simulations for Lennard-Jones particles in the near-critical region varying the system size, the overall system density, and the cutoff distance. For a system with N = 5500 particles, potential truncation at 8σ and analytical tail corrections, an extrapolation of GEMC simulation data at temperatures in the range from 1.27 to 1.305 yields T(c) = 1.3128 ± 0.0016, ρ(c) = 0.316 ± 0.004, and p(c) = 0.1274 ± 0.0013 in excellent agreement with the thermodynamic limit determined by Potoff and Panagiotopoulos [J. Chem. Phys. 109, 10914 (1998)] using grand canonical Monte Carlo simulations and finite-size scaling. Critical properties estimated using GEMC simulations with different overall system densities (0.296 ≤ ρ(t) ≤ 0.336) agree to within the statistical uncertainties. For simulations with tail corrections, data obtained using r(cut) = 3.5σ yield T(c) and p(c) that are higher by 0.2% and 1.4% than simulations with r(cut) = 5 and 8σ but still with overlapping 95% confidence intervals. In contrast, GEMC simulations with a truncated and shifted potential show that r(cut) = 8σ is insufficient to obtain accurate results. Additional GEMC simulations for hard-core square-well particles with various ranges of the attractive well and for n-decane molecules represented by the TraPPE force field yield data that support the trends observed for Lennard-Jones particles. The finite-size dependence of the critical properties obtained from GEMC simulations is significantly smaller than those from grand-canonical ensemble simulations. Thus, when resources are not available for a rigorous finite-size scaling study, GEMC simulations provide a straightforward route to determine fairly accurate critical properties using relatively small system sizes.

We have employed two-to-one mapping scheme to develop three coarse-grained (CG) water models, namely, 1-, 2-, and 3-site CG models. Here, for the first time, particle swarm optimization (PSO) and gradient descent methods were coupled to optimize the force-field parameters of the CG models to reproduce the density, self-diffusion coefficient, and dielectric constant of real water at 300 K. The CG MD simulations of these new models conducted with various timesteps, for different system sizes, and at a range of different temperatures are able to predict the density, self-diffusion coefficient, dielectric constant, surface tension, heat of vaporization, hydration free energy, and isothermal compressibility of real water with excellent accuracy. The 1-site model is ∼3 and ∼4.5 times computationally more efficient than 2- and 3-site models, respectively. To utilize the speed of 1-site model and electrostatic interactions offered by 2- and 3-site models, CG MD simulations of 1:1 combination of 1- and 2-/3-site models were performed at 300 K. These mixture simulations could also predict the properties of real water with good accuracy. Two new CG models of benzene, consisting of beads with and without partial charges, were developed. All three water models showed good capacity to solvate these benzene models.

The effects of solid-fluid interactions on the vapor-liquid phase diagram, coexistence density, relative volatility and vaporization enthalpy have been investigated for confined binary systems of CO2-CH4, CO2-N2 and CH4-N2. The Gibbs ensemble Monte Carlo (GEMC) simulation results indicate that the confinement and the solid-fluid interaction have significant influences on the vapor-liquid equilibrium properties. The confinement and the strength of the solid-fluid interaction make the p-xi phase diagram move to higher pressure regions. They also make the two-phase region become narrower for each binary mixture. The strength of the solid-fluid interactions can cause increases in the coexistence liquid and vapor densities, and cause the decrease of the relative volatility and the vaporization enthalpy for the systems studied. As the pore width is decreased, the two-phase region of the binary mixture becomes narrower.

TraPPE-UA and OPLS-AA force fields were used to investigate the applications of the Gibbs ensemble Monte Carlo (GEMC) method and the Gibbs-Duhem integration (GDI) method in predicting the vapor-liquid equilibrium properties. It was found that the GDI method was much faster than the GEMC method if all-atom force field was employed. Based on the calculation results, the two methods could be viewed as complementary to the problem of vapor liquid coexistence. The calculated liquid densities, enthalpies of vaporization, critical temperatures, and critical densities by the two methods were very close to each other for a given force field. If the force field caused errors in the calculated heats of vaporization, the calculated vapor pressures and densities by the two methods were clearly different, which also led to difference in the predicted critical pressure.

Gibbs ensemble Monte Carlo (GEMC) simulations in the isochoric–isothermal (NVT) ensemble were used to simulate vapour–liquid–liquid equilibrium (VLLE) for binary n-hexane–water and ethane–ethanol mixtures. The GEMC simulation of binary VLLE data proved to be extremely difficult and that is probably the reason why the open literature is so sparse with simulations for these types of systems. The results presented in this paper are to our knowledge the first successful binary three-phase GEMC simulations of non-idealised fluid systems. This paper also shows that the isobaric–isothermal (NPT) ensemble is unsuitable for the simulation of phase equilibria of binary three-phase systems.

Molecular simulation of the impact of surface roughness on carbon dioxide adsorption in organic-rich shales

The Gibbs ensemble Monte Carlo (GEMC) method is a versatile approach for the prediction of fluid phase equilibria from particle-based simulations. For a one-component system, a GEMC simulation utilises two separate simulation boxes for the vapour and the liquid phases and a significant fraction of the computational effort is expended on special trial moves that transfer (swap) particles and exchange volume between the two boxes. The user needs to specify the frequency of swap and volume moves and the overall volume that controls the phase ratio. In this study, the efficiency of GEMC simulation protocols that yield three different frequencies of accepted swap and volume moves and three different phase ratios is assessed for the computation of the saturated vapour pressure and liquid density of n-octane and water at three reduced temperatures. Differences in the simulation efficiency of up to an order of magnitude are observed, and recommendations are made for suitable GEMC simulation protocols.

We present evidence that the transition between organic and third phases, which can be observed in the plutonium uranium reduction extraction (PUREX) process at high metal loading, is an unusual transition between two isotropic bicontinuous microemulsion phases. As this system contains so many components, however, we have been seeking first to investigate the properties of a simpler system, namely, the related metal-free, quaternary water/n-dodecane/nitric acid/tributyl phosphate (TBP) system. This quaternary system has been shown to exhibit, under appropriate conditions, three coexisting phases: a light organic phase, an aqueous phase, and the so-called third phase. In the current work, we focused on the coexistence of the light organic phase with the third phase. Using Gibbs ensemble Monte Carlo (GEMC) simulations, we found coexistence of a phase rich in nitric acid and dilute in n-dodecane (the third phase) with a phase more dilute in nitric acid but rich in n-dodecane (the light organic phase). The compositions and densities of these two coexisting phases determined using the simulations were in good agreement with those determined experimentally. Because such systems are generally dense and the molecules involved are not simple, the particle exchange rate in their GEMC simulations can be rather low. To test whether a system having a composition between those of the observed third and organic phases is indeed unstable with respect to phase separation, we used the Bennett acceptance ratio method to calculate the Gibbs energies of the homogeneous phase and the weighted average of the two coexisting phases, where the compositions of these phases were taken both from experimental results and from the results of the GEMC simulations. Both demixed states were determined to have statistically significant lower Gibbs energies than the uniform, mixed phase, providing confirmation that the GEMC simulations correctly predicted the phase separation. Snapshots from the simulations and a cluster analysis of the organic and third phases revealed structures akin to bicontinuous microemulsion phases, with the polar species residing within a mesh and with the surface of the mesh formed by amphiphilic TBP molecules. The nonpolar n-dodecane molecules were observed in these snapshots to be outside this mesh. The only large-scale structural differences observed between the two phases were the dimensions of the mesh. Evidence for the correctness of these structures was provided by the results of small-angle X-ray scattering (SAXS) studies, where the profiles obtained for both the organic and third phases agreed well with those calculated from simulations. Finally, we looked at the microscopic structures of the two phases. In the organic phase, the basic motif was observed to be one nitric acid molecule hydrogen-bonded to a TBP molecule. In the third phase, the most common structure was that of the hydrogen-bonded TBP-HNO3-HNO3 chain. A cluster analysis provided evidence for TBP forming an extended, connected network in both phases. Studies of the effects of metal ions on these systems will be presented elsewhere.

This article presents a methodology for checking the existence of the azeotrope and computing its composition, density, and pressure at a given temperature by integrating chemical engineering insights with molecular simulation principles. Liquid-vapor equilibrium points are computed by molecular simulations using the Gibbs ensemble Monte Carlo (GEMC) method at constant volume. The appearance of the azeotropic point is marked by a shift of the equilibrium constant from one side of the unity to the other. After each GEMC simulation, an identity change move is derived in the grand canonical ensemble to progress towards the azeotrope along the equilibrium curve. The effectiveness of the proposed methodology is successfully tested for several binary Lennard-Jones mixtures reported in the literature.

The phase equilibrium of mixtures of Yukawa and charged Yukawa particles is studied by means of Gibbs ensemble Monte Carlo (GEMC) simulation method and the mean spherical approximation (MSA). The strength of the Coulomb energy compared to that of the Yukawa attraction is characterized by a coupling constant. For low coupling constants a classical vapor--liquid phase separation appears with a good agreement between GEMC and the MSA. For high coupling constant, a phase separation between a salt poor and a salt rich phase occurs that resembles the phase equilibrium behavior of the solvent primitive model.

Extensive Gibbs ensemble Monte Carlo (GEMC) simulations of a rigid molecule model of C60, characterized by a deeply attractive short-ranged interaction potential, are performed with the aim to establish the effect of the system size on the existence and location of the liquid–vapor binodal line and of its critical point. The results obtained with N=300, 600, and 1500 particles indicate that the position and the overall shape of the binodal is only minorly influenced by finite size effects. The estimated critical temperature and density at the different N fall in the ranges 1920–1940 K, and 0.4–0.45 nm−3, respectively. The results are discussed by making reference to previous studies of finite size effects in GEMC simulations. The GEMC predictions are also compared with previous computer simulation and theoretical calculations for the same model fluid. The agreement is on the whole satisfactory for both the liquid–vapor coexistence line and the critical point parameters. On the basis of previously determined freezing lines of C60, and of the actual binodal line, different estimates of the location of the triple point are also made. Triple point temperatures are found, in any case, definitely lower (by at least 150 K) than the critical temperatures, thus confirming the existence of a relatively narrow liquid phase region in the phase diagram, as predicted in previous molecular dynamics and theoretical works. The existence of such a liquid phase for the adopted model potential is discussed and assessed in the more general framework of liquid–vapor coexistence conditions in fluids interacting through short-ranged forces. The possibility to get the liquid phase of “real life” C60, hitherto not observed experimentally, is also discussed in connection with recent high temperature experimental results on fullerite samples.

The phase diagrams of hard-core Yukawa mixtures (HCYM), constituted of equal sized hard spheres interacting through an attractive Yukawa tail, are determined by means of Gibbs Ensemble Monte Carlo (GEMC) simulations, Semi-grand Canonical Monte Carlo (SGCMC) simulations, and through the modified hypernetted-chain (MHNC) theory. Freezing lines are obtained according to an approach recently proposed by Giaquinta and co-workers [Physica A 187, 145 (1992); Phys Rev. A 45, 6966 (1992)] in which an analysis of multiparticle contributions to the excess entropy, Δs, is performed, with the determination of the Δs=0 locus. Liquid–vapor coexistence, determined through GEMC simulations, turns out to be favored when the strength ratio ν of unlike to like particle interaction, is close to 1. For lower ν’s, liquid–vapor coexistence is favored at low densities, and liquid–liquid coexistence, determined through SGCMC simulations, at high densities. The liquid–vapor binodal shifts downward in temperature and flattens when ν decreases, with a decrease of the critical temperature. At ν=0.9 a triple point can be identified from the intersection of the freezing line with the binodal line; at ν=0.7, instead, the binodal ends on the line of liquid–liquid (consolute) critical points, the intersection of the two lines thus identifying the “crossover” density and temperature between the two equilibrium regimes which correspond to the critical end point of the mixture. We find that, for not too high densities, consolute equilibrium can be also explored through GEMC simulations; the results for liquid–liquid coexistence obtained through this method and SGCMC simulations compare quite satisfactorily with each other. The trend of the liquid–vapor binodal to disappear for relatively weak unlike interactions is discussed in connection with the disappearance of liquid–vapor equilibrium which occurs in one component hard-core Yukawa fluids characterized by very short ranged attractive forces. The latter behavior has been conjectured to be relevant for the onset of crystallization in protein solutions; the implications of the present results, which are obtained in the context of a two component, albeit rough, modelization of a realistic solution, are discussed. In agreement with similar results obtained by Giaquinta et al., we finally find that the Δs=0 locus not only brings the signature of the freezing transition, but also of structural rearrangements preluding to other phase equilibria; in fact, the Δs=0 line turns out to be coincident to a high accuracy with the line of consolute critical points and with the gas branches of the liquid–vapor binodals.

A study on the vapor–liquid equilibria of vinyl acrylate-acetic acid and vinyl acetate-acetic acid with GEMC method

The vapor–liquid equilibrium (VLE) for the system acetic acid, acetaldehyde, and crotonaldehyde was investigated by combining the Gibbs ensemble Monte Carlo (GEMC) method and activity coefficient models prediction. The VLE of above pure components and their binary systems were first studied using GEMC method. Based on the strong association characteristics of acetic acid, the revised TraPPE-UA-D force field was developed with quantum chemistry methods for the VLE calculation of the acetic acid system. The original and mixed TraPPE-UA force field was introduced for the VLE calculation of the acetaldehyde and crotonaldehyde system, respectively. The results showed that the predictions agreed fairly well with the average relative deviations of the saturated densities smaller than 4.60%. Moreover, the binary VLE diagram was calculated using the above force fields. The simulation results of the mole composition were in good agreement with the experimental values, which verified the high simulation accuracy of ...