Large-Scale Refinement of Metal−Organic Framework Structures

Nov 15, 2016 - Atomistic simulations of metal–organic frameworks (MOFs) depend on ...... Smit , B. The materials genome in action: identifying the p...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/cm

Large-Scale Refinement of Metal−Organic Framework Structures Using Density Functional Theory Dalar Nazarian,† Jeffrey S. Camp,† Yongchul G. Chung,‡,§ Randall Q. Snurr,§ and David S. Sholl*,† †

School of Chemical & Biomolecular Engineering, Georgia Institute of Technology, Atlanta, Georgia 30332, United States School of Chemical and Biomolecular Engineering, Pusan National University, Busan 609-735, South Korea § Department of Chemical and Biological Engineering, Northwestern University, Evanston, Illinois 60208, United States ‡

S Supporting Information *

ABSTRACT: Efforts to computationally characterize large numbers of nanoporous materials often rely on databases of experimentally resolved crystal structures. The accuracy of experimental crystal structures used in such calculations has a significant impact on the reliability of the results. In this work, we report structures optimized using periodic density functional theory (DFT) for more than 800 experimentally synthesized metal−organic frameworks (MOFs). Many MOFs changed significantly upon structural optimization, particularly materials that were crystallographically resolved in their solvated form. For each MOF, we simulated the adsorption of CH4 and CO2 using grand canonical Monte Carlo both before and after DFT optimization. The DFT optimization has a large impact on simulated gas adsorption in some cases. For example, CO2 loading at 1 bar changed by more than 25% in over 25% of the MOFs we considered.

A

upon desolvation;10 second, there may be structural defects in experimental samples or trapped solvent molecules from incomplete activation11 and third, there can be variations between even high-quality MOF samples produced using various synthesis conditions. The variance between experimentally reported structures may be observed between different representations of the same MOF deposited in the CSD. There are at least 50 entries for CuBTC (HKUST-1) in the CSD reported by research groups using different synthesis conditions, activation conditions, and XRD refinement procedures.3 The a direction lattice constants of these structures vary between 18.34 and 18.66 Å (see Figure S1), and there are small differences in textural descriptors such as the largest cavity diameters. Because different structures may be of interest to different researchers, these “duplicate” structures were retained in the original CoRE MOF database. This adds the complication of having to select the highest quality or most applicable representation of duplicated structures when using the CoRE MOF database. Small differences in the crystal structures can have important consequences in the modeling of adsorption properties of MOFs. For example, Lawler and co-workers found that GCMC-simulated isotherms and heats of adsorption for Kr and Xe varied significantly between 12 different experimentally resolved HKUST-1 structures.12 Considerable variance has also

tomistic simulations of metal−organic frameworks (MOFs) depend on accurate crystallographic structural data obtained from diffraction using X-rays or neutrons.1 More than 90% of MOF structures reported in the Cambridge Structural Database (CSD)2 were resolved using single-crystal X-ray diffraction (XRD). In most cases, MOF structures are resolved in their as-synthesized form with solvent molecules in their nanopores even though most applications of MOFs envision the use of activated materials for which these solvent species have been removed. Other MOF structures deposited in the CSD include molecular guests from in situ crystallographic studies of adsorption.1 These guest molecules must be removed prior to computer simulations of fully activated structures. In addition, researchers must select a single representation of any partially occupied or disordered atoms in experimentally refined MOF structures. In previous work, we performed these preparation steps on over 5,000 experimental MOF structures from the CSD.3 The result of that work, the ComputationReady Experimental MOF (CoRE MOF) database, eliminates a significant initial hurdle to high-throughput atomistic simulations of MOFs.3,4 The CoRE MOF database has been used to screen MOF structures for utility in methane storage and purification,3,5 CO2 separations,6,7 CH4/H2 separations,8 and Xe/Kr separations9 using grand-canonical Monte Carlo (GCMC) simulations of adsorption with classical potentials. A key assumption made in each of these GCMC studies is that experimental MOF samples are well represented by the structures that have been completely desolvated in silico. However, this assumption may be problematic for three reasons: first, some MOFs undergo large phase transitions or even completely collapse © 2016 American Chemical Society

Special Issue: Computational Design of Functional Materials Received: October 4, 2016 Revised: November 15, 2016 Published: November 15, 2016 2521

DOI: 10.1021/acs.chemmater.6b04226 Chem. Mater. 2017, 29, 2521−2528

Article

Chemistry of Materials

sets for all elements except La. The basis set used for La is provided in the Supporting Information. Spin polarization was considered for all calculations. Geometry optimization of the system was performed using a BFGS optimizer, allowing for full atomic and cell relaxation without symmetry restrictions until the largest force on atoms reached less than 0.0003 hartree/Bohr. Optimization was attempted for MOFs in the CoRE MOF database without lanthanide or actinide metals because basis sets for these atom types are missing from CP2K. As an exception, MOFs with lanthanum were included in the attempted optimization (see Supporting Information). Of the 2,612 structures, calculations for 879 structures converged successfully within a reasonable computation time (see Supporting Information for details.) Structure Analysis. The changes in geometry associated with optimization were analyzed for each successfully optimized structure. Structural parameters considered include unit cell parameters (a, b, c, α, β, and γ), unit cell volume, helium void fraction, largest cavity diameter (LCD), and pore-limiting diameter (PLD).28 The void fraction was computed with the Widom particle insertion method using a helium atom probe.29 The LCD, PLD, and pore size distributions were calculated using the Zeo++ geometry analysis package.30 Methane and CO2 Adsorption. Classical grand canonical Monte Carlo (GCMC) simulations of methane adsorption were conducted using the RASPA 2.0 code31 for 465 structures. These 465 structures have readily available DDEC point charges to use for GCMC calculations. Methane adsorption was simulated at 65 bar and 298 K. The Peng−Robinson equation of state32 was used to calculate the fugacity values necessary to impose equilibrium between the simulation cell and the external gas reservoir. Methane−methane and methane−framework interactions were modeled with a LennardJones (LJ) 12-6 potential using the Lorentz−Berthelot mixing rules. The LJ parameters for all framework atoms were obtained from the Universal Force Field (UFF).33 The LJ parameters for methane (ε/kB = 148.0 K; σ = 3.73 Å) were obtained from the TraPPE force-field, modeled as a single sphere with one LJ interaction site.34 All LJ interaction potentials were truncated at 12.8 Å. For the minimum image convention to be satisfied, all simulation cells were replicated to at least 25.6 Å along each axis. All GCMC simulations included a 2,500 cycle equilibration period followed by a 2,500 cycle production period.35 GCMC simulations included random insertion, deletion, translation, and reinsertion moves with equal probabilities. GCMC calculations of CO2 adsorption were conducted for structures with available charges in the MOF DDEC charge database.4 Simulations were performed for 465 structures at 1 bar and 298 K. DDEC point charges from the experimentally observed structures were directly mapped onto the optimized structures. Tests for multiple structures showed that the DFT optimization did not lead to changes in atomic coordination, and the changes in DDEC point charges associated with relaxation were negligible (see Figure S2). The LJ parameters for CO2 were obtained from the TraPPE force field modeled as a 3-site molecule as described in the Supporting Information.36 CO2 framework LJ parameters were defined using Lorentz−Berthelot combining rules in the same way as our calculations for methane adsorption.

been reported in simulated adsorption capacities for CH4 at 298 K and 65 bar in 11 different MIL-53(Al) structures from the CoRE MOF database.3 Structural uncertainties in crystals can potentially be resolved through computational energy minimization to obtain optimized structures. These efforts typically result in small changes relative to the initial experimental structure and are referred to by Catlow and co-workers as “structural refinements”.13 Computational structural refinements are distinct from crystal structure prediction in that the experimental atomic positions are used as a starting point for optimization. Density functional theory (DFT) is the most accurate computational approach that can be feasibly applied to periodic crystal structural refinement.14 A range of studies have used DFT to refine and analyze large sets of crystal structures. For example, Sharma et al. conducted a search for novel dielectric polymeric materials by screening 1D repeat units for high dielectric constants using DFT.15 Chandrasekhar et al. used DFT methods to predict solubility, diffusivity, and permeability of hydrogen in intermetallic membrane materials.16 Armiento et al. performed a large-scale DFT study of the ABO3 chemical space in the perovskite crystal structure to identify promising piezoelectric materials.17 Yan et al. screened hundreds of transition metal oxides for photocatalytic materials for water splitting.18 Nicholson et al. used DFT-based methods to assess this stability of metal hydride systems by predicting thermodynamic properties.19 The Materials Project has applied DFT+U methods using a combination of sophisticated high-throughput infrastructure and crowd-sourcing to optimize thousands of structures obtained from the Inorganic Crystal Structure Database and predict a range of material properties.20,21 Because of the relatively high computational cost of such calculations, DFT-based structure refinements of MOFs have previously only been applied to a small number of structures.14 To date, all high-throughput screening studies of MOFs have used experimentally observed crystal structures or collections of hypothetical structures generated in silico. In this work, we systematically refine a diverse set of over 800 MOF structures using periodic DFT. For each starting and geometry-optimized structure, we calculate adsorption properties of CH4 and CO2 at conditions relevant to adsorbed natural gas storage and CO2 capture from flue gas, respectively. The goals of this work are three-fold. First, we evaluate the extent to which MOFs change structurally after energy minimization. Second, we evaluate how DFT optimization changes GCMC-predicted gas adsorption properties relative to the starting CoRE MOF structure. Finally, for structures with multiple representations in the CoRE MOF database, such as HKUST-1, we evaluate whether these structures converge to a self-consistent representation after DFT optimization. The DFT-optimized structures in this work are publicly available at http://dx.doi.org/10.11578/1118280 as a supplement to the CoRE MOF database.





METHODS

RESULTS AND DISCUSSION

Structure Refinement. To compare the extent and type of geometry changes associated with each MOF structure after energy minimization, we considered changes in four structural parameters: the unit cell volume, the helium void fraction, the unit cell angles, and largest cavity diameter (LCD). A majority of structures showed less than 10% change in the unit cell volume and the void fraction and 2° change in unit cell angles after optimization (see Figure S3). We also find that most structures experience small changes in the pore diameter where more than 90% of structures experience less than 1 Å change in

Structure Refinement. All DFT calculations were carried out using the Gaussian plane-wave (GPW) computational package CP2K 2.622 on the Argonne National Laboratories supercomputer MIRA. On the basis of the results of our previous benchmarking study,23 we chose to use the Gordecker, Teter, and Hutter dual-space pseudopotentials (GTH)24 with the PBE-D325,26 functional. Dispersion coefficients used in PBE-D3 are geometry dependent and are adjusted on the basis of the local coordination number around the atoms of interest. After a series of convergence tests, a plane-wave energy cutoff of 800 Ry was chosen. We have used the double-ζ valence polarized (DZVP)27 basis 2522

DOI: 10.1021/acs.chemmater.6b04226 Chem. Mater. 2017, 29, 2521−2528

Article

Chemistry of Materials

Figure 1. Box and whisker diagrams for change in volume, change in LCD, change in cell angles, and change in void fraction for structures after DFT optimization. Each diagram distinguishes results between structures with and without solvent in the CSD version of the structure. The box represents the first quartile (Q1) and third quartile (Q3) with the center line representing the data median. Going outward, the set of horizontal lines represent the Q1 − 1.5 × IQR and Q3 + 1.5 × IQR, where IQR is the interquartile range Q3−Q1. The points outside these bounds represent the outliers in each data set.

Figure 2. (left) Representation of JUC-63 as found in the CSD containing both bound (purple) and free (green) DMF molecules. (middle) Representation of JUC-63 discussed in synthesis literature and found in the CoRE MOF database. (right) Representation of JUC-63 after DFT minimization, showing that the MOF changes drastically after activation.

bound solvent and those with only free solvent. We find only small differences in the impact of solvent type on changes in structural parameters (see Figure S4). We find that some MOFs show substantial changes in their structure following DFT minimization. For example, JUC-63 (CSD REFCODE: OFODAP)37 shows a substantial change in structure after geometry optimization. Figure 2(a) shows a representation of the MOF as found in the CSD. This structure includes free DMF molecules inside the MOF pores as well as DMF molecules bound to the Cd metal center. Figure 2(b) shows the CoRE MOF database version of OFODAP in which all DMF molecules present during experimental crystal structure refinement have been removed in silico. When the desolvated structure is energy minimized, the structure deforms substantially and decreases in volume (see Figure 2(c)). This DFT-optimized structure is 44.4% smaller in volume than the original structure with an LCD 2 Å smaller than the experimentally reported structure. The observation that a MOF structure can change significantly upon solvent removal is not surprising. Indeed, one of the major advances in early work on MOFs was the discovery of materials that did not undergo pore collapse after solvent removal. The results above, however, are a reminder that not all materials are stable upon solvent removal. Our results have important implications for efforts to use highthroughput computational methods to assess properties of

the LCD. We find little correlation between changes in volume and void fraction. Impact of Residual Solvent. Although most MOFs in the CoRE MOF database show only small structural changes between their experimental and DFT-optimized representations, a subset of MOFs exist with substantial changes in their structures after energy minimization. Examination of these materials indicates that this subset is comprised mainly of structures that contained residual solvent in the structure during crystal structure refinement. In the process of constructing the CoRE MOF database, Chung et al. identified MOFs with and without solvent in the crystallographic structure information provided in the CSD.3 We used this information to divide the 879 structures from our DFT calculations into these two categories. When the structural parameters discussed above (overall change in unit cell volume, helium void fraction, unit cell angle, and LCD) are analyzed for each group, we find a stark difference between the groups on the impact of optimization on geometry. Figure 1 shows a box and whisker diagram for each structural parameter with and without solvent in the CSD version of the structure. Experimental MOF structures that included solvent during structural analysis display a larger number and magnitude of outliers. MOF structures with solvent in the CSD representation can further be divided into two categories: those with at least one 2523

DOI: 10.1021/acs.chemmater.6b04226 Chem. Mater. 2017, 29, 2521−2528

Article

Chemistry of Materials

3% increase in the unit cell volume and a 0.047 change in void fraction after energy minimization. Even with this apparently small change in geometric properties, the computed methane uptake increases from 8.48 to 25.97 cm3 of adsorbate (STP)/ cm3 of MOF. In another example, NASCIV,39 a Cd-based MOF, shows an 8% change in the unit cell volume and 0.04 change in the helium void fraction after DFT optimization but an increase in CO2 uptake from 9.28 to 18.26 cm3 of adsorbate (STP)/cm3 of MOF. These examples highlight the observation that adsorption properties can depend in a complex manner on the free energy landscape of the adsorbate molecule within the structure. Figure 4 shows the relationship between change in adsorbate

MOFs because they represent an important refinement to the original CoRE MOF database, which did not consider the potential structural effects due to solvent removal. Adsorption Properties. GCMC simulations of methane uptake at 65 bar and CO2 uptake at 1 bar were performed for the experimentally reported and DFT-optimized representations of the 465 structures for which DDEC point charges are available.23 For consistency, all charges, including charges for MOFs with significant structural changes, were assigned based on the CoRE MOF DDEC Point Charge database.4 The framework atomic positions were held fixed during all GCMC calculations. The simulations show that methane adsorption at 65 bar changed by more than 25% for 22% of structures after structural optimization, whereas CO2 uptake at 1 bar changed by more than 25% for 28% of structures (Figure S6). Among the 465 MOFs, the percent difference in adsorbate loading after geometry optimization is uncorrelated with the magnitude of uptake in the original CoRE structure. Qualitatively, a 5% difference in simulated uptake between the initial and refined structures might be considered small. For CO2, more than 82% of the structures we examined showed changes larger than a 5% threshold. This is a significant observation for efforts to examine the properties of MOFs using computational highthroughput methods and indicates that calculated results that do not use structurally refined structures may be subject to considerable uncertainty. Figure 3 compares the changes associated with DFT optimization of MOF structures for CO2 and CH4 uptake.

Figure 4. Percent difference in CO2 (1 bar and 298 K) and methane (65 bar and 298 K) uptake between the structure found in the CoRE MOF database and DFT-optimized structure as a function of the original structure’s LCD.

uptake and the pore diameter of the original structure. The results show that almost all structures with a significant difference in uptake have pore sizes less than 5 Å, indicating a strong adsorbate interaction with the pore walls. It is important to note, however, that not all MOFs with an LCD of less than 5 Å behave this way. Of the 154 structures with an LCD of less than 5 Å, 60 and 56 structures showed a difference in computed CO2 and methane uptake of less than 25% after DFT optimization of the structure, respectively. Figure 4 shows that, for structures with 6 Å or larger LCDs, methane uptake is less sensitive to structural change. This is not true for CO2 uptake, where structures with LCD of 12 Å can still undergo up to 100% change in uptake upon structural refinement. Among the 465 structures were the six different representations of MIL-53(Al)-lp found in the CoRE MOF database. The standard deviation of the unit cell volume of these structures is reduced from 609 to 162 Å3 following DFT-based geometry optimization. The standard deviation in LCD is 0.18 Å after minimization, compared to 0.49 Å in the original structures. The existence of the range of structures even after DFT optimization highlights the need to investigate multiple crystallographic structures of the MOF of interest when available. To investigate the effect of these structural changes on adsorption properties, we calculated methane adsorption isotherms for the six structures using both experimentally observed and DFT-optimized structures (see Figure 5). The results show that using DFT-optimized structures in GCMC

Figure 3. Percent difference in CO2 (1 bar and 298 K) and methane (65 bar and 298 K) uptake and change in void fraction between structures found in the CoRE MOF database and DFT-optimized structures.

These results show that MOFs with significant changes in the helium void fraction after DFT optimization tend to exhibit larger differences in GCMC-predicted adsorption properties between experimental and optimized structures. For structures with minor changes in geometry, we find only a weak correlation between change in void fraction and uptake. As shown in Figure 3, some structures experience a change in the helium void fraction of less than 0.05 but experience a more than 150% increase in the calculated CO2 or methane uptake. For example, HUHJAW,38 a Cd and Cl based MOF, shows a 2524

DOI: 10.1021/acs.chemmater.6b04226 Chem. Mater. 2017, 29, 2521−2528

Article

Chemistry of Materials

in some cases, lead to changes in predicted adsorption properties relative to those in a rigid structure.42 For the sensitivity of CO2 uptake in CICYIX to be examined further, a 400 fs ab initio molecular dynamics (AIMD) calculation was carried out in the empty MOF to obtain different representations of the MOF under thermal motion. This approach, a “snapshot method”, does not require specification of force field parameters for the MOF degrees of freedom.43,44 Further details on AIMD calculations are given in the Supporting Information. GCMC calculations for CO2 adsorption were performed for four snapshots taken from the AIMD simulations at 2, 260, 310, and 360 fs. The snapshot at 2 fs represents a minor change in atomic positions of the structure relative to the DFT-minimized structure. As expected, we find that these minor changes do not impact the CO2 uptake (see Figure 7(b)). For more significant thermal displacements in the

Figure 5. Methane adsorption in experimentally resolved MIL-53(Al) structures (orange) and their DFT-minimized counterparts (blue).

calculations leads to more consistent methane adsorption isotherms, although there is still some variation in predicted uptake among DFT-optimized structures. This variation may be more prevalent in flexible MOFs that exhibit “breathing” modes because these modes correspond to deformations in the crystal structure that are associated with small changes in total energy.40 For example, we find that, although there are many representations of HKUST-1 in the CoRE MOF database, the differences among the HKUST-1 structures found in the database are much smaller than for the MIL-53 structures found in the database. Although DFT optimization of HKUST1 reduces the differences among structures, the impact of DFT optimization on the geometric properties is small and has negligible impact on the methane adsorption properties of methane (see Figure S7). For some MOFs, we find that the geometry optimization leads to an increase in methane uptake but a decrease in CO2 uptake. One such example is CICYIX,41 a Cd-based MOF shown in Figure 6. During geometry optimization, the unit cell

Figure 7. (a) Heat of adsorption (b) and uptake of CO2 at 1 bar and 298 K among AIMD snapshots of CICYIX. (c) Framework of CICYIX (shown in gray) and adsorbed CO2 molecules (shown in red) during GCMC calculations at 1 bar and 298 K. (d) Pore size distribution of each framework (calculated using Zeo++).30 Results show that choice of framework can drastically affect calculated results.

Figure 6. Framework of CICYIX (shown in gray) and adsorbed CO2 molecules (shown in red) during GCMC calculations at 1 bar and 298 K for (a) the structure of CICYIX found in the CoRE MOF database and (b) the structure after energy minimization.

atomic positions and changes in the unit cell to be evaluated, three snapshots were obtained at 50 fs intervals starting at the arbitrarily chosen 260 fs AIMD run. The structures chosen in this way show less than 1° change in the unit cell angles and 0.5 Å change in the unit cell lengths. All reported CO2 uptakes and the heats of adsorption for CICYIX were calculated by averaging and calculating the standard error among 10 GCMC calculations at 1 bar and 298 K. Adsorption studies in these snapshots show CO2 uptakes ranging from 20 to 45 cm3 (STP) /cm3, i.e., results that are very different from the negligible adsorption predicted in the DFT-minimized structure. As shown in Figure 7(d), the DFT optimization of the structure results in a shift to larger pores. After 2 fs of AIMD in the DFT-optimized structure, we find an additional increase in pore size. After 260 fs of AIMD, we find a wider range of pore sizes. This structure includes pore sizes observed in both the

of CICYIX transitions from an orthorhombic structure to a triclinic structure. The LCD of the structure increases by 5%, and framework density increases by 3%. This results in an increase from 0.22 to 15.87 cm3 (STP)/cm3 in methane uptake. CO2 uptake drops from a relatively high uptake of 154.40 cm3 (STP)/cm3 to almost no uptake of 3.29 cm3 (STP)/cm3. As shown in Figure 6, in the experimental structure, CO2 adsorbs in a highly specific location in every pore during a GCMC simulation. That is, there are specific pockets within the MOF pore that are ideal for CO2 adsorption. Internal Degrees of Freedom and Their Impact on Adsorption Properties. All of the GCMC simulations described above treat the MOF framework as rigid, so they neglect possible impacts from internal motions of the MOF or coupling between the MOF and adsorbates. Recent calculations have indicated that thermal motions of MOF frameworks may, 2525

DOI: 10.1021/acs.chemmater.6b04226 Chem. Mater. 2017, 29, 2521−2528

Article

Chemistry of Materials

results indicate that adsorption properties can be dependent in a complex manner on the free energy landscape of the adsorbate molecule within the structure. We have also shown that almost all of the structures with a significant difference in uptake have small pore sizes because of the importance of adsorbate interaction with the pore walls. We further examined one MOF, CICYIX, where the adsorbate−wall interaction makes the computed adsorption properties of CO2 highly sensitive to the changes in the framework structure. These results indicate that, for MOFs with pore sizes close to those of the adsorbate in question, framework dynamics due to thermal motions can significantly impact computationally calculated adsorption properties. In structures of this type, performing GCMC in a rigid DFTrelaxed structure may not be sufficient to accurately predict the adsorption properties of the real material. Our results show a striking example of this situation in a specific MOF, but it is not possible from our results to date to predict how common this situation is. Although such structures require more computationally demanding simulations, incorporation of framework flexibility in adsorption calculations could enable more accurate prediction of adsorption properties, which in turn can enable accelerated design of the MOFs that are highly CO 2 selective.45−47 Our results have significant implications for the MOF simulation community. So far, framework flexibility is typically ignored in the simulation of adsorption for MOFs. We have considered two phenomena that may be relevant when performing adsorption calculations in MOFs. First, because adsorption calculations typically treat the adsorbing framework as rigid, care should be taken to use the most physically accurate structure possible. Our results show that if the MOF structure was obtained from experiments with residual solvent present in the pores, the structure should ideally be DFT optimized. Even if the structure does not change significantly after energy minimization, GCMC calculations in the new DFT-optimized structure may result in significantly different uptake, especially for polar molecules where the free energy of adsorption inside the MOF changes in a complex manner. We feel that whenever a DFT-optimized structure of the kind we have reported here is available for a particular MOF, that structure should be used in adsorption calculations preferentially over a purely experimental structure. Even if this statement is too strong, performing calculations that examine the change in properties between the experimental and DFToptimized structure should become a routine approach. This will provide more consistency to the adsorption simulations performed in the MOF community. For these reasons, we anticipate that our database of refined structures will have important implications for such high-throughput MOF screening efforts and recommend that the use of structurally refined materials be preferred in studies of this kind whenever possible. We also considered whether fluctuations in MOF crystal structures due to thermal vibrations can have a significant impact on simulated adsorption uptake. For MOFs with pore sizes comparable to the size of the adsorbate molecule, small fluctuations in the pore size can drastically change the adsorption properties in at least some cases. In this regard, the consideration of framework vibrations and related effects appears to be important for reliable prediction of adsorption. Further development of computationally efficient approaches that incorporate these phenomena will significantly advance the

original structure and the DFT-optimized structure. Upon visualization of adsorption (see Figure 7(c)), we observe a diverse range of CO2 adsorption sites leading to a different amount of CO2 uptake and heat of adsorption. The significant variation in adsorption properties observed among these CICYIX structures highlights the sensitivity of the MOF’s CO2 adsorption properties to structure. Our results for adsorption of CO2 in CICYIX indicate that, for MOFs with pore sizes close to those of the adsorbate in question, framework dynamics associated with thermal motions of framework atoms can significantly impact the calculated adsorption properties. Small fluctuations in the pore size can drastically change uptake properties. For such situations, consideration of framework flexibility may be essential for a reliable predication of adsorption. It is currently difficult to speculate on how prevalent this kind of sensitivity is because we have only examined this internal flexibility effect in one material. As mentioned above, Gee et al. have demonstrated that framework flexibility is crucial for accurately predicting selectivity for o-xylene/ethylbenzene (oX/eb) in MIL-47.42 In their study, Gee et al. showed that GCMC with a rigid DFToptimized structure for MIL-47 predicts a selectivity that is an order of magnitude larger than the experimentally predicted value. For this discrepancy to be resolved, the authors captured flexibility in MIL-47 using a “snapshot method”, which led to a more accurate prediction of selectivity closer in magnitude to the experimental result. Similar to the methods used to study CICYIX above, Gee et al. obtained structural snapshots from a fully flexible molecular dynamics simulation of an empty framework of MIL-47. They used the snapshots as input to a GCMC simulation of adsorption where the framework atoms are held fixed during the simulation. In both cases, the calculations do not include possible effects associated with coupling between the adsorbed molecules and flexible degrees of freedom in the porous material. Taken together, these results indicate that further development of computationally inexpensive methods to assess the impact of internal flexibility on adsorption in MOFs may be important for future improvements in the accuracy of computational predictions of these properties.



SUMMARY We have prepared a diverse set of 879 DFT-optimized MOF structures that is publicly available at http://dx.doi.org/10. 11578/1118280 as a supplement to the CoRE MOF database. We assessed the change to the experimentally refined structures upon DFT optimization and found that a majority of the structures undergo less than 10% change in structural parameters such as the pore size, the unit cell length and angles, the unit cell volume, and the helium void fraction. Moreover, we showed that the MOFs whose crystal structures were solved with solvents in the pores undergo a larger structural change during the geometry optimization. We used the large set of DFT-optimized structures to assess the effect of structural changes on adsorption properties of the MOFs. By studying uptake of CH4 and CO2 before and after structural relaxation, we have shown that, for the majority of MOFs, methane uptake is not sensitive to small structural changes. However, we find that, for 82% of MOFs (382 structures), CO2 uptake changed by more than 5% upon minimization of the structure. We found a weak correlation between change in structure and change in uptake. These 2526

DOI: 10.1021/acs.chemmater.6b04226 Chem. Mater. 2017, 29, 2521−2528

Article

Chemistry of Materials

(7) Li, S.; Chung, Y. G.; Snurr, R. Q. High-Throughput Screening of Metal−Organic Frameworks for CO2 Capture in the Presence of Water. Langmuir 2016, 32, 10368−10376. (8) Erucar, I.; Keskin, S. Computational assessment of MOF membranes for CH4/H2 separations. J. Membr. Sci. 2016, 514, 313− 321. (9) Simon, C. M.; Mercado, R.; Schnell, S. K.; Smit, B.; Haranczyk, M. What are the best materials to separate a xenon/krypton mixture? Chem. Mater. 2015, 27, 4459−4475. (10) Coudert, F. X.; Boutin, A.; Jeffroy, M.; Mellot-Draznieks, C.; Fuchs, A. H. Thermodynamic Methods and Models to Study Flexible Metal-Organic Frameworks. ChemPhysChem 2011, 12, 247−258. (11) Sholl, D. S.; Lively, R. P. Defects in Metal-Organic Frameworks: Challenge or Opportunity? J. Phys. Chem. Lett. 2015, 6, 3437−3444. (12) Lawler, K. V.; Hulvey, Z.; Forster, P. M. On the importance of a precise crystal structure for simulating gas adsorption in nanoporous materials. Phys. Chem. Chem. Phys. 2015, 17, 18904−18907. (13) Woodley, S. M.; Catlow, R. Crystal structure prediction from first principles. Nat. Mater. 2008, 7, 937−946. (14) Odoh, S. O.; Cramer, C. J.; Truhlar, D. G.; Gagliardi, L. Quantum-Chemical Characterization of the Properties and Reactivities of Metal-Organic Frameworks. Chem. Rev. 2015, 115, 6051−6111. (15) Sharma, V.; Wang, C.; Lorenzini, R. G.; Ma, R.; Zhu, Q.; Sinkovits, D. W.; Pilania, G.; Oganov, A. R.; Kumar, S.; Sotzing, G. A.; Boggs, S. A.; Ramprasad, R. Rational design of all organic polymer dielectrics. Nat. Commun. 2014, 5, 4845. (16) Chandrasekhar, N.; Sholl, D. S. Large-Scale Computational Screening of Binary Intermetallics for Membrane-Based Hydrogen Separation. J. Phys. Chem. C 2015, 119, 26319−26326. (17) Armiento, R.; Kozinsky, B.; Fornari, M.; Ceder, G. Screening for high-performance piezoelectrics using high-throughput density functional theory. Phys. Rev. B: Condens. Matter Mater. Phys. 2011, 84, 014103. (18) Yan, F.; Zhang, X.; Yu, Y. G.; Yu, L.; Nagaraja, A.; Mason, T. O.; Zunger, A. Design and discovery of a novel half-Heusler transparent hole conductor made of all-metallic heavy elements. Nat. Commun. 2015, 6, 7308. (19) Nicholson, K. M.; Chandrasekhar, N.; Sholl, D. S. Powered by DFT: Screening Methods That Accelerate Materials Development for Hydrogen in Metals Applications. Acc. Chem. Res. 2014, 47, 3275− 3283. (20) Jain, A.; Ong, S. P.; Hautier, G.; Chen, W.; Richards, W. D.; Dacek, S.; Cholia, S.; Gunter, D.; Skinner, D.; Ceder, G.; Persson, K. A. Commentary: The Materials Project: A materials genome approach to accelerating materials innovation. APL Mater. 2013, 1, 011002. (21) Jain, A.; Hautier, G.; Moore, C. J.; Ong, S. P.; Fischer, C. C.; Mueller, T.; Persson, K. A.; Ceder, G. A high-throughput infrastructure for density functional theory calculations. Comput. Mater. Sci. 2011, 50, 2295−2310. (22) Hutter, J.; Iannuzzi, M.; Schiffmann, F.; VandeVondele, J. CP2K: atomistic simulations of condensed matter systems. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2014, 4, 15−25. (23) Nazarian, D.; Ganesh, P.; Sholl, D. S. Benchmarking density functional theory predictions of framework structures and properties in a chemically diverse test set of metal-organic frameworks. J. Mater. Chem. A 2015, 3, 22432−22440. (24) Goedecker, S.; Teter, M.; Hutter, J. Separable dual-space Gaussian pseudopotentials. Phys. Rev. B: Condens. Matter Mater. Phys. 1996, 54, 1703−1710. (25) Grimme, S. Semiempirical GGA-type density functional constructed with a long-range dispersion correction. J. Comput. Chem. 2006, 27, 1787−1799. (26) 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. (27) VandeVondele, J.; Hutter, J. Gaussian basis sets for accurate calculations on molecular systems in gas and condensed phases. J. Chem. Phys. 2007, 127, 114105.

quantitative use of computational screening of MOFs for adsorption and other applications.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.chemmater.6b04226. Computational details and analysis of atomic charges (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Randall Q. Snurr: 0000-0003-2925-9246 David S. Sholl: 0000-0002-2771-9168 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work has been supported by the Department of Energy Nanoporous Materials Genome Center, supported by the U.S. Department of Energy, Office of Basic Energy Sciences, Division of Chemical Sciences, Geosciences and Biosciences under Award DEFG02-12ER16362 and used resources of the Argonne Leadership Computing Facility, which is a DOE Office of Science User Facility supported under Contract DEAC02-06CH11357. This research was supported in part through the computational resources and staff contributions provided for the Quest high-performance computing facility (Project Allocation: p20663) at Northwestern University, which is jointly supported by the Office of the Provost, the Office for Research, and Northwestern University Information Technology. Y.G.C was supported by the Pusan National University Research Grant 2016. We gratefully acknowledge Christopher Knight for assistance with running calculations on MIRA and Cory Simon for helpful discussions.



REFERENCES

(1) Carrington, E. J.; Vitorica-Yrezabal, I. J.; Brammer, L. Crystallographic studies of gas sorption in metal-organic frameworks. Acta Crystallogr., Sect. B: Struct. Sci., Cryst. Eng. Mater. 2014, 70, 404− 422. (2) Allen, F. H. The Cambridge Structural Database: a quarter of a million crystal structures and rising. Acta Crystallogr., Sect. B: Struct. Sci. 2002, 58, 380−388. (3) Chung, Y. G.; Camp, J.; Haranczyk, M.; Sikora, B. J.; Bury, W.; Krungleviciute, V.; Yildirim, T.; Farha, O. K.; Sholl, D. S.; Snurr, R. Q. Computation-Ready, Experimental Metal-Organic Frameworks: A Tool To Enable High-Throughput Screening of Nanoporous Crystals. Chem. Mater. 2014, 26, 6185−6192. (4) Nazarian, D.; Camp, J. S.; Sholl, D. S. A Comprehensive Set of High-Quality Point Charges for Simulations of Metal-Organic Frameworks. Chem. Mater. 2016, 28, 785−793. (5) Simon, C. M.; Kim, J.; Gomez-Gualdron, D. A.; Camp, J. S.; Chung, Y. G.; Martin, R. L.; Mercado, R.; Deem, M. W.; Gunter, D.; Haranczyk, M.; Sholl, D. S.; Snurr, R. Q.; Smit, B. The materials genome in action: identifying the performance limits for methane storage. Energy Environ. Sci. 2015, 8, 1190−1199. (6) Qiao, Z.; Zhang, K.; Jiang, J. In silico screening of 4764 computation-ready, experimental metal−organic frameworks for CO2 separation. J. Mater. Chem. A 2016, 4, 2105−2114. 2527

DOI: 10.1021/acs.chemmater.6b04226 Chem. Mater. 2017, 29, 2521−2528

Article

Chemistry of Materials (28) Haldoupis, E.; Nair, S.; Sholl, D. S. Efficient Calculation of Diffusion Limitations in Metal Organic Framework Materials: A Tool for Identifying Materials for Kinetic Separations. J. Am. Chem. Soc. 2010, 132, 7528−7539. (29) Widom, B. Some topics in the theory of fluids. J. Chem. Phys. 1963, 39, 2808−2812. (30) Willems, T. F.; Rycroft, C.; Kazi, M.; Meza, J. C.; Haranczyk, M. Algorithms and tools for high-throughput geometry-based analysis of crystalline porous materials. Microporous Mesoporous Mater. 2012, 149, 134−141. (31) Dubbeldam, D.; Calero, S.; Ellis, D. E.; Snurr, R. Q. RASPA: molecular simulation software for adsorption and diffusion in flexible nanoporous materials. Mol. Simul. 2016, 42, 81−101. (32) Peng, D.; Robinson, D. B. A New Two-Constant Equation of State. Ind. Eng. Chem. Fundam. 1976, 15, 59−64. (33) Rappé, A. K.; Casewit, C. J.; Colwell, K.; Goddard, W.; Skiff, W. UFF, a full periodic table force field for molecular mechanics and molecular dynamics simulations. J. Am. Chem. Soc. 1992, 114, 10024− 10035. (34) Martin, M. G.; Siepmann, J. I. Transferable potentials for phase equilibria. 1. United-atom description of n-alkanes. J. Phys. Chem. B 1998, 102, 2569−2577. (35) Wilmer, C. E.; Farha, O. K.; Bae, Y. S.; Hupp, J. T.; Snurr, R. Q. Structure-property relationships of porous materials for carbon dioxide separation and capture. Energy Environ. Sci. 2012, 5, 9849−9856. (36) Potoff, J. J.; Siepmann, J. I. Vapor-liquid equilibria of mixtures containing alkanes, carbon dioxide, and nitrogen. AIChE J. 2001, 47, 1676−1682. (37) Xue, M.; Zhu, G. S.; Li, Y. X.; Zhao, X. J.; Jin, Z.; Kang, E.; Qiu, S. L. Structure, hydrogen storage, and luminescence properties of three 3D metal-organic frameworks with NbO and PtS topologies. Cryst. Growth Des. 2008, 8, 2478−2483. (38) Tong, X. L.; Wang, D. Z.; Hu, T. L.; Song, W. C.; Tao, Y.; Bu, X. H. Zinc and Cadmium Coordination Polymers with Bis(tetrazole) Ligands Bearing Flexible Spacers: Synthesis, Crystal Structures, and Properties. Cryst. Growth Des. 2009, 9, 2280−2286. (39) Zhong, D.-C.; Meng, M.; Deng, J.-H.; Luo, X.-Z.; Xie, Y.-R. Synthesis, crystal structure, and luminescent property of a threedimensional (3day) microporous Cd (II) metal-organic framework (MOF) based on mixed ligands. Inorg. Chem. Commun. 2011, 14, 1952−1956. (40) Sarkisov, L.; Martin, R. L.; Haranczyk, M.; Smit, B. On the Flexibility of Metal-Organic Frameworks. J. Am. Chem. Soc. 2014, 136, 2228−2231. (41) Yasodha, V.; Govindarajan, S.; Low, J. N.; Glidewell, C. Cationic, neutral and anionic metal(II) complexes derived from 4-oxo4H-pyran-2,6-dicarboxylic acid (chelidonic acid). Acta Crystallogr., Sect. C: Cryst. Struct. Commun. 2007, 63, M207−M215. (42) Gee, J. A.; Sholl, D. S. Effect of Framework Flexibility on C8 Aromatic Adsorption at High Loadings in Metal−Organic Frameworks. J. Phys. Chem. C 2015, 120, 370−376. (43) Haldoupis, E.; Watanabe, T.; Nair, S.; Sholl, D. S. Quantifying Large Effects of Framework Flexibility on Diffusion in MOFs: CH4 and CO2 in ZIF-8. ChemPhysChem 2012, 13, 3449−3452. (44) Watanabe, T.; Sholl, D. S. Accelerating applications of metal− organic frameworks for gas adsorption and separation by computational screening of materials. Langmuir 2012, 28, 14114−14128. (45) Deria, P.; Li, S.; Zhang, H.; Snurr, R. Q.; Hupp, J. T.; Farha, O. K. A MOF platform for incorporation of complementary organic motifs for CO2 binding. Chem. Commun. 2015, 51, 12478−12481. (46) Li, J. R.; Yu, J. M.; Lu, W. G.; Sun, L. B.; Sculley, J.; Balbuena, P. B.; Zhou, H. C. Porous materials with pre-designed single-molecule traps for CO2 selective adsorption. Nat. Commun. 2013, 4, 1538. (47) Llewellyn, P. L.; Bourrelly, S.; Serre, C.; Vimont, A.; Daturi, M.; Hamon, L.; De Weireld, G.; Chang, J. S.; Hong, D. Y.; Hwang, Y. K.; Jhung, S. H.; Ferey, G. High uptakes of CO2 and CH4 in mesoporous metal-organic frameworks MIL-100 and MIL-101. Langmuir 2008, 24, 7245−7250.

2528

DOI: 10.1021/acs.chemmater.6b04226 Chem. Mater. 2017, 29, 2521−2528