Improved Predictive Tools for Structural Properties of Metal-Organic Frameworks.

The accurate determination of structural parameters is necessary to understand the electronic and magnetic properties of metal–organic frameworks (MOFs) and is a first step toward accurate calculations of electronic structure and function for separations and catalysis. Theoretical structural determination of metal-organic frameworks is particularly challenging because they involve ionic, covalent, and noncovalent interactions, which must be treated in a balanced fashion. Here, we apply a diverse group of local exchange-correlation functionals (PBE, PBEsol, PBE-D2, PBE-D3, vdW-DF2, SOGGA, MN15-L, revM06-L, SCAN, and revTPSS) to a broad test set of MOFs to seek the most accurate functionals to study various structural aspects of porous solids, in particular to study lattice constants, unit cell volume, two types of pore size characteristics, bond lengths, bond angles, and torsional angles). The recently developed meta functionals revM06-L and SCAN, without adding any molecular mechanics terms, are able to predict more accurate structures than previously recommended functionals, both those without molecular mechanics terms (PBE, PBEsol, vdW-DF2, and revTPSS) and those with them (PBE-D2 and PBE-D3). To provide a broader test, these two functionals are also tested for lattice constants and band gaps of unary, binary, and ternary semiconductors.

Kohn-Sham density functional theory (KS-DFT) is the most popular way to calculate structural parameters of solids [15][16][17][18][19]. The accuracy of a KS-DFT calculation of a system depends on the choice of exchange-correlation (XC) functional. Our study builds on previous work by Sholl and co-workers [20], who studied a set of 12 MOFs containing various metal atoms. Of six functionals tested, they found that the PBE-D2, PBE-D3, and vdW-DF2 functionals are most promising for the accurate determination of the structural parameters of the MOFs (lattice parameters, unit cell volume, bonded parameters, and pore descriptors). Because improved functionals have become available since the work of Sholl et al., the present work considers the same set of 12 MOFs (listed in Table 1) studied by Sholl et al. [20] and two other MOFs (UiO-66 and MOF-5 or IRMOF-1, which were added to increase the diversity of the test set) with ten functionals-the three that they found to perform best and seven others. We

•
GA functionals: gradient approximations in which the energy density depends on the local spin-specific electron densities and their gradients. In particular, we consider PBE [21,22], PBEsol [23], and SOGGA [24] in this group. • Meta functionals in which the energy density also depends on the spin-specific local kinetic energy densities. In particular, we consider MN15-L [25], revM06-L [26], SCAN [27], and revTPSS [28] in this group. • GA+MM functionals in which a GA is augmented with molecular mechanics terms. In particular, we consider PBE-D2 and PBE-D3 in this group [29,30]. • A functional with nonlocal correlation, vdW-DF2 [31], which has local exchange and nonlocal correlation.
Because nonlocal Hartree-Fock exchange is very expensive in plane wave codes, most calculations on solids are restricted to using only local exchange (sometimes called semilocal exchange) in which the exchange energy density at a point is approximated in terms of functions evaluated at that point, and both the work of Sholl and coworkers and the present work consider only functionals with local exchange. All functionals considered here except vdW-DF also have local correlation, but an efficient algorithm [32] is available for vdW-DF so the cost is not impractically higher than that of the local functionals. In particular, we found that the relative costs for the calculations reported here with the different kinds of functionals are approximately: GA, 1; GA+MM, 1.2, meta, 1.5-2, and vdW-DF2, 3.5). a Chemical formula of the unit cell of the MOF; b The reference code associated with each structure in the Cambridge Structural Database (CSD) and CoRE MOF database [47]; c Experimental reference. Table 1 includes the reference code from the Cambridge Structural Database (CSD) and CoRE MOF database [47] of the crystal structure associated with each MOF in the test set. There are nine different metal centers among the twelve MOFs; this includes two lanthanides (Dy and Sm). We considered two different oxidation states of Cu (Cu +1 and Cu +2 ) and Fe (Fe +2 and Fe +3 ) (Table 1), and both magnetic and non-magnetic MOFs are considered in the study.
We also note the work of Tran et al. [48,49], who studied lattice constants, cohesive energies, and bulk moduli of strongly and weakly bound solids and found the best behavior for rev-VDW-DF2 and various meta functionals. The reader is directed to their studies for information complementary to the present study.

Computational Methods
We performed all the density functional theory calculations with the latest version (5.4.4) of the Vienna Ab initio Simulation Package (VASP) [50,51] with periodic boundary conditions. The core-valence electron interactions were described by using the projector augmented wave (PAW) method [52]. The kinetic energy cutoff was 520 eV for all calculations (the same as used by Sholl and coworkers [20]). All the lattice constants, cell shapes, and atomic coordinates were fully relaxed during the optimizations. The C 6 and R 0 parameters for PBE-D2 method for elements H-Xe are the original values [29] where available, but for Dy and Sm we used C 6 = 33.98 Jnm 6 /mol and R 0 = 2.226 Å, as used by Sholl and co-workers [20]. We set a convergence criterion of 10 −5 eV for electronic energy minimization and a convergence criterion of 0.0002 eV/Å on the force on each atom. Monkhorst-Pack k-point grids [53] were used as specified in Table S1 (tables with a prefix S are in the Supplementary Materials).
The Cu 3 H 4 C 10 O 10 (MURCEH) MOF requires additional considerations. We used a 2 × 1 × 1 supercell because the unit cell contains an odd number of metal centers (3 Cu centers); the supercell is used to make an even number of metal centers (6 Cu centers) in order to study the magnetic ground state. All the structural parameters of MURCEH are calculated based on the 2 × 1 × 1 supercell.
We performed spin-polarized KS-DFT calculations to determine the magnetic ground states. The initial guesses of the number of unpaired electrons on metal centers and the magnetic ground states are listed in Table S1. The magnetic properties of the MOFs are determined based on the number of unpaired electrons in the metal center of the ground state and are listed in Table S1. We used the PBEsol and PBE-D3 functionals to verify the magnetic ground states, and the magnetic ground states of the MOFs determined in our calculations agree well with the experimental results, as well as with the work of Sholl and co-workers [20].
Using the Zeo++ software [54][55][56], we calculated two kinds of pore sizes, in particular the pore limiting diameter (PLD) and the largest cavity diameter (LCD). The PLD is defined as the diameter of the largest sphere that can move through the MOF. The LCD is defined as the diameter of the largest sphere that can fit in the pore of the MOFs. According to the Zeo++ software website [57], the LCD and PLD are designated as "Di" and "Dif".
To assess the performance of the functionals on the various properties, we calculated the mean signed error (MSE), mean unsigned error (MUE), and mean unsigned percentage error (MUPE) for each combination of functional and property. The errors of the structural parameters (lattice parameters, lattice angles, bond lengths, bond angles, and torsional angles) are defined as the deviations between experimental results and the results calculated by DFT.
We also studied unary, binary, and ternary semiconductors with two of the functionals (revM06-L and SCAN). For these calculations, we used a kinetic energy cutoff of 520 eV, a 12 × 12 × 12 Monkhorst-Pack k-point mesh, and the same convergence criteria and force criteria as for the MOFs.

Lattice Constants and Unit Cell Volume
The first and most important structural parameters of a solid are the lattice constants. The three lattice constants are designated as a, b, and c (see Figure 1). To ensure that all the crystals are of high quality, the experimental lattice constants of the 14 MOFs are from solvent-free and disorder-free single-crystal X-ray diffraction (XRD) data at room temperature (290-310 K) with an R-factor (the measure of the agreement between the experimental X-ray diffraction data and the crystallographic model, also known as residual factor, reliability factor, R-value, or RWork) of less than 10 [58]. Figure 2a shows the error statistics of the various functionals for predicting lattice constants. The MSEs show that on average, PBE, PBE-D3, vdW-DF2, and revTPSS overestimate the lattice constants, and PBEsol, PBE-D2, SOGGA, MN15-L, revM06-L, and SCAN underestimate them. The SCAN (0.06 Å) and revM06−L (0.07 Å) functionals are the most accurate (i.e., have the smallest MUE). The PBE functional has the largest MUE, which is consistent with what is usually found for PBE for other kinds of solids.
Full details given in the Supplementary Materials show that the maximum deviation of the PBE functional is for Ag 4 C 12 Cl 4 O 8 MOF (RORQOE) (a = 0.50 Å, b = 0.28 and c = 0.98 Å), which also agrees well with the results of Sholl and co-workers [20]. However, SCAN and revM06-L have errors in the lattice constants of this MOF that are 5-24 times smaller. SCAN and revM06-L are very accurate for the lattice constants of Cd12H48C72N72O48 MOF (GUPCUQ01); the deviation is only 0.01 Å for both functionals (Table S2). We find that SOGGA and PBEsol are almost equally good for predicting the lattice constants of the MOFs, with MUEs of 0.11 Å for PBEsol and 0.12 Å for SOGGA, although PBEsol has an empirical parameter adjusted to lattice constants, and SOGGA does not. The revM06-L and SCAN functionals are even better and have the lowest MUEs of the functionals studied here.
We also studied the changes in lattice angles of the unit cells of the test sets of MOFs. We found that there is no effect of the functionals on the lattice angles for seven of the MOFs, but there are some changes in lattice angles for the rest of the MOFs. Based on those changes we have calculated the MSE, MUE, and MUPE (given in Table S8) statistics, and these values show that SCAN and revM06-L perform best for this property; both have an MUPE of only 0.3%. To ensure that all the crystals are of high quality, the experimental lattice constants of the 14 MOFs are from solvent-free and disorder-free single-crystal X-ray diffraction (XRD) data at room temperature (290-310 K) with an R-factor (the measure of the agreement between the experimental X-ray diffraction data and the crystallographic model, also known as residual factor, reliability factor, R-value, or R Work ) of less than 10 [58].  (Table S2). We find that SOGGA and PBEsol are almost equally good for predicting the lattice constants of the MOFs, with MUEs of 0.11 Å for PBEsol and 0.12 Å for SOGGA, although PBEsol has an empirical parameter adjusted to lattice constants, and SOGGA does not. The revM06-L and SCAN functionals are even better and have the lowest MUEs of the functionals studied here.
We also studied the changes in lattice angles of the unit cells of the test sets of MOFs. We found that there is no effect of the functionals on the lattice angles for seven of the MOFs, but there are some changes in lattice angles for the rest of the MOFs. Based on those changes we have calculated the MSE, MUE, and MUPE (given in Table S8) statistics, and these values show that SCAN and revM06-L perform best for this property; both have an MUPE of only 0.3%.  Table S3.

Pore Size
The accurate prediction of pore diameters is critical for understanding separations and catalytic selectivity. Here, we have studied two different kinds of pore diameters of the MOFs; the LCD and PLD nomenclature (explained in Section 2) is taken from the work of Haldoupis and coworkers [59].
We found that PBE, PBEsol, SOGGA, and revTPSS overestimate both measures of pore size (   Table S3.

Pore Size
The accurate prediction of pore diameters is critical for understanding separations and catalytic selectivity. Here, we have studied two different kinds of pore diameters of the MOFs; the LCD and PLD nomenclature (explained in Section 2) is taken from the work of Haldoupis and coworkers [59].
We found that PBE, PBEsol, SOGGA, and revTPSS overestimate both measures of pore size (

Bond Lengths
Bond length is the distance between two bonded atoms. We only considered the bond lengths of the inorganic nodes. We considered 60 bond lengths (Table S5) (Table  S5). The overestimated bond lengths by PBE and vdW-DF2 are consistent with their overestimation of the lattice constants (Table S2). The Cu−N bond in Cu8H8C8N12Cl8 (QEJZUB01) is highly overestimated by 0.23 Å by PBE, as also found by Sholl et al. [20]. The Fe−O bond in Fe4H4C4O12 (DEMLIR) also has large errors for some functionals: 0.23 and 0.18 Å by MN15-L and SOGGA respectively, and this large deviation of the Fe−O bond lengths is consistent with the high lattice constant errors for this system (Table S2). Overall, revM06-L has the smallest MUE (0.026 Å) and SCAN has the second smallest error (0.027 Å).

Bond Lengths
Bond length is the distance between two bonded atoms. We only considered the bond lengths of the inorganic nodes. We considered 60 bond lengths (Table S5) (Table S5). The overestimated bond lengths by PBE and vdW-DF2 are consistent with their overestimation of the lattice constants (Table S2). The Cu-N bond in Cu 8 H 8 C 8 N 12 C l8 (QEJZUB01) is highly overestimated by 0.23 Å by PBE, as also found by Sholl et al. [20]. The Fe-O bond in Fe 4 H 4 C 4 O 12 (DEMLIR) also has large errors for some functionals: 0.23 and 0.18 Å by MN15-L and SOGGA respectively, and this large deviation of the Fe-O bond lengths is consistent with the high lattice constant errors for this system (Table S2). Overall, revM06-L has the smallest MUE (0.026 Å) and SCAN has the second smallest error (0.027 Å).

Broader Test on Semiconductors
The results presented above show that the revM06-L and SCAN functionals give the most accurate results for the structural parameters of MOFs. Since the performance of these two functionals for structural properties of MOFs is similar, we carried out a broader test of these functionals for important crystalline properties, in particular for lattice parameters and band gaps of ten semiconductors, namely Si, Ge, SiC, β-GaN, CdS, CdSe, ZnS, ZnO, TiO 2 -anatase, and TiO 2 -rutile. The results are in Figure 6 and Table S8. This comparison shows that revM06-L and SCAN both are accurate to better than 0.7% for predicting the lattice constants of the semiconductors, but revM06-L is much more accurate for predicting the band gaps of the semiconductors.
Molecules 2020, 25, x FOR PEER REVIEW 9 of 13 functionals for important crystalline properties, in particular for lattice parameters and band gaps of ten semiconductors, namely Si, Ge, SiC, β-GaN, CdS, CdSe, ZnS, ZnO, TiO2-anatase, and TiO2rutile. The results are in Figure 6 and Table S8. This comparison shows that revM06-L and SCAN both are accurate to better than 0.7% for predicting the lattice constants of the semiconductors, but revM06-L is much more accurate for predicting the band gaps of the semiconductors. Previous work in the literature has shown that the band gaps of semiconductors are most accurately determined by hybrid functionals (e.g., HSE06 [60]) and by high-local-exchange functionals (e.g., HLE17 [61,62]). For large systems like MOFs, HLE17 is as accurate as HSE06 with almost 100 times less computational cost [61,62], and for this reason we focus here on local functionals. The previous study showed that HLE17 is far more accurate for predicting the band gaps of unary and binary semiconductors (Si, Ge, SiC, β-GaN, CdS, CdSe, ZnS, and ZnO) than local functionals like PBE, PBE+U with U set to 4.0 eV, PBEsol, B3LYP, TPSS, GAM, and HCTH/407, but the lattice constants are not accurate [61]. A similar trend was found in our recent paper where HLE17 gave more accurate band gaps for TiO2 (both anatase and rutile) than other local functionals like PBE, PBEsol, and PBE+U (U = 4.0 eV) [62], but HLE17 did not produce good structural parameters. In contrast, the present study shows that the revM06-L and SCAN local functionals are both able to produce good structural parameters. There is, however, still a trade-off in that neither of these functionals produces band gaps as accurate as those calculated by HLE17. The key finding of the present section, though, is that the trade-off is less severe with revM06-L than with SCAN, because revM06-L yields more accurate band gaps than SCAN.

Conclusions for MOFs
We studied the structural parameters-lattice parameters, unit cell volume, bond lengths, bond angles, torsional angles, and pore sizes (largest cavity diameters and pore limiting diameters) of a diverse set of 14 MOFs with nine local functionals and one nonlocal functional. The mean unsigned percentage errors (MUPEs) for the various properties are summarized in Table 2. The last column of the table averages these MUPEs to get an average MUPE (AMUPE) to provide a more global summary of the present findings. The three functionals previously [20] found to be most accurate have AMUPEs of 3.0-3.5%. Five of the functionals introduced in the present study have smaller AMUPEs than that, namely SCAN (1.5%), revM06-L (1.7%), PBEsol (2.4%), SOGGA (2.5%), and revTPSS (2.5%). Thus, these functionals, especially SCAN and revM06-L, are highly recommended for future work on the computational design of MOFs. Previous work in the literature has shown that the band gaps of semiconductors are most accurately determined by hybrid functionals (e.g., HSE06 [60]) and by high-local-exchange functionals (e.g., HLE17 [61,62]). For large systems like MOFs, HLE17 is as accurate as HSE06 with almost 100 times less computational cost [61,62], and for this reason we focus here on local functionals. The previous study showed that HLE17 is far more accurate for predicting the band gaps of unary and binary semiconductors (Si, Ge, SiC, β-GaN, CdS, CdSe, ZnS, and ZnO) than local functionals like PBE, PBE+U with U set to 4.0 eV, PBEsol, B3LYP, TPSS, GAM, and HCTH/407, but the lattice constants are not accurate [61]. A similar trend was found in our recent paper where HLE17 gave more accurate band gaps for TiO 2 (both anatase and rutile) than other local functionals like PBE, PBEsol, and PBE+U (U = 4.0 eV) [62], but HLE17 did not produce good structural parameters. In contrast, the present study shows that the revM06-L and SCAN local functionals are both able to produce good structural parameters. There is, however, still a trade-off in that neither of these functionals produces band gaps as accurate as those calculated by HLE17. The key finding of the present section, though, is that the trade-off is less severe with revM06-L than with SCAN, because revM06-L yields more accurate band gaps than SCAN.

Conclusions for MOFs
We studied the structural parameters-lattice parameters, unit cell volume, bond lengths, bond angles, torsional angles, and pore sizes (largest cavity diameters and pore limiting diameters) of a diverse set of 14 MOFs with nine local functionals and one nonlocal functional. The mean unsigned percentage errors (MUPEs) for the various properties are summarized in Table 2. The last column of the table averages these MUPEs to get an average MUPE (AMUPE) to provide a more global summary of the present findings. The three functionals previously [20] found to be most accurate have AMUPEs of 3.0-3.5%. Five of the functionals introduced in the present study have smaller AMUPEs than that, namely SCAN (1.5%), revM06-L (1.7%), PBEsol (2.4%), SOGGA (2.5%), and revTPSS (2.5%). Thus, these functionals, especially SCAN and revM06-L, are highly recommended for future work on the computational design of MOFs.