We present a new molecular simulation code, Brick-CFCMC, for performing Monte Carlo simulations using state-of-the-art simulation techniques. The Continuous Fractional Component (CFC) method is implemented for simulations in the NVT/NPT ensembles, the Gibbs Ensemble, the Grand-Canonical Ensemble, and the Reaction Ensemble. Molecule transfers are facilitated by the use of fractional molecules which significantly improve the efficiency of the simulations. With the CFC method, one can obtain phase equilibria and properties such as chemical potentials and partial molar enthalpies/volumes directly from a single simulation. It is possible to combine trial moves from different ensembles. This enables simulations of phase equilibria in a system where also a chemical reaction takes place. We demonstrate the applicability of our software by investigating the esterification of methanol with acetic acid in a two-phase system.

@en

Computer simulations of adsorption of aromatics in zeolites are typically performed using rigid zeolite frameworks. However, adsorption isotherms for aromatics are very sensitive to small differences in the atomic positions of the zeolite (Chem. Phys. Lett., 1999, 308, 155-159). This article studies the effect of framework flexibility on the adsorption of aromatics in MFI-type zeolites computed by grand-canonical Monte Carlo simulations. New experimental data of adsorption of ethylbenzene in a MFI-type zeolite at 353 K is presented. The adsorption of n-heptane, ethylbenzene, and xylene isomers is computed in three MFI-type zeolite structures. It is observed that the intraframework interactions in flexible framework models induce small but important changes in the atom positions of the zeolite and hence in the adsorption isotherms. Framework flexibility is differently "rigid": flexible force fields produce a zeolite structure that vibrates around a new equilibrium configuration with limited capacity to accommodate to a bulky guest molecule. The vibration of the zeolite atoms only plays a role at high loadings, and the adsorption is mainly dependent on the average positions of the atoms. The simulations show that models for framework flexibility should not be blindly applied to zeolites and a general reconsideration of the parametrization schemes for such models is needed.

@en

Partial molar properties are of fundamental importance for understanding properties of non-ideal mixtures. Josephson and co-workers (Mol. Phys. 2019, 117, 3589–3602) used least squares multiple linear regression to obtain partial molar properties in open constant-pressure ensembles. Assuming composition-independent partial molar properties for the narrow composition range encountered throughout simulation trajectories, we rigorously prove the equivalence of two approaches for computing thermodynamic derivatives in open ensembles of an n-component system: (1) multiple linear regression, and (2) thermodynamic fluctuations. Multiple linear regression provides a conceptually simple and computationally efficient way of computing thermodynamic derivatives for multicomponent systems. We show that in the reaction ensemble, the reaction enthalpy can be computed directly by simple multiple linear regression of the enthalpy as a function of the number of reactant molecules. Non-linear regression and a Gaussian process model taking into account the compositional dependence of partial molar properties further support that multiple linear regression captures the correct physics.

@en

We introduce an alternative method to perform free energy calculations for mixtures at multiple temperatures and pressures from a single simulation, by combining umbrella sampling and the continuous fractional component Monte Carlo method. One can perform a simulation of a mixture at a certain pressure and temperature and accurately compute the chemical potential at other pressures and temperatures close to the simulation conditions. This method has the following advantages: (1) Accurate estimates of the chemical potential as a function of pressure and temperature are obtained from a single state simulation without additional postprocessing. This can potentially reduce the number of simulations of a system for free energy calculations for a specific temperature and/or pressure range. (2) Partial molar volumes and enthalpies are obtained directly from the estimated chemical potentials. We tested our method for a Lennard-Jones system, aqueous mixtures of methanol at T = 298 K and P = 1 bar, and a mixture of ammonia, nitrogen, and hydrogen at T = 573 K and P = 800 bar. For pure methanol (N = 410 molecules), we observed that the estimated chemical potentials from umbrella sampling are in excellent agreement with the reference values obtained from independent simulations, for ΔT = ±15 K and ΔP = 100 bar (with respect to the simulated system). For larger systems, this range becomes smaller because the relative fluctuations of energy and volume become smaller. Without sufficient overlap, the performance of the method will become poor especially for nonlinear variations of the chemical potential.

@en

We systematically study how the degree of framework flexibility affects the adsorption and diffusion of aromatics in MFI-type zeolites as computed by Monte Carlo simulations. It is observed that as the framework is more flexible, the zeolite structure is inherently changed. We have found that framework flexibility has a significant effect on the adsorption of aromatics in MFI-type zeolites, especially at high pressure. Framework flexibility allows the zeolite framework to accommodate to the presence of guest aromatic molecules. For very flexible zeolite frameworks, loadings up to two times larger than that in a rigid zeolite framework are obtained at a given pressure. We assessed the "flexible snapshot"method, which captures framework flexibility using independent snapshots of the framework. We have found that this method only works well when the loadings are low. This suggests that the effect of the guest molecules on the zeolite framework is important. Framework flexibility lowers the free-energy barriers between low energy states, increasing the rate of diffusion of aromatics in the straight channel of MFI-type zeolites for many orders of magnitude compared to a rigid zeolite framework. The simulations show that framework flexibility should not be neglected and that it significantly affects the diffusion and adsorption properties of aromatics in an MFI-type zeolite.

@en

The CFCMC simulation methodology considers an expanded ensemble to solve the problem of low insertion/deletion acceptance probabilities in open ensembles. It allows for a direct calculation of the chemical potential by binning of the coupling parameter λ and using the probabilities p(λ = 0) and p(λ = 1), which require extrapolation. Here, we show that this extrapolation leads to systematic errors when the distribution p(λ) is steep. We propose an alternative binning scheme which improves the accuracy of computed chemical potentials. We also investigate the use of multiple fractional molecules needed in simulations of multiple components, and show that these fractional molecules are very weakly correlated and that calculations of chemical potentials are not affected. The statistics of Boltzmann averages in systems with multiple fractional molecules is shown to be poor. Good agreement is found between CFCMC averages (uncorrected for the bias) and Boltzmann averages when the number of fractional molecules is less than 1% of the total number of all molecules. We found that, in dense systems, biased averages have a smaller uncertainty compared to Boltzmann averages.@en

The combination of the TraPPE and OPLS/2016 force fields with five water models, TIP3P, SPC/E, OPC, TIP4P/2005 and TIP4P/EW was used to compute mixing enthalpies, excess chemical potentials and activity coefficients of water and methanol. Excess chemical potentials and activity coefficients were computed in an expanded version of the NPT ensemble. We found the best agreement between experimental data for all the computed properties of water–methanol mixtures for the TIP4P/2005-TraPPE force fields. The performance of the spherical cutoff methods in MC and MD simulations was compared to the Ewald summation. The radial distribution functions obtained from the Ewald summation and the Damped-Shifted Force (DSF) method were in excellent agreement. Numerical artifacts appeared at the cutoff radius when the original Wolf method was used to calculate the electrostatic interactions. The calculated excess mixing enthalpies, excess chemical potentials, and activity coefficients of water and methanol obtained from the Wolf method were in good agreement with the DSF method. Our simulation results show that the numerical artifacts of the original Wolf method have little effect for energy calculations in aqueous methanol mixtures.@en

A mathematical model is used to predict adsorption isotherms from experimentally measured breakthrough curves. Using this approach, by performing only breakthrough experiments for a mixture of two (or more) components, one can obtain pure component adsorption isotherms up to the pressure of the experiment. As a case study, the adsorption of an equimolar mixture of CO₂ and CH₄ in zeolite ITQ-29 is investigated. Pure component linear adsorption isotherms for CO₂ and CH₄ are predicted by fitting the theoretical breakthrough curves to the experimental ones. Henry coefficients obtained from our approach are in excellent agreement with those measured experimentally. A similar procedure is applied to predict the complete Langmuir adsorption isotherm from breakthrough curves at high pressures. The resulting adsorption isotherms are in very good agreement with those measured experimentally. In our model for transient adsorption, mass transfer from the gas phase to the adsorbed phase is considered using the Linear Driving Force model and dispersion of the gas phase in the packed bed is taken into account. IAST is used to compute the equilibrium loadings for a mixture of gases. The influence of the dispersion coefficient and the effective mass transfer coefficient on the shape of breakthrough curves is investigated and discussed. Rough estimations of these values are sufficient to predict adsorption isotherms from breakthrough curves.

@en

A new macOS software package, iRASPA, for visualisation and editing of materials is presented. iRASPA is a document-based app that manages multiple documents with each document containing a unique set of data that is stored in a file located either in the application sandbox or in iCloud drive. The latter allows collaboration on a shared document (on High Sierra). A document contains a gallery of projects that show off the main features, a CloudKit-based access to the CoRE MOF database (approximately 8000 structures), and local projects of the user. Each project contains a scene of one or more structures that can initially be read from CIF, PDB or XYZ-files, or made from scratch. Main features of iRASPA are: structure creation and editing, pictures and movies, ambient occlusion and high-dynamic range rendering, collage of structures, (transparent) adsorption surfaces, cell replicas and supercells, symmetry operations like space group and primitive cell detection, screening of structures using user-defined predicates, and GPU-computation of helium void fraction and surface areas in a matter of seconds. Leveraging the latest graphics technologies like Metal, iRASPA can render hundreds of thousands of atoms (including ambient occlusion) with stunning performance. The software is freely available from the Mac App Store.

@en

Chemical potentials of coexisting gas and liquid phases for water, methanol, hydrogen sulphide and carbon dioxide for the temperature range (Formula presented.) K to (Formula presented.) K are computed using two different methodologies: (1) Widom’s test particle insertion (WTPI) method in the conventional Gibbs Ensemble (GE), and (2) the Continuous Fractional Component Gibbs Ensemble Monte Carlo (CFCGE MC) method. It is shown that the WTPI method fails to accurately compute the chemical potentials of water and methanol in the liquid phase at low temperatures, while accurate chemical potentials in the liquid phase are computed using CFCGE MC method. For the CFCGE MC method, the statistical uncertainty for computed chemical potentials of water and methanol in the liquid phase are considerably smaller compared to the WTPI method. For the water models considered in this study (SPC, TIP3P-EW, TIP4P-EW, TIP5P-EW), computed excess chemical potentials based on three-site models are in better agreement with the chemical potentials computed from an empirical equation of state from the NIST database. For water, orientational biasing is applied during test particle insertion to check whether certain orientations of test particle are energetically unfavourable. A two-dimensional Overlapping Distribution Method (ODM) in the NVT ensemble is derived for this purpose. It is shown that failure of the WTPI method for systems with a strong hydrogen bonding network does not depend on orientation of the test molecule in that system. For all systems in this study, the WTPI method breaks down when the void fraction of the system drops below approximately 0.50.@en

The effect of confinement on the equilibrium reactive system containing nitrogen dioxide and dinitrogen tetroxide is studied by molecular simulation and the reactive Monte Carlo (RxMC) approach. The bulk-phase reaction was successfully reproduced and five all-silica zeolites (i.e. FAU, FER, MFI, MOR, and TON) with different topologies were selected to study their adoption behavior. Dinitrogen tetroxide showed a stronger affinity than nitrogen dioxide in all the zeolites due to size effects, but exclusive adsorption sites in MOR allowed the adsorption of nitrogen dioxide with no competition at these sites. From the study of the adsorption isotherms and isobars of the reacting mixture, confinement enhanced the formation of dimers over the full range of pressure and temperature, finding the largest deviations from bulk fractions at low temperature and high pressure. The channel size and shape of the zeolite have a noticeable influence on the dinitrogen tetroxide formation, being more important in MFI, closely followed by TON and MOR, and finally FER and FAU. Preferential adsorption sites in MOR lead to an unusually strong selective adsorption towards nitrogen dioxide, demonstrating that the topological structure has a crucial influence on the composition of the mixture and must be carefully considered in systems containing nitrogen dioxide.

@en

Natural gas, synthesis gas, and flue gas typically contain a large number of impurities (e.g., acidic gases), which should be removed to avoid environmental and technological problems, and to meet customer specifications. One approach is to use physical solvents to remove the acidic gases. If no experimental data are available, the solubility data required for designing the sweetening process can be obtained from molecular simulations. Here, Monte Carlo (MC) simulations are used to compute the solubility of the gas molecules, i.e., carbonyl sulfide, carbon disulfide, sulfur dioxide, hydrogen sulfide, methyl mercaptan, carbon dioxide, and methane in the commercial solvents tetraethylene-glycol-dimethyl-ether (Selexol), n-methyl-2-pyrrolidone, propylene carbonate, methanol (Rectisol), and the ionic liquid 1-butyl-3-methylimidazolium bis(trifluoromethylsulfonyl)imide ([bmim][Tf₂N]). Henry coefficients of the gases in the investigated solvents are obtained from the computed solubilities. The ratio of Henry coefficients is used to compute ideal selectivities of the solvents. The solubilites and selectivities computed from MC simulations are compared with available experimental data. Some guidelines are provided to remove acidic gases using the investigated solvents. Rectisol is the best solvent for acid gas removal, but it should be used at low temperatures. Selexol and the ionic liquid have similar selectivity of sulfur compounds with respect to methane and may be used at elevated pressures and temperatures since both have low vapor pressures. The solubility of carbon disulfide, sulfur dioxide, and methyl mercaptan in these solvents is the highest. Hence, these components can be removed easily prior to hydrogen sulfide, carbonyl sulfide, and carbon dioxide in a pre-absorber.

@en

It is shown that ensemble averages computed in the Gibbs Ensemble with Continuous Fractional Component Monte Carlo (CFCMC GE) are different from those computed in the conventional Gibbs Ensemble (GE). However, it is possible to compute averages corresponding to the conventional GE while performing simulations in the CFCMC GE. In this way, one can benefit from the nice features of CFCMC GE (e.g. more efficient particle exchange) and at the same time compute the ensemble averages that correspond to the conventional GE. As a case study, the equilibrium pressure and densities of the systems of 256 and 512 LJ particles at different reduced temperatures ((Formula presented.)) are computed in the conventional GE and CFCMC GE. The validity of the expressions derived for computation of the thermodynamic pressure and densities corresponding to the conventional GE and computed in the CFCMC GE is examined numerically. The thermodynamic pressure in the conventional GE and CFCMC GE typically differs by at most 4%. It is shown that a very good estimate of the average pressure and densities corresponding to the conventional GE can be obtained by performing simulation in CFCMC GE and ignoring the contributions of the fractional molecule. It is also shown that the fractional molecule does not have an influence on the structure of the liquid, even for very small system sizes (e.g. 40 particles). The approach used here to compute the equilibrium pressure and densities of the conventional GE using the CFCMC GE can be easily extended to other thermodynamic properties and other ensembles.

@en

The influence of product shape selectivity on the bifunctional conversion of n-C₇ by zeolite catalysts is investigated. Three different zeolite catalysts with different pore sizes (MFI-type, MEL-type, and BEA-type zeolites) have been investigated experimentally. For all three catalysts, n-C₇ is isomerized to monobranched isomers which are further isomerized into dibranched isomers, and these dibranched molecules are converted into cracking products. More dibranched isomers and less cracking products are produced by BEA-type zeolite compared to MFI-type and MEL-type zeolites and clear differences are observed in the distribution of dibranched isomers produced by different catalysts. Molecular simulation is used to compute the adsorption isotherms and free energy barriers for diffusion of dibranched isomers in MFI-type, MEL-type, and BEA-type zeolites. Combining simulation results and experimental observations, it is shown that product shape selectivity can explain the distribution of dibranched molecules while transition state shape selectivity fails to do so. For the medium-pore zeolites (MFI-type and MEL-type zeolites), free energy barriers for diffusion of dibranched molecules are significant. For MFI-type and MEL-type zeolites, the dibranched molecule that has to overcome lower diffusion barrier is produced with a higher yield and the distribution of dimethylpentane molecules is determined by their diffusion rate. It is shown that there is almost no free energy barrier for the diffusion of any of these molecules in BEA-type zeolite. As BEA-type zeolite imposes no free energy barrier for diffusion of any of dibranched isomers, the distribution of dibranched isomers is very close to the equilibrium distribution in the gas phase. Due to the limited mobility of dimethylpentanes within the pores of MFI-type and MEL-type zeolites, most of the dimethylpentane molecules are trapped inside the zeolite and undergo consecutive cracking. Dimethylpentane molecules diffuse sufficiently fast in the large pores of BEA-type zeolite and transfer to the gas phase, before consecutive reaction converts these molecules into cracking products. Moreover, the effect of the MFI-type crystal size on the production of dibranched isomers is investigated. The yield of dibranched isomers reduces by increasing the size of the crystal and larger part of dibranched isomers are cracked as the crystal size of MFI-type is increased.

@en

An increase in demand for energy efficient processes for the separation of saturated and unsaturated light hydrocarbons mixtures drives the need of noncryogenic processes. The adsorptive separation using Metal-Organic Frameworks with coordinatively unsaturated metal sites may provide a cost-effective alternative due to the strong binding of the metal cation with the unsaturated hydrocarbons. Since experiments on adsorption equilibrium of gas mixtures are challenging, we propose classical force field based simulations to analyse the ability of MOF-74 with different metal substitutions for the separation of C2 and C3 olefin/paraffin binary mixtures. We parametrized the force field by fitting to available experimental single-component adsorption isotherms of ethane, ethene, propane, and propene in M–MOF-74 (M=Co, Fe, Mn, and Ni). The force field was validated for a variety of temperatures ranged from 273 K to 353 K. We then conducted Monte Carlo simulations in the Grand-Canonical ensemble to elucidate the adsorption mechanisms of the saturated/unsaturated hydrocarbon mixtures, at 318 K and 353 K. We computed the adsorption isotherms, and from these the adsorption selectivity, and addressed the variations of MOF properties with different metal cations. Fe-based MOF-74 appears the best option for both ethane/ethene and propane/propene separation applications. This finding partly agrees with previous work based on the Ideal Adsorbed Solution Theory.@en

A new formulation of the Reaction Ensemble Monte Carlo technique (RxMC) combined with the Continuous Fractional Component Monte Carlo method is presented. This method is denoted by serial Rx/CFC. The key ingredient is that fractional molecules of either reactants or reaction products are present and that chemical reactions always involve fractional molecules. Serial Rx/CFC has the following advantages compared to other approaches: (1) One directly obtains chemical potentials of all reactants and reaction products. Obtained chemical potentials can be used directly as an independent check to ensure that chemical equilibrium is achieved. (2) Independent biasing is applied to the fractional molecules of reactants and reaction products. Therefore, the efficiency of the algorithm is significantly increased, compared to the other approaches. (3) Changes in the maximum scaling parameter of intermolecular interactions can be chosen differently for reactants and reaction products. (4) The number of fractional molecules is reduced. As a proof of principle, our method is tested for Lennard-Jones systems at various pressures and for various chemical reactions. Excellent agreement was found both for average densities and equilibrium mixture compositions computed using serial Rx/CFC, RxMC/CFCMC previously introduced by Rosch and Maginn (Journal of Chemical Theory and Computation, 2011, 7, 269–279), and the conventional RxMC approach. The serial Rx/CFC approach is also tested for the reaction of ammonia synthesis at various temperatures and pressures. Excellent agreement was found between results obtained from serial Rx/CFC, experimental results from literature, and thermodynamic modeling using the Peng–Robinson equation of state. The efficiency of reaction trial moves is improved by a factor of 2 to 3 (depending on the system) compared to the RxMC/CFCMC formulation by Rosch and Maginn.@en

Computing bubble-points of multicomponent mixtures using Monte Carlo simulations is a non-trivial task. A new method is used to compute gas compositions from a known temperature, bubble-point pressure, and liquid composition. Monte Carlo simulations are used to calculate the bubble-points of carbon dioxide (CO2) and methane (CH4) mixtures in the ionic liquids (ILs) 1-n-butyl-3-methylimidazolium bis(trifluoromethylsulfonyl)imide [bmim][Tf2N] and 1-ethyl-3-methylimidazolium diethylphosphate [emim][dep]. The Continuous Fractional Component Monte Carlo (CFCMC) method in the osmotic ensemble has been used to compute the solubility of CO2/CH4 gas mixtures at different temperatures (T), pressures (P), and gas compositions (yi). The effect of T, P, and yi on the real CO2/CH4 selectivity (i.e., the selectivity of CO2 in the presence of CH4) is investigated. The real selectivity will differ from the ideal selectivity, which is defined as the ratio of the Henry's constants, if the solubility of CO2 is influenced by the presence of CH4. The computed real selectivities are compared with the experimentally obtained real and ideal selectivities. The real CO2/CH4 selectivity decreases with increasing temperature and pressure, while the gas phase composition has a minor effect. The real selectivity is approximately identical to the ideal selectivity for relatively low pressures and low solute concentrations in the liquid phase. The real selectivity deviates from the ideal selectivity as the solute concentration in the liquid phase increases.@en