Accurate Evaluation of the Dispersion Energy in the Simulation of Gas

Subscriber access provided by UB + Fachbibliothek Chemie | (FU-Bibliothekssystem)

Article

Accurate Evaluation of the Dispersion Energy in the Simulation of Gas Adsorption into Porous Zeolites. Alberto Fraccarollo, Lorenzo Canti, Leonardo Marchese, and Maurizio Cossi J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.6b01021 • Publication Date (Web): 27 Feb 2017 Downloaded from http://pubs.acs.org on March 1, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Theory and Computation is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Accurate Evaluation of the Dispersion Energy in the Simulation of Gas Adsorption into Porous Zeolites. Alberto Fraccarollo, Lorenzo Canti, Leonardo Marchese, and Maurizio Cossi∗ Dipartimento di Scienze e Innovazione Tecnologica (DISIT), Università del Piemonte Orientale, via T. Michel 11, Alessandria, Italy E-mail: [email protected]

Abstract The force fields used to simulate the gas adsorption in porous materials are strongly dominated by the van der Waals terms: here we discuss the delicate problem to estimate these terms accurately, analyzing the effect of different models. To this end, we simulated the physisorption of CH4 , CO2 and Ar into various Al-free microporous zeolites (ITQ-29, SSZ-13 and silicalite-1), comparing the theoretical results with accurate experimental isotherms. The van der Waals (vdW) terms in the force fields were parameterized against the free gas densities and high level quantum mechanical (QM) calculations, comparing different methods to evaluate the dispersion energies. In particular, MP2 and DFT with semiempirical corrections, with suitable basis sets, were chosen to approximate the best QM calculations; either Lennard-Jones or Morse expressions were used to include the vdW terms in the force fields. The comparison of the simulated and experimental isotherms revealed that a strong interplay exists between the definition of the dispersion energies and the functional form used in the force field; these results are fairly general and reproducible, at least for the systems considered here. On this basis, the reliability of different models can be discussed and a recipe can be provided to obtain accurate simulated adsorption isotherms.

1

ACS Paragon Plus Environment

1 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation Introduction

Nanoporous materials are of great interest for applications in molecular separation, 1–3 gas storage, 4–6 heterogeneous catalysis. 7–9 The performance of such materials in all the fields of application is critically affected by their porous structure, 10 i. e. the specific surface area, the micro and mesoporous volumes and the pore size distribution, as well as by their physical and chemical stability. Alongside many experimental techniques, computational modeling can be of invaluable help during the design of materials with tailored structure and properties, supporting the physico-chemical characterization, providing structure-properties relationships, and rationalizing the host-guest interactions. 11–16 Most of the modeling activities are based on Monte Carlo or molecular mechanics simulations, 17–24 and rely on classical force fields (FF), whose quality depends on the functional forms and the set of parameters defining the various energy terms. 25–29 This work is focused on the modeling of gas adsorption inside porous materials: this process is used either to characterize the porous structure and volume (in this case the adsorbate is typically N2 at 77 K or Ar at 87 K, though CO2 at room temperature is sometimes used too), or for gas storage at high pressure (involving methane, other hydrocarbons or hydrogen as energy carriers, or CO2 ), or even for gas separation at lower pressures. The adsorption isotherms are usually simulated by Grand Canonical Monte Carlo (GCMC) calculations, 30 and the results depend mostly on the host-guest non-bonded interactions. Even if in some cases the deformation of the material structure during the adsorption can be taken into account also, in the following we will consider rigid systems only. Then the relevant FF terms are van der Waals and in some cases (typically for CO2 ) the Coulomb charge-charge contributions. To obtain a greater accuracy, these terms are often fitted for the specific systems under study against experimental data or quantum mechanical (QM) calculations on model systems. Another possible approach is to derive the full FF from QM calculations: an efficient authomated protocol has been proposed for this purpose, 31 and another procedure has been described recently for condensed phase simulations. 27,32,33 Other suitable force fields able to describe liquids as well as gas phases have been derived from ab initio calculations, also including many-body and polarization effects. 34,35 2


Page 2 of 43

Page 3 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


When QM calculations are used to model the van der Waals (vdW) host-guest interactions, a delicate question arises about the level of theory needed to include the elusive dispersion energy, which is in this case the dominant contribution. In past decades, the dispersion energies were mostly estimated with post-HartreeFock methods, typically second order Möller-Plesset (MP2) or Coupled Clusters (CC). The latter, extended to single, double and (partially) triple excitations, i.e. CCSD(T), and applied with a large basis set (ideally a complete basis set, for systems small enough) is the state-of-the-art QM method for this calculation, 36 but it is also very expensive computationally. Quite recently, a different strategy has been proposed, based on semiempirical expressions which provide the dispersion energy as an additive contribution to DFT calculations, thus avoiding the heavy QM methods mentioned above. 37 These contributions rely on sets of atom-atom parameters, optimized for many density functionals: the most popular expression, proposed by S. Grimme, 38–41 is considered very reliable in the recent literature, 42–45 and is implemented in several QM codes. Using the tag proposed by Grimme and applied e.g. in Gaussian09 implementation, 46 we shall refer to this method as DFT-D3. Clearly, the use of atom-atom pair potentials cannot account for many-body interactions, yet affecting the physisorption processes: in this sense, effective pair potentials are needed, able to reproduce the experimental behaviour. It is interesting that such effective potentials can be optimized on the basis of pair interactions only and still reproduce the observed adsorption isotherms accurately; even more remarkably, the potentials are highly transferable to similar sorbates and sorbents, as shown below. On the other hand, the same approach based on pair potentials is adopted in many computational studies of physisorption processes, 47–54 though more complex models, with different potentials including many-body polarizations and inhomogeneity effects, have also been described and used recently. 55–59 In the present paper we compare the MP2 and the DFT-D3 calculation of dispersion energies, to obtain reliable force fields for the gas adsorption modeling; the interplay of the chosen method and the functional form of the non-bonded term (i.e. Morse or Lennard-Jones functions) is also carefully investigated. Our goal is twofold: i) provide an insight on the most reliable and convenient strat3



egy to parameterize a FF for this kind of application, and ii) discuss to what extent 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the different definitions of van der Waals interactions are related to experimentally measurable data. The paper is organized as follows. First we consider in detail the adsorption of CH4 , CO2 and Ar into the microporous all-silica zeolite ITQ-29, recently synthesized and characterized: 60,61 the crystalline structure of this solid, along with the adsorption isotherms of the mentioned gas have been accurately determined, and they were used as the experimental benchmarks in this study. Since this material is highly crystalline, and the X-Ray structure was used in the simulations, we assume that the agreement or disagreement between theoretical and experimental adsorptions is directly related to the FF quality. Then, to test the transferability of our results, the same FF are used to simulate the adsorption of CH4 and CO2 in other highly crystalline zeolites, i.e. SSZ-13 and silicalite-1.

2

Methods

QM calculations were performed with Gaussian09 package, 46 at DFT, MP2 and CCSD(T) 62,63 levels. In DFT calculations the hybrid functional B3LYP was used, supplemented by the semiempirical dispersion correction proposed by Grimme with the D3 set of parameters 40 (referred to as DFT-D3 in the following). The following Dunning’s correlation consistent basis sets, with double-, triple- and quadruple-ζ valence shells and possibly added diffuse functions, were used throughout the paper: cc-pVDZ (DZ), aug-cc-pVDZ (aug-DZ), cc-pVTZ (TZ), cc-pVQZ (QZ). All the QM interaction energies have been corrected for the basis set superposition error (BSSE) with counterpoise (CP) method. 64 Classical energies were computed with purposely modified force fields: the bonding parameters were taken from Dreiding FF 13,65 without changes, while for vdW nonbonded interactions we used either Lennard-Jones (LJ) 6-12,

E(r) = D0

"

R0 r

12

−2

R0 r

4


6 #

(1)

Page 4 of 43

Page 5 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

or Morse,


r y E(r) = D0 X − 2X ; X = exp 1− 2 R0

2

(2)

equations, where r is the atom-atom distance, R0 and D0 the equilibrium distance and the corresponding energy, and y a dimensionless factor determining the stiffness of the potential. These parameters were optimized as detailed below, and their values are collected in the Supporting Information (SI) file. For CO2 the electrostatic interaction term was also included, with the atomic charges proposed in ref. 21: +0.48222 for C and −0.24111 for O; since no atomic charges were assigned to zeolite atoms, the electrostatic term affects only the gas-gas interactions. Molecular mechanics calculations of single point energies were performed with the Forcite module provided by Materials Studio 6.0 package; the Sorption module of the same package was used to perform GCMC simulations of the adsorption isotherms. GCMC calculations included 106 equilibration and 5 × 106 production steps; Ewald summation method was used for electrostatic contribuiotns, while for the other non-bonded terms the default value of the cutoff for vdW interactions (15.5 Å) was used when not otherwise stated. Three types of moves were applied in GCMC calculations: i) translation and rotation steps, handled with the usual Metropolis algorithm; 66 ii) insertion with random orientation; iii) removal. The probabilities were 0.2 for translations and rotations, 0.02 for insertion and deletion moves. Gas fugacities (directly related to the GCMC chemical potential) were converted to pressures by means of the suitable fugacity coefficients, which were computed for all the gas of interest as described e.g. in ref. 50. It is important to stress that in this approach the fugacity coefficients are obtained from the experimental compression factors, thus assuming that the pure gas simulations match the experimental densities well over the whole pressure range, as discussed in detail below; note however that the distinction between fugacity and pressure is relevant for rather high pressures only. The crystallographic structure (CIF) of the all-silica zeolite ITQ-29 was taken from ref. 60, all the other structures were found in the database of the International Zeolite Association. The experimental adsorption isotherms of Ar, CH4 and CO2 in ITQ-29 were kindly provided by prof. Fernando Rey of the Universitat Politècnica de València, CO2 adsorption isotherms in SSZ-13 and silicalite-1 were taken from 5



refs. 55 and 67, respectively, and that of CH4 in silicalite-1 from ref. 68. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

As detailed in the following, various FF were fitted against QM calculations involving model fragments of ITQ-29, cut out of its periodic structure: two such fragments were defined, with stoichiometries Si8 O8 H16 (fragment I) and Si12 O15 H20 (fragment II). Fragment I features a 4-T ring, each silicon atom forming a further Si-O-Si bridge, while II comprises a 8-T ring fused to two 4-T rings, as shown in Figure 1; all the Si free valencies were saturated by hydrogen atoms. Figure 1

3

Results and discussion

3.1

Gas-gas interactions

The parameters in equations 1 and 2 were optimized for Ar, C, H and O atoms in order to fit the experimental densities 69 of argon at 87 K, methane at 283 K and carbon dioxide at 298 K (the corresponding Dreiding atom types are Ar, C_3, C_1, H_ and O_2). We choosed to fit the FF against the experimental densities instead of QM calculations to concentrate as much as possible on the host-guest interactions, as detailed in the following. The simulated and experimental results are compared in Figure 2: both LJ and Morse functions provide gas densities in excellent agreement with the experiment in all the considered range of pressures. Figure 2 However, the FF validation for Ar and CH4 could be questioned as in this range of temperature and pressure the experimental curves are too close to the ideal gas behavior, calling for a further check of the fitted parameters. Then we extended the pressure range in the simulations: at 87 K, argon condensates at 0.99 bar and this phase transition was correctly reproduced by both FF for pressure between 1 and 1.05 bar, as shown in Figure S1-a in the SI file. To test the sensitivity of the simulated curves to the size of FF parameters, we repeated the calculation with parameters increased by 20% and reduced by 5%, with the results reported in Figure 6


Page 6 of 43

Page 7 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


S1-b,c: clearly, the quality of the results is strongly dependent on the parameter size. As for methane at 283 K, the simulated density is in excellent agreement with the experiment up to 200 bar, with both FF; also in this case, increasing the parameters by 20% leads to completely wrong simulated curves, as illustrated by Figure S2 in the SI file. We conclude that the fitted parameters are significant and adequate to describe gas-gas interactions also outside the ideal gas region.

3.2

QM benchmarks

A preliminary test was performed about the reliability of the various QM methods for the calculation of host-guest interactions: the addition energies of one molecule of CH4 and one of CO2 with a minimal silica model, SiH3 − O − SiH3 , were scanned with respect to the intermolecular distance at B3LYP, B3LYP-D3 and MP2 levels with DZ basis set, and compared with the results of CCSD(T) calculations with DZ, TZ and QZ bases. As shown in Figure 3, the CC calculation is quite sensitive to the quality of the basis set: among the variational methods used here, the CCSD(T)/QZ is the most reliable one, 36 and it was taken as a benchmark. The simple DFT interaction energy is way too small both for methane and for carbon dioxide: on the other hand, almost all the correlation energy (one could say, in this context, the dispersion energy) is recovered by DFT-D3, which agrees very well with the benchmark. MP2 performs better than DFT, but the energies at this level are still far from the benchmark. Figure 3 The basis set dependence of MP2 and DFT-D3 energies, compared to the CC/QZ reference, is illustrated in Figure 4. DFT-D3 is affected only slightly by the basis set, the DZ result being anyway in the best agreement with the benchmark; MP2 energies are much more sensitive to the basis set, the agreement with the reference improving markedly with TZ and even more with aug-DZ. The deviations of all the addition energies with respect to the QM benchmark are also illustrated in Figure S1 in the Supporting Information. All the DFT-D3 energies are slightly too large, in contrast with the other methods, and the smaller basis set (DZ) is in better agreement with the reference, especially with CO2 , where the error is very small; on 7



the other hand, MP2 results improve with larger basis sets, obtaining with aug-DZ 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

a very small deviation for CH4 . Interestinlgy, both DFT-D3 and MP2 with the suitable basis sets appear more reliable than CC/DZ and better or comparable to CC/TZ and CC/aug-DZ, though much less expensive computationally. Based on the above results, in the following we will use either DFT-D3/DZ or MP2/aug-DZ methods to compute QM host-guest energies: for methane, the two methods approximate equally well the reference curve, while for CO2 DFT-D3/DZ is quite closer to CC/QZ. Figure 4

3.3

CH4 adsorption in ITQ-29

To simulate CH4 adsorption into ITQ-29 the first step was the parameterization of the vdW terms for Si3 and O_3 atom types, keeping the same parameters for methane atoms as in the gas phase. To do that, we considered one of the zeolite models described above (fragment I, see Figure 1): the addition energy of one methane molecule with this cluster was scanned along a trajectory “perpendicular” to the fragment, as depicted in Figure S2 in the SI, at MP2/aug-DZ and DFT-D3/DZ levels, and these energy scans were fitted with LJ and Morse functions. It is noteworthy that the fitting procedures tend to weigh the points around the potential mimina more than the interactions at large distance, which become rapidly very small: however, as discussed in detail below, in some cases the cumulative effect of long-range interactions is not negligible. On the other hand, the accurate modeling of physisorption processes relies on the additive effects of short and long range interactions, so that the whole potential function has to be carefully optimized: we will investigate what fitting strategy takes this point in better account. Four FF were thus obtained, which will be indicated in the following as LJ-MP2, LJ-D3, Mo-MP2 and Mo-D3, in obvious notation: all the fitted parameters are listed in the SI. The QM and classical energy scans are shown in Figure 5: it is worth noting that the DFT-D3 curve, as well as the FF fitted on it, are markedly more attractive than the MP2 counterparts, as expected after the results of the 8


Page 8 of 43

Page 9 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


previous section. Moreover, the LJ function fits the MP2 scan better than Morse in the medium distance range (approximately 4-8 Å), while both expressions perform equally well with the DFT-D3 curve: MP2 interaction energy, though less attractive in the minimum, rises less steeply than DFT-D3 and LJ function seems to be flexible enough to follow both trends, unlike Morse function. Figure 5 The FF so obtained were used in GCMC calculations to simulate the adsorption isotherms of CH4 in ITQ-29 at 283 K, which are compared to the experimental isotherm in Figure 6. Figure 6 Quite unexpectedly, LJ-MP2 and Mo-D3 isotherms are very similar and in excellent agreement with the experiment, despite the marked difference in the well depth of their interaction curves (Figure 5). On the other hand, LJ-D3 predicts a strikingly large gas uptake, while Mo-MP2 strongly underestimates it. Figure S3 in SI reports LJ-MP2 and Mo-D3 curves only, along with the experimental isotherm, to better evaluate their agreement in a tighter scale. Apparently, when LJ function is used for the vdW terms, a less attractive potential (at least around the interaction energy minimum) is needed than with Morse expression: to obtain a simulated isotherm very close to the experiment, the FF can be fitted either to MP2 or to DFT-D3 interaction curves, provided the suitable functional form is used. A related behaviour was observed in ref. 25 for methane adsorbed on carbon nanotubes: the GCMC simulations performed with LJ and Morse functions provided quite different gas uptakes, though both FF had been fitted on the same MP2 energy scans. In that work, it was also pointed out that simple DFT scans cannot be used for the FF fitting, since most of the dispersion energy is missing, and also that the size of the basis set in MP2 calculations is crucial in determining which fit agrees better with the benchmark. It is worth analyzing in greater detail how two FF fitted to quite different potential curves, like MP2 and DFT-D3 in Figure 5, can provide such similar GCMC simulated

9



isotherms. In fact, the host-guest interaction energy in the model fragment is much 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

smaller in LJ-MP2 than in Mo-D3: at some point, passing from the model fragment to the actual periodic crystal, this difference must be compensated. To investigate this matter, a configuration was extracted from one of the GCMC simulations, comprising one CH4 molecule in the ITQ-29 unit cell, and the single point energy of this periodic structure was computed with Forcite: the LJ-MP2 and Mo-D3 addition energies (referred to isolated methane and empty unit cell) are −2.069 and −1.950 kcal/mol, respectively. Note that in this particular arrangement LJ-MP2 provides an energy slightly larger in absolute value, but of course the opposite occurs in some of the other configurations which contribute to the GCMC simulation. Then a finite cluster was cut out from this configuration, with the same stoichiometry as model fragment I, leaving the methane molecule in its position (this is called Model 1 in the following). Next, Model 1 was “inflated” adding more and more zeolite atoms, until the whole unit cell was rebuilt (these are the Models 2 to 5); then larger models were formed, placing the methane molecule in n × n × n zeolite supercells, with n = 2, 3, 4, 6, corresponding to Models 6 to 9. All the Models defined here are depicted in Figure S4 in the SI. The formation energies of all the Models were computed with single point (nonperiodic) Forcite calculations, using LJ-MP2 and Mo-D3 (the energies are referred to the isolated CH4 and zeolite fragments): we guess that passing from Model 1 to Model 9, all the relevant host-guest interactions affecting the methane molecule are included. To obtain the periodic formation energies we also have to consider guest-guest interactions: they were estimated by computing the formation energies of n × n × n “supercells” formed only by the central methane and its n3 − 1 periodic images (in this case the reference is the isolated CH4 plus the whole set of images); these supercells are shown in Figure S5 in the SI. The results are illustrated in Figures 7 (host-guest) and 8 (guest-guest contribution). Figures 7, 8 Clearly, the host-guest contribution is largely dominant: the guest-guest energies are indeed negligible, noting the different y scales in Figures 7 and 8. This does not exclude, of course, that guest-guest interactions play a significant role at higher 10


Page 10 of 43

Page 11 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


pressures, when more than one gas molecule are comprised in the unit cell. As expected, for Model 1 Mo-D3 (based on a more attractive potential curve) provides a formation energy quite larger, in absolute value, than LJ-MP2; as the Model size increases, the formation energies become more and more negative (since more atoms are added to the clusters), but keeping their relative order, until the complete unit cell (Model 5) is formed. On the other hand, when supercells are formed Mo-D3 energies change very little (with 2 × 2 × 2 cell, Model 6, the periodic energy is already practically reached), while with LJ-MP2 they keep growing in absolute value, eventually overpassing Mo-D3 values and converging to the periodic energy at 6 × 6 × 6 cell (Model 9). It can be added that LJ-MP2 guest-guest energies (Figure 8), though negligible compared to the previous contribution, are always more negative than Mo-D3 values, which are close to zero for all the supercells. These results show that LJ-MP2 is much more sensitive to long range interactions than its competitor: as far as the atoms in the unit cell are concerned, Mo-D3 provides more negative formation energies, but when farther atoms are added in the supercells this FF quickly gets to the asymptotic value. Such atoms, on the contrary, keep improving the LJ-MP2 energies: their contributions, though very small individually, sum up to a significant fraction of the periodic value. Note that, since distant atoms had to be accounted for, it was necessary to increase the vdW energy cufoff up to 60 Å (the Forcite default value being 15.5 Å). It has already been observed that Morse function leads to “harder” potentials than LJ: 70,71 when the two functions are fitted to reproduce the same minimum energy, the latter always provides more negative long range interactions, even if this effect is not dramatic usually, unless a large number of atoms contribute to the interaction as in our case. On the other hand, Han, Goddard III and coworkers 21,29 claimed that Morse function is preferable over LJ for the simulation of CO2 adsorption in metal-organic frameworks (MOF) and CH4 in covalent frameworks (COF), since LJ potentials are judged too stiff in the short range region and long range interactions, where the 1/r6 character is important, are considered not relevant. The results discussed above, however, show that in the present case the inclusion of long range interactions can change the order of LJ and Morse energies and modify substantially the simulated isotherms. In fact we find here that DFT-D3/DZ energies are attractive enough to compensate the smaller interaction range of Morse function, but too 11



attractive when LJ is used (see LJ-D3 isotherm in Figure 6); the opposite is true for 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

MP2/aug-DZ interaction energies. The results discussed above shed some light on the wisest strategy to build an effective atom-atom potential able to reproduce the experimental behaviour, suggesting that the interplay between the model function interaction range and the depth of the potential curves plays an important role.

3.4

CO2 adsorption in ITQ-29

The same approach described above for methane was applied to simulate carbon dioxide adsorption: a cluster was defined, with a CO2 molecule moving above fragment I (see Figure S6-a in the SI), and LJ and Morse parameters for the zeolite atoms were fitted on DFT-D3/DZ and MP2/aug-DZ energy scans. The vdW parameters and atomic charges for CO2 atoms were left as determined previously for gas-gas interactions; the optimized parameters for all the FF are reported in the SI. The QM and the classical energy scans are compared in Figure 9: similarly to what found for methane, MP2 provides interaction energies markedly smaller than DFTD3, the well minimum being around 2 and 3.4 kcal/mol, respectively. The four FF behave accordingly: also in this case, Morse curves tend to rise more steeply than LJ, so that long range interactions are likely to be more important in FF containing LJ terms. Figure 9 To check the sensitivity of the fitted parameters to the orientation of the CO2 molecule with respect to the zeolite fragment, we repeated the QM and FF calculations also along the trajectory drawn in Figure S6-b (SI): as shown in Figure S7 (SI) all the FF fitted above reproduce very well the QM calculations also for this orientation. Then the adsorption isotherms were simulated with GCMC, and the FF results compared to the experimental isotherm in Figure 10. LJ-MP2 and Mo-D3 are again in very good agreement with the experimental findings, even if in this case the isotherms are slightly less overlapped than for methane; also for CO2 LJ-D3 adsorption is strongly overestimated, while Mo-MP2 uptakes are way too low. Then the same conclusion as for methane holds for CO2 too: the FF can be parameterized 12


Page 12 of 43

Page 13 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


either on DFT-D3 or on MP2 curves, if the proper functional form is selected for the vdW terms. Figure 10

3.5

Ar adsorption and pore size analysis

As noted above, the adsorption of argon at 87 K is often used to characterize the fine structure of microporous solids 72–74 (we recall that, according to the IUPAC definition, 75 micro- and mesopores have diameters below 2 nm or comprised between 2 and 50 nm, respectively). An accurate description of the pore size distribution (PSD), as well as of the cumulative pore volume (CPV), can be obtained by analyzing the adsorption isotherms with computational tools recently developed and based on DFT parameterization. 74,76 The same tools can also be applied to GCMC simulated isotherms, treating them as their experimental counterparts: so the PSD of the model and the laboratory material can be directly compared. We resorted to the procedure described above to optimize the vdW parameters for the zeolite atoms, fitting the QM energy scans of Ar on the same fragment used above (values reported in the SI): on the basis of the previous sections findings, we optimized LJ-MP2 and Mo-D3 FF only. The corresponding adsorption isotherms are reported in Figure 11 along with the experimental isotherm: at such a low temperature the zeolite pores are saturated very quickly, so that the curves are better compared in the logarithmic plot reported as inset in Figure 11. Figure 11 The PSD and CPV obtained from the experimental and simulated isotherms are illustrated in Figure 12. Despite the difficulty of reproducing the very steep adsorption at extremely low pressures, the computed pore structure is in remarkable agreement with the experimental one: all the methods find a narrow family of micropores around 8.5 − 9 Å, though the two FF slightly overstimate the number of pores. The CPV are in good agreement also (Figure 12, inset). Figure 12

13


3.6 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Effect of the model size

In the previous sections we have seen that LJ-MP2 and Mo-D3 are able to reproduce the adsorption of CH4 , CO2 and Ar in ITQ-29 very well, and we have discussed how the short- and long range interactions balance in the two FF, so that different interaction curves provide very similar GCMC adsorption isotherms. A question still remains open, whether the size of the fragment used to model the zeolite walls in the FF parameterization affects the quality of the results. To clarify this point, we resorted to a larger model (fragment II, described in section “Methods”) and computed its interaction energy with CH4 along a perpendicular path as done with fragment I: the adduct is depicted in Figure S8 in the SI. Then the vdW parameters were re-optimized with this model: LJ and Morse functions were fitted against MP2/aug-DZ and DFT-D3/DZ calculations, respectively, obtaining two new FF indicated as LJ-MP2-II and Mo-D3-II. The energy scans at QM level are reported in Figure 13, along with the Forcite energies obtained with LJ-MP2-II and Mo-D3-II; we also report for comparison the energies from LJ-MP2 and Mo-D3. Note that fragment II contains a large 8membered ring, which can be crossed by the methane molecule (in fact, the gas diffuses into the zeolite through these windows): in this case then the scan begins at distance zero, with the molecule in the same plane as the ring atoms (see Figure S8, SI). Figure 13 It is clear that LJ-MP2 reproduces also the scan with fragment II quite well, and the Forcite curves obtained with LJ-MP2 and LJ-MP2-II are very close to each other (the latter being slightly more accurate for short distances); on the other hand, Mo-D3 reproduces the QM curve quite poorly, with a minimum at higher energy and markedly shifted: as a consequence, the Mo-D3 and Mo-D3-II profiles are very different, the latter being definitely more attractive around the minimum. This different behaviour reflects on the GCMC simulated isotherms, as shown in Figure 14: LJ-MP2-II performs as its homologue parameterized on fragment I, reproducing the experimental isotherm very well, whereas Mo-D3-II overestimates the gas uptake largely. The conclusion is that LJ parameters optimized on the MP2 14


Page 14 of 43

Page 15 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


scan are satisfactorily independent on the cluster size, and thus more reliable than Morse parameters fitted on DFT-D3: apparently, the good performance of the latter FF in the simulations reported above is due to the favorable choice of the model fragment. Figure 14 One could wonder if the failure of Mo-D3 to reproduce the scan on fragment II (leading to the new parameterization, Mo-D3-II, which proved to be too attractive) is due to the rigidity of the Morse expression or rather to the use of the semiempirical D3 parameters. To get some insight on this question, we cross-checked the functions and the potentials, scanning the energies on fragment II also with Mo-MP2 and LJ-D3 (i.e. the other FF obtained above for fragment I), with the results shown in Figures S9 and S10 in the SI. We find that both LJ-MP2 and Mo-MP2 fit the fragment II MP2 curve quite well, reproducing the well depth and the position of the minimum: in the mediumrange region, as observed above with fragment I, Morse function rises too steeply, and we know that this feature ultimately causes Mo-MP2 to underestimate the gas uptake. Anyway, if we had chosed fragment II to perform the initial FF fitting, we would have obtained the same, or very similar parameters both for LJ-MP2 and for Mo-MP2. On the other hand, LJ-D3 curve performs better than Mo-D3, since it reproduces the DFT-D3 well depth in fragment II, but the minimum distance is still too large: in general, also LJ-D3 curve seems to be too loose compared to QM; using fragment II, we would likely obtain more attractive parameters also for LJ-D3. Then, the qualitative conclusion is that Morse function confirms to be too rigid, compared to LJ, but at least a part of the non-transferability of the FF passing from fragment I to II is due to the use of DFT-D3 reference data.

3.7

Gas adsorption in other zeolites

In the previous sections the adsorption of various gases in ITQ-29 was discussed in detail, and the different behavior of the FF was analyzed depending on the vdW 15



term functional form and the QM fitting data. One could wonder to what extent 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

these results depend on the particular system considered: then we simulated the adsorption of CH4 and CO2 in other all-Si crystalline zeolites, to investigate the effect of different structures and pore sizes. In all the simulations the FF described above were used without changes: since the building units of all the host systems are the same, we expect that the atom-atom interactions are reproduced reasonably well also in these zeolites (a discussion on the transferability of force fields among different zeolite frameworks can be found e.g. in ref. 77). We considered the following host-guest systems: methane 68,78 at 303 K and carbon dioxide 67,78 at 308 K in silicalite-1 (an all-Si MFI-type zeolite with a 7 Å wide channel network), and carbon dioxide 55 at 303 K in SSZ-13 (a chabazite-like framework with pores of about 8 Å width). The simulated isotherms are compared to the experimental data in Figure 15; note that for CO2 @SSZ-13 the pressure range is limited to 1 bar. Figure 15 One can see that the trends discussed above for ITQ-29 are reproduced very well in all the other systems also: the isotherms simulated with LJ-MP2 and Mo-D3 are very similar to each other and in good agreement with the experiment, while MoMP2 severely underestimates and LJ-D3 overestimates the gas uptake. For CO2 in SSZ-13 the Monte Carlo results are less close to the experimental isotherm, though the trend of the different FF is still very clear: note that ref. 55, besides providing the measured uptakes, also presents a simulation with a purposely optimized FF in better agreement with the experiment. However we remind that our purpose is not to parameterize an optimal FF for all the systems, bur rather to investigate the effect of the dispersion energy models. From this point of view, all these simulations confirm what discussed above for gases in ITQ-29: for all the considered adsorbates, also in different zeolites, there is a strong correlation between the functional form of vdW terms and the QM level used to parameterize them, and satisfactory results can be obtained from MP2 and DFTD3 benchmarks only if LJ or Morse functions are used, respectively. Furthermore, we tested the two FF optimized on the larger fragment, i.e. LJ-MP2-II and MoD3-II, by repeating the simulation for CH4 in silicalite-1 at 303 K. The resulting 16


Page 16 of 43

Page 17 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


isotherms are illustrated in Figure S11 in the SI file: as already found for ITQ-29, also in this zeolite LJ-MP2-II performs as well as its homologue optimized on the small fragment, while the gas uptake is strongly overestimated by Mo-D3-II. It is also interesting to compare the present results with the simulations recently obtained for CO2 in ITQ-29 and silicalite-1 with an empirical Lennard-Jones potential, 79 optimized against experimental data 61,68 (and combined with electrostatic terms based on DFT-derived atomic charges). Clearly the simulated adsorptions in ref. 79 are in excellent agreement with the experiment, since the FF parameters were fitted to these data, but they also agree very well with our results: for instance, the CO2 uptake at 6 bar in ITQ-29 is 97.6 cm3 /g for the empirical potential and 98.0 and 95.1 cm3 /g for our LJ-MP2 and Mo-D3 FF, respectively; in silicalite-1 at 6 bar the empirical potential uptake is 59.0 cm3 /g, compared to 59.0 cm3 /g for LJ-MP2 and 60.2 cm3 /g for Mo-D3. Then LJ-MP2 and Mo-D3, though obtained by fitting the quantum mechanical pair interaction energies, provide results fully comparable to those obtained by a FF directly optimized on the experimental data, confirming their reliability as effective pair potentials.

4

Conclusions

We have discussed how the dispersion/repulsion interactions should be accounted for in the simulation of gas adsorption in porous materials (various all-Si zeolites); in particular, what QM level is most suited for the fitting of the vdW terms in the force fields. On the basis of high-level CCSD(T) benchmark calculations, two QM levels were tested for the FF fitting: DFT with semiempirical corrections (D3) and MP2, with suitable basis sets. Then we proposed an (indirect) experimental validation: the FF fitted to the different QM models are used to simulate the adsorption isotherms of various gases in zeolites of known structure. Comparing the simulation results to accurate experimental isotherms obtained for the same systems, we could evaluate the fitness of the FF and the reliability of the underlying dispersion/repulsion models. The case of methane adsorbed in ITQ-29 has been analyzed in great detail, also 17



to clarify how the interactions energies evaluated on small model systems, and used 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

to fit the vdW terms, sum up to produce the adsorption energies in the periodic frameworks. Then the same approach was used to describe the adsorption of CO2 and Ar in ITQ-29; the results were generalized, applying the same FF to simulate the adsorption of CH4 , and CO2 in other zeolites with different framework structures (i.e. silicalite-1 and SSZ-13). A strong relationship was found between the QM level at which the dispersion/repulsion energy is evaluated and the functional form (namely, Lennard-Jones or Morse) used to include this term in the FF: the best results, using a small fragment model of the zeolite, are obtained either with MP2/aug-cc-pVDZ or with B3LYP/cc-pVDZ plus Grimme’s D3 correction, provided that LJ is used in the former and Morse in the latter case. The other combinations lead to strong over- or underestimations of the gas uptake. Besides providing a recipe for a reliable fitting of the vdW terms, thus obtaining useful effective potentials for the modeling of the sorption phenomena, such a correlation is interesting to discuss the physical content of the dispersion energy models. In fact, DFT-D3 method usually provides more attractive potential curves, which are in good agreement with the CCSD(T) results: however, the resulting vdW terms are too attractive if a 1/r6 trend is used for the long range interactions, as in LJ expression. On the other hand, fitting LJ terms on the less attractive MP2 potential allows to reproduce the experimental adsorption; Morse function, rising more steeply with the interatomic distance and thus less sensitive to distant atoms, provides good results when fitted on the DFT-D3 curve. If a larger zeolite model is used in the fitting, however, only the LJ function optimized on MP2 curves keeps the good agreement with the experiment, while the Morse function fitted on DFT-D3 leads to overestimate the gas adsorption strongly: then it seems that the good performance of Mo-D3 with the smaller fragment was fortuitous, and that DFT-D3 potential tends to be generally too attractive, at least for the systems examined here.

18


Page 18 of 43

Page 19 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

of Chemical Theory and Computation SupportingJournal Information

The Supporting Information is available free of charge the Internet at http://pubs.acs.org. Optimized FF parameters; paths for molecules/zeolite fragments scans; models used in Forcite calculations; interaction energy curves for methane and large zeolite fragment; adsorption isotherms of methane in silicalite-1 (with LJ-MP2-II and Mo-D3-II).

Acknowledgements Prof. Fernando Rey (Universitat Politècnica de València) is gratefully thanked for providing the adsorption isotherms of methane, CO2 and argon in ITQ-29, which were used as experimental references. The financial support by SOL Group and by Compagnia di San Paolo (CSP2014-HEPYCHEM project) is acknowledged.

References (1) Guerra, G.; Daniel, C.; Rizzo, P.; Tarallo, O. Advanced Materials Based on Polymer Cocrystalline Forms. J. Polym. Sci., Part B: Polym. Phys. 2012, 50, 305–322. (2) Kim, W.-G.; Nair, S. Membranes from Nanoporous 1D and 2D Materials: a Review of Opportunities, Developments, and Challenges. Chem. Eng. Sci. 2013, 104, 908–924. (3) Ouyang, Y.; Shi, H.; Fu, R.; Wu, D. Highly Monodisperse Microporous Polymeric and Carbonaceous Nanospheres with Multifunctional Properties. Sci. Rep. 2013, 3 . (4) Matsuda, R. Design and Synthesis of Porous Coordination Polymers Showing Unique Guest Adsorption Behaviors. Bull. Chem. Soc. Japan 2013, 86, 1117– 1131.

19



(5) Broom, D.; Thomas, K. Gas Adsorption By Nanoporous Materials: Future 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Applications and Experimental Challenges. MRS Bull. 2013, 38, 412–421. (6) Casco, M.; Martínez-Escandell, M.; Gadea-Ramos, E.; Kaneko, K.; SilvestreAlbero, J.; Rodríguez-Reinoso, F. High-pressure Methane Storage in Porous Materials: Are Carbon Materials in the Pole Position? Chem. Mater. 2015, 27, 959–964. (7) Zheng, Y.; Liu, J.; Liang, J.; Jaroniec, M.; Qiao, S. Graphitic Carbon Nitride Materials: Controllable Synthesis and Applications in Fuel Cells and Photocatalysis. Energy Environ. Sci. 2012, 5, 6717–6731. (8) Rose, M. Nanoporous Polymers: Bridging the Gap Between Molecular and Solid Catalysts? ChemCatChem 2014, 6, 1166–1182. (9) Pal, N.; Bhaumik, A. Mesoporous Materials: Versatile Supports in Heterogeneous Catalysis for Liquid Phase Catalytic Transformations. RSC Adv. 2015, 5, 24363–24391. (10) Díaz, U.; Corma, A. Ordered Covalent Organic Frameworks, COFs and PAFs. from Preparation to Application. Coord. Chem. Rev. 2016, 311, 85–124. (11) Jelfs, K.; Cooper, A. Molecular Simulations to Understand and to Design Porous Organic Molecules. Curr. Opin. Solid State Mater. Sci. 2013, 17, 19– 30. (12) Krishna, R. Diffusion in Porous Crystalline Materials. Chem. Soc. Rev. 2012, 41, 3099–3118. (13) Düren, T.; Bae, Y.-S.; Snurr, R. Using Molecular Simulation to Characterise Metal-organic Frameworks for Adsorption Applications. Chem. Soc. Rev. 2009, 38, 1237–1247. (14) Haldoupis, E.; Nair, S.; Sholl, D. Pore Size Analysis of >250 000 Hypothetical Zeolites. Phys. Chem. Chem. Phys. 2011, 13, 5053–5060. (15) Wilmer, C.; Leaf, M.; Lee, C.; Farha, O.; Hauser, B.; Hupp, J.; Snurr, R. Largescale Screening of Hypothetical Metal-organic Frameworks. Nature Chem. 2012, 4, 83–89. 20


Page 20 of 43

Page 21 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(16) Abbott, L. J.; Colina, C. M. Atomistic Structure Generation and Gas Adsorption Simulations of Microporous Polymer Networks. Macromolecules 2011, 44, 4511–4519. (17) Jiang, S.; Jelfs, K.; Holden, D.; Hasell, T.; Chong, S.; Haranczyk, M.; Trewin, A.; Cooper, A. Molecular Dynamics Simulations of Gas Selectivity in Amorphous Porous Molecular Solids. J. Am. Chem. Soc. 2013, 135, 17818– 17830. (18) Mendoza-Cortes, J.; Pascal, T.; Goddard, W. Design of Covalent Organic Frameworks for Methane Storage. J. Phys. Chem. A 2011, 115, 13852–13857. (19) Trewin, A.; Willock, D. J.; Cooper, A. I. Atomistic Simulation of Micropore Structure, Surface Area, and Gas Sorption Properties for Amorphous Microporous Polymer Networks. J. Phys. Chem. C 2008, 112, 20549–20559. (20) Ohba, T.; Yamamoto, S.; Takase, A.; Yudasaka, M.; Iijima, S. Evaluation of Carbon Nanopores Using Large Molecular Probes in Grand Canonical Monte Carlo Simulations and Experiments. Carbon 2015, 88, 133–138. (21) Han, S.; Kim, D.; Jung, D.; Cho, S.; Choi, S.-H.; Jung, Y. Accurate Ab Initiobased Force Field for Predictive CO 2 Uptake Simulations in MOFs and ZIFs: Development and Applications for MTV-MOFs. J. Phys. Chem. C 2012, 116, 20254–20261. (22) Zwijnenburg, M.; Cheng, G.; McDonald, T.; Jelfs, K.; Jiang, J.-X.; Ren, S.; Hasell, T.; Blanc, F.; Cooper, A.; Adams, D. Shedding Light on Structureproperty Relationships for Conjugated Microporous Polymers: the Importance of Rings and Strain. Macromolecules 2013, 46, 7696–7704. (23) Garberoglio, G.; Pugno, N.; Taioli, S. Gas Adsorption and Separation in Realistic and Idealized Frameworks of Organic Pillared Graphene: a Comparative Study. J. Phys. Chem. C 2015, 119, 1980–1987. (24) Collell, J.; Galliero, G.; Gouth, F.; Montel, F.; Pujol, M.; Ungerer, P.; Yiannourakou, M. Molecular Simulation and Modelisation of Methane/ethane Mixtures Adsorption Onto a Microporous Molecular Model of Kerogen Under Typical Reservoir Conditions. Micropor. Mesopor. Mat. 2014, 197, 194–203. 21



(25) Li, K.; Cao, D. The Bridge Between First-principles Calculations and Grand 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Canonical Monte Carlo Simulations: Morse and Lennard-Jones Force Fields. Mol. Simul. 2010, 36, 1157–1163. (26) Lan, J.; Coa, D.; Wang, W. High Uptakes of Methane in Li-doped 3D Covalent Organic Frameworks. Langmuir 2010, 26, 220–226. (27) Cacelli, I.; Cinacchi, G.; Prampolini, G.; Tani, A. Computer Simulation of Solid and Liquid Benzene with An Atomistic Interaction Potential Derived from Ab Initio Calculations. J. Am. Chem. Soc. 2004, 126, 14278–14286. (28) Fu, C.-F.; Tian, S. A Comparative Study for Molecular Dynamics Simulations of Liquid Benzene. J. Chem. Theory Comput. 2011, 7, 2240–2252. (29) Mendoza-Cortés, J.; Han, S.; Furukawa, H.; Yaghi, O.; Goddard III, W. Adsorption Mechanism and Uptake of Methane in Covalent Organic Frameworks: Theory and Experiment. J. Phys. Chem. A 2010, 114, 10824–10833. (30) Ward, B.; Getman, R. Molecular Simulations of Physical and Chemical Adsorption Under Gas and Liquid Environments Using Force Field- and Quantum Mechanics-based Methods. Mol. Simul. 2014, 40, 678–689. (31) Grimme, S. A General Quantum Mechanically Derived Force Field (QMDFF) for Molecules and Condensed Phase Simulations. J. Chem. Theory Comput. 2014, 10, 4497–4514. (32) Cacelli, I.; Cimoli, A.; Livotto, P.; Prampolini, G. An Automated Approach for the Parameterization of Accurate Intermolecular Force-fields: Pyridine as a Case Study. J. Comput. Chem. 2012, 33, 1055–1067. (33) Prampolini, G.; Livotto, P.; Cacelli, I. Accuracy of Quantum Mechanically Derived Force-fields Parameterized from Dispersion-corrected DFT Data: the Benzene Dimer as a Prototype for Aromatic Interactions. J. Chem. Theory Comput. 2015, 11, 5182–5196. (34) Schmidt, J.; Yu, K.; McDaniel, J. Transferable Next-Generation Force Fields from Simple Liquids to Complex Materials. Acc. Chem. Res. 2015, 48, 548– 556.

22


Page 22 of 43

Page 23 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(35) Paesani, F. Getting the Right Answers for the Right Reasons: Toward Predictive Molecular Simulations of Water with Many-Body Potential Energy Functions. Acc. Chem. Res. 2016, 49, 1844–1851. (36) øRezáč, J.; Hobza, P. Describing Noncovalent Interactions Beyond the Common Approximations: How Accurate Is the "gold Standard," CCSD(T) At the Complete Basis Set Limit? J. Chem. Theory Comput. 2013, 9, 2151–2155. (37) Brandenburg, J.; Hochheim, M.; Bredow, T.; Grimme, S. Low-cost Quantum Chemical Methods for Noncovalent Interactions. J. Phys. Chem. Lett. 2014, 5, 4275–4284. (38) Schwabe, T.; Grimme, S. Towards Chemical Accuracy for the Thermodynamics of Large Molecules: New Hybrid Density Functionals Including Non-local Correlation Effects. Phys. Chem. Chem. Phys. 2006, 8, 4398–4401. (39) Grimme, S. Semiempirical GGA-type Density Functional Constructed with a Long-Range Dispersion Correction. J. Comput. Chem. 2006, 27, 1787–1799. (40) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (dft-d) for the 94 Elements H-pu. J. Chem. Phys. 2010, 132, 154104. (41) Grimme, S. Density Functional Theory with London Dispersion Corrections. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2011, 1, 211–228. (42) Reimers, J.; Ford, M.; Goerigk, L. Problems, Successes and Challenges for the Application of Dispersion-corrected Density-functional Theory Combined with Dispersion-based Implicit Solvent Models to Large-scale Hydrophobic Self-assembly and Polymorphism. Mol. Simul. 2016, 42, 494–510. (43) Banerjee, P.; Chandrakumar, K.; Das, G. Exploring Adsorption and Desorption Characteristics of Molecular Hydrogen on Neutral and Charged Mg Nanoclusters: a First Principles Study. Chem. Phys. 2016, 469-470, 123–131. (44) Ricca, C.; Ringuedé, A.; Cassir, M.; Adamo, C.; Labat, F. A Comprehensive DFT Investigation of Bulk and Low-index Surfaces of ZrO2 Polymorphs. J. Comput. Chem. 2015, 36, 9–21.

23



(45) Davis, J.; Baletto, F.; Johnston, R. The Effect of Dispersion Correction on the 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Adsorption of CO on Metallic Nanoparticles. J. Phys. Chem. A 2015, 119, 9703–9709. (46) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, Ö.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09 revision a.1. Gaussian Inc. Wallingford CT, 2009. (47) Lee, S.-J.; Bae, Y.-S. Review of Molecular Simulations of Methane Storage in Metal-organic Frameworks. J. Nanosci. Nanotech. 2016, 16, 4284–4290. (48) Keskin, S.; Sholl, D. Efficient Methods for Screening of Metal Organic Framework Membranes for Gas Separations Using Atomically Detailed Models. Langmuir 2009, 25, 11786–11795. (49) ColÃşn, Y.; Snurr, R. High-throughput Computational Screening of Metalorganic Frameworks. Chem. Soc. Rev. 2014, 43, 5735–5749. (50) Cossi, M.; Gatti, G.; Canti, L.; Tei, L.; Errahali, M.; Marchese, L. Theoretical Prediction of High Pressure Methane Adsorption in Porous Aromatic Frameworks (PAFs). Langmuir 2012, 28, 14405–14414. (51) Fraccarollo, A.; Canti, L.; Marchese, L.; Cossi, M. Monte Carlo Modeling of Carbon Dioxide Adsorption in Porous Aromatic Frameworks. Langmuir 2014, 30, 4147–4156. 24


Page 24 of 43

Page 25 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(52) Babarao, R.; Dai, S.; Jiang, D.-E. Functionalizing Porous Aromatic Frameworks with Polar Organic Groups for High-capacity and Selective CO 2 Separation: a Molecular Simulation Study. Langmuir 2011, 27, 3451–3460. (53) Yang, Z.; Peng, X.; Cao, D. Carbon Dioxide Capture By PAFs and An Efficient Strategy to Fast Screen Porous Materials for Gas Separation. J. Phys. Chem. C 2013, 117, 8353–8364. (54) Cheng, J.; Yuan, X.; Zhao, L.; Huang, D.; Zhao, M.; Dai, L.; Ding, R. GCMC Simulation of Hydrogen Physisorption on Carbon Nanotubes and Nanotube Arrays. Carbon 2004, 42, 2019–2024. (55) Pham, T.; Xiong, R.; Sandler, S.; Lobo, R. Experimental and Computational Studies on the Adsorption of CO2 and N2 on Pure Silica Zeolites. Micropor. Mesopor. Mat. 2014, 185, 157–166. (56) Pham, T.; Forrest, K.; Chen, K.-J.; Kumar, A.; Zaworotko, M.; Space, B. Theoretical Investigations of CO2 and H2 Sorption in Robust Molecular Porous Materials. Langmuir 2016, 32, 11492–11505. (57) Cioce, C.; McLaughlin, K.; Belof, J.; Space, B. A Polarizable and Transferable PHAST N2 Potential for Use in Materials Simulation. J. Chem. Theory Comput. 2013, 9, 5550–5557. (58) Borycz, J.; Tiana, D.; Haldoupis, E.; Sung, J. C.; Farha, O. K.; Siepmann, J. I.; Gagliardi, L. CO2 Adsorption in M-IRMOF-10 (M = Mg, Ca, Fe, Cu, Zn, Ge, Sr, Cd, Sn, Ba). J. Phys. Chem. C 2016, 120, 12819–12830. (59) Shah, M. S.; Tsapatsis, M.; Siepmann, J. I. Monte Carlo Simulations Probing the Adsorptive Separation of Hydrogen Sulfide/Methane Mixtures Using AllSilica Zeolites. Langmuir 2015, 31, 12268–12278. (60) Corma, A.; Rey, F.; Rius, J.; Sabater, M.; Valencia, S. Supramolecular Selfassembled Molecules as Organic Directing Agent for Synthesis of Zeolites. Nature 2004, 431, 287–290. (61) Palomino, M.; Corma, A.; Rey, F.; Valencia, S. New Insights on C02-methane Separation Using LTA Zeolites with Different Si/Al Ratios and a First Comparison with MOFs. Langmuir 2010, 26, 1910–1917. 25



(62) Scuseria, G.; Janssen, C.; Schaefer III, H. An Efficient Reformulation of the 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Closed-shell Coupled Cluster Single and Double Excitation (CCSD) Equations. J. Chem. Phys. 1988, 89, 7382–7387. (63) Pople, J.; Head-Gordon, M.; Raghavachari, K. Quadratic Configuration Interaction. a General Technique for Determining Electron Correlation Energies. J. Chem. Phys. 1987, 87, 5968–5975. (64) Boys, S.; Bernardi, F. Calculation of Small Molecular Interactions by Differences of Separate Total Energies - Some Procedures with Reduced Errors. Mol. Phys. 1970, 19, 553–560. (65) Mayo, S.; Olafson, B.; Goddard III, W. DREIDING: a Generic Force Field for Molecular Simulations. J. Phys. Chem. 1990, 94, 8897–8909. (66) Allen, M.; Tildesley, T. Computer Simulation of Liquids; Clarendon Press, Oxford, 1987. (67) Sun, M.; Shah, D.; Xu, H.; Talu, O. Adsorption Equilibria of C1 to C4 Alkanes, CO2, and SF6 on Silicalite. J. Phys. Chem. B 1998, 102, 1466–1473. (68) Zhu, W.; Hrabanek, P.; Gora, L.; Kapteijn, F.; Moulijn, J. Role of Adsorption in the Permeation of CH4 and CO2 Through a Silicalite-1 Membrane. Ind. Eng. Chem. Res. 2006, 45, 767–776. (69) Lemmon, E. W.; McLinden, M. O.; Friend, D. G. In NIST Chemistry WebBook, NIST Standard Reference Database Number 69 ; J., L. P., G., M. W., Eds.; National Institute of Standards and Technology, 2011. (70) Lim, T.-C. The Relationship between Lennard-Jones (12-6) and Morse Potential Functions. Z. Naturforsch. 2003, 58, 615–617. (71) Lim, T.-C. Alignment of Buckingham Parameters to Generalized LennardJones Potential Functions. Z. Naturforsch. 2009, 64, 200–204. (72) Pastore, H.; De Oliveira, E.; Superti, G.; Gatti, G.; Marchese, L. Reaction At Interfaces: the Silicoaluminophosphate Molecular Sieve CAL-1. J. Phys. Chem. C 2007, 111, 3116–3129.

26


Page 26 of 43

Page 27 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(73) Thommes, M. Chapter 15. Textural Characterization of Zeolites and Ordered Mesoporous Materials By Physical Adsorption. Stud. Surf. Sci. Catal. 2007, 168 . (74) Thommes, M.; Cychosz, K. A. Physical Adsorption Characterization of Nanoporous Materials: Progress and Challenges. Adsorption 2014, 20, 233– 250. (75) McCusker, L.; Liebau, F.; Engelhardt, G. Nomenclature of Structural and Compositional Characteristics of Ordered Microporous and Mesoporous Materials with Inorganic Hosts: (IUPAC Recommendations 2001). Pure Appl. Chem. 2001, 73, 381–394. (76) Landers, J.; Gor, G.; Neimark, A. Density Functional Theory Methods for Characterization of Porous Materials. Colloids Surf. A 2013, 437, 3–32. (77) Liu, B.; Smit, B.; Calero, S. Evaluation of a New Force Field for Describing the Adsorption Behavior of Alkanes in Various Pure Silica Zeolites. J. Phys. Chem. B 2006, 110, 20166–20171. (78) Guliants, V.; Huth, A. Force Fields for Classical Atomistic Simulations of Small Gas Molecules in Silicalite-1 for Energy-related Gas Separations At High Temperatures. J. Porous Mater. 2013, 20, 741–751. (79) Fischer, M.; Bell, R. Influence of Zeolite Topology on CO2/N2 Separation Behavior: Force-field Simulations Using a DFT-derived Charge Model. J. Phys. Chem. C 2012, 116, 26449–26463.

27



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

I

II Figure 1: ITQ-29 model clusters used for FF parameterization: I, Si8 O8 H16 ; II, Si12 O15 H20 .

28


Page 28 of 43

Page 29 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Ar, 87 K

CH4 , 283 K

CO2 , 298 K

Figure 2: Measured and simulated gas densities: black squares, experimental; red circles, Lennard-Jones; blue triangles, Morse.

29



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 3: Interaction energies (CP corrected for BSSE) of CH4 and CO2 with SiH3 − O − SiH3 model at different QM levels.

30


Page 30 of 43

Page 31 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 4: Interaction energies (CP corrected for BSSE) of CH4 and CO2 with SiH3 − O − SiH3 model at different QM levels: effect of the basis set on MP2 and DFT-D3.

31



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 5: Interaction energy scan of CH4 with the ITQ-29 fragment: black squares, QM calculations; blue triangles, fitted FF with Morse function; red circles, fitted FF with LJ function; open symbols, fit to DFT-D3; filled symbols, fit to MP2.

32


Page 32 of 43

Page 33 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 6: Simulated and experimental adsorption isotherms of methane in ITQ-29 at 283 K; black squares, experimental; open red circles, LJ-D3; filled red circles, LJMP2; open blue triangles, Mo-D3; filled blue triangles, Mo-MP2. Note that LJ-MP2 and Mo-D3 curves overlap almost completely.

33



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 7: Formation energies of clusters (Model 1 to 9, see Figure S3 in SI) comprising one CH4 molecule in zeolite models of increasing size; open blue triangles, Mo-D3; filled red circles, LJ-MP2; dotted blue line, Mo-D3 periodic formation energy; solid red line, LJ-MP2 periodic energy.

34


Page 34 of 43

Page 35 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 8: Formation energies of clusters formed by one CH4 molecule and its images in n × n × n supercells (see Figure S4 in SI); open blue triangles, Mo-D3; filled red circles, LJ-MP2.

35



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 9: Interaction energy scan of CO2 with the ITQ-29 fragment I: black squares, QM calculations; blue triangles, fitted FF with Morse function; red circles, fitted FF with LJ function; open symbols, fit to DFT-D3; filled symbols, fit to MP2.

36


Page 36 of 43

Page 37 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 10: Simulated and experimental adsorption isotherms of carbon dioxide in ITQ-29 at 298 K; black squares, experimental; open red circles, LJ-D3; filled red circles, LJ-MP2; open blue triangles, Mo-D3; filled blue triangles, Mo-MP2.

37



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 11: Simulated and experimental adsorption isotherms of argon in ITQ-29 at 87 K; black squares, experimental; filled red circles, LJ-MP2; open blue triangles, Mo-D3. Inset: logarithmic plot of the low pressure region.

38


Page 38 of 43

Page 39 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 12: Pore size distribution (PSD) and cumulative pore volume (CPV, inset) of ITQ-29, obtained from the Ar adsorption isotherms at 87 K by NLDFT analysis with parameters for Ar and silica, and spherical pore shape (ref. 74). Black squares, data from the experimental isotherm; filled red circles, LJ-MP2; open blue triangles, Mo-D3.

39



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 13: Interaction energy scan of CH4 with the ITQ-29 fragment II: black squares, QM calculations; filled red circles, LJ-MP2 (previously fitted on fragment I); purple stars, LJ-MP2-II (fitted on this fragment); open blue triangles, Mo-D3 (previously fitted on fragment I); open green diamonds, Mo-D3-II (fitted on this fragment).

40


Page 40 of 43

Page 41 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 14: Simulated and experimental adsorption isotherms of methane in ITQ-29 at 283 K; black squares, experimental; purple stars, LJ-MP2-II; open green diamonds, Mo-D3-II.

41



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

CH4 in silicalite-1, 303 K

CO2 in silicalite-1, 308 K

CO2 in SSZ-13, 303 K

Figure 15: Simulated and experimental adsorption isotherms of methane and carbon dioxide in various zeolites: black squares, experimental; open red circles, LJ-D3; filled red circles, LJ-MP2; open blue triangles, Mo-D3; filled blue triangles, MoMP2.

42


Page 42 of 43

Page 43 of 43

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


TOC graphic.

43


Accurate Evaluation of the Dispersion Energy in the Simulation of Gas

Recommend Documents