UvA-DARE (Digital Academic Repository) Understanding and solving disorder in the substitution pattern of amino functionalized MIL-47(V)

Electronic energies and elastic constants of four amino functionalized MIL-47(V) supercells were computed using plane wave density functional theory to determine the in ﬂ uence of the substituent positions on the organic linker. An inverse relationship between the ab initio energies and the elastic constants was found, indicating that the high electronic stability correlates with high mechanical stability. Torsion in all supercells was induced upon substitution, which caused strain in the NH 2 -MIL-47(V) supercell. The combined e ﬀ ect of the substituent bulkiness and the induced torsion reduced the pore volume of the NH 2 -MIL-47(V) structures by >7% and the surface area by >14% with respect to MIL-47(V). This reduction was con ﬁ rmed by lower saturation capacities of methane, CO 2 and benzene. When unfavourable substituent positions are chosen, large torsions caused a further reduction of the saturation capacity. Di ﬀ er-ences in surface area, pore volume and saturation capacity illustrate the importance of choosing the correct NH 2 -MIL-47(V) supercell.


Introduction
Metal-organic frameworks (MOFs), also known as coordination polymers, are a class of nanoporous materials consisting of metal clusters interconnected with organic linkers. Due to their large surface area, high porosity and enormous number of synthesizable topologies, MOFs have promising applications in hydrogen, 1,2 CO 2 3 and methane 4 storage, adsorptive separation 5-7 and catalysis. 8,9 In contrast to other porous materials, such as zeolites, the pore characteristics of MOFs can be easily functionalized by substituting one or more hydrogens of the organic linker for various functional groups. [10][11][12] This can be achieved by prefunctionalizing the organic linkers or by post-synthetic modification. Functionalization of MOFs can, for example, heavily affect the breathing properties as was shown by Devic et al. 10 for the narrow porelarge pore transition of functionalized MIL-53(Fe). An excellent example of a functionalized MOF was recently shown by Hu et al. 13 They constructed UTSA-100, an amino derivative of tetrazol-1,3-benzenedicarboxylic acid, that efficiently removes acetylene at room temperature from ethylene/acetylene mixtures containing only 1% acetylene. Another interesting example is NH 2 -MIL-88B(Fe) which improved the performance of multiple applications with respect to MIL-88B-(Fe), such as acting as a better photocatalyst for the reduction of Cr(VI) 14 and having a higher saturation capacity for ethanol, without losing the flexibility. 15 For NH 2 -MIL-53(Fe) it was found to be an excellent solid base catalyst for the transesterification of triglycerides with methanol 16 and Serra-Crespo et al. 17 showed that it can be used in a vacuum pressure swing adsorption process for natural gas upgrading.
As described above, amino functionalized MOFs can provide improved performance for various applications. From single X-ray structures it is possible to determine the appropriate unit cell of the functionalized MOF; the question remains however, how these substituted linkers are distributed with respect to each other? This is important for modeling these materials, since it is expected that different substitution patterns, resulting in linker-pair combinations, might lead to different adsorptive properties. For mono-substituted linkers there are numerous degrees of freedom with respect to the substituent positions, which are discussed below.
In this study, we address the substitution pattern of monoamino functionalized MIL-47(V). MIL-47(V) 18 is a vanadiumbased metallic cluster connected to BDC (= 1,4-benzenedicarboxylate) linkers and MIL-47(V) has the formula [V IV (O)BDC]. In contrast to the MIL-53 family, MIL-47(V) is missing the μ 2 -OH group due to thermal activation and does not show significant breathing. 19,20 The unit cell of MIL-47(V) is displayed in Fig. 1a.
Consider two MIL-47(V) unit cells along the a-axis that form a supercell; this results in adjacent BDC linkers, as shown in Fig. 1b. Only the left/right ortho/meta positions of the aromatic ring of the BDC-linker are available for mono-substitution since the ipso and para positions are connected to the carboxylate. To our knowledge, it is currently unclear if the substituent positions on both BDC-linkers are identical or that they differ. We expect that the most stable structures are those with low energetics and high mechanical stability.
By comparing the ab initio energies and the elastic constants of the different linker-pair combinations, we can determine the effect that linker substitution has on the overall supercell. Furthermore, we computed single component adsorption isotherms (methane, CO 2 and benzene) and the corresponding heat of adsorptions for the various isomeric supercells to illustrate the difference in adsorption properties (surface area and saturation capacity) due to these different linker-pair combinations.

Methodology
Approach Two MIL-47(V) unit cells, see Fig. 1b, each containing four BDC-linkers, were considered as a supercell to resemble two independent adjacent BDC-linkers. If only two adjacent BDClinkers are substituted, this is called partial substitution which shows the local effect of two adjacent linkers upon substitution. Full substitution implies that all BDC-linkers are substituted and will reflect the overall effect on the supercell. The proposed substitution patterns are shown in Fig. 2. A and B refer to the first and second BDC-linker, numbers indicate the positions of the aromatic ring with 1 = top ipso position and going clockwise until 6 = left ortho position.

Calculation of relative ab initio energies
To calculate the relative ab initio energies of the isomeric supercells, the geometry and lattice parameters of the structures were optimized using the Vienna Ab Initio Software Package 21,22 (VASP) with the Grimme dispersion corrected 23 PBE 24 exchange-correlation functional. A plane wave basis was employed with an energy density cut-off of 600 eV in order to achieve high accuracy, needed for the computation of the stress tensor. Core-electrons were described by the projected augmented wave method 25,26 and integration in k-space was performed over the Γ-point using the Monkhorst-Pack sampling scheme. 27 A Gaussian electron smearing of 0.05 eV was applied to converge the electronic structure. The geometric convergence criteria of the ionic forces were set to 10 −3 eV Å −1 .

Calculation of elastic constants
The mechanical and elastic properties of various MOFs have been studied intensively [28][29][30] and are essential to determine the stability of the supercells. A MOF subjected to a certain strain ε can be described by Hooke's law (within the elastic limit) as 31 with σ i being the stress tensor and C ij being the 6 × 6 elastic constant matrix in Voigt notation. ‡ For a low symmetry orthorhombic unit cell, such as MIL-47(V), the elastic constant matrix has the form with nine independent elastic constants. Dots indicate symmetry equivalent values. From the optimized structures, the elastic constants were computed using the same settings as described above. For the computation of the Hessian matrix, numerical derivatives were calculated using four displacements with the default value of 0.015 Å each. 32 VASP computes the elastic constants according to The first term of the r.h.s. of eqn (3), known as the clamped ion term (or the Born term) is calculated by the stress-strain relationship 33 for which six finite distortions (3 normal, 3 shear) of the lattice are applied, resulting in strained supercells. For each strain E, E + e1,…, E + e6, the corresponding stresses are computed ab initio. The clamped ionic elastic constants are then obtained by the least squares fitting procedure as described by Le Page and Saxe. 33 The Born expression is generalized by including the relaxation effects which are taken into account in the second term. The contributions due to relaxation of the elastic constants are calculated by multiplying the internal strain tensor Λ mj with the inverted Hessian matrix H −1 . Here, Λ mj is the internal strain tensor 34 given by with U is the energy per undeformed unit cell volume V 0 and r m being the displacement of the atom from the equilibrium position. By subtracting the elastic moduli from ionic relaxation from the clamped ion term, we obtain the elastic constants.
It should be noted that by substituting the BDC-linkers, the symmetry of the supercell is destroyed reducing the space group to P1. Therefore, calculation of the elastic constants is subjected to numerical inaccuracies. For this reason, we also computed the elastic constants for MIL-47(V) in the P1 space group to assess trends in the elastic constants.

Calculation of structural and adsorption properties
For the optimized structures we computed the pore volume, accessible surface area and single component adsorption isotherms using the molecular software package RASPA-2.0. 35 The void fraction is computed using the Widom particle insertion method 36 with helium as a probe at room temperature. Multiplying the helium void fraction by the total volume gives the pore volume.
The accessible surface area is calculated by the Shrake and Rupley method 37 in which a probe atom is used to 'roll' over the surface. It has been shown that the accessible surface area is a valid characteristic for describing the actual surface area of a MOF. 38 Values are close to surface areas obtained by the frequently used BET method. [38][39][40] Single component adsorption isotherms for methane, CO 2 and benzene were computed using Grand Canonical Monte Carlo 41 simulations. Here, the temperature (T ), chemical potential (μ) and volume (V) are fixed whereas the number of molecules are allowed to fluctuate. 42,43 We adopted the Con-figurational Bias Monte Carlo (CBMC) approach [44][45][46] where growth of the adsorbate molecules is biased towards favourable configurations. The amount of initialization cycles was set to 1000 and number of production cycles to 30 000. Fugacity coefficients and excess loadings were computed using the Peng-Robinson equation 47 of state. For statistical purposes, the simulation is divided into five intervals. For each interval the standard deviation within a 95% confidence interval was determined. The adsorbent interactions were described using the generic DREIDING 48 force field and adsorbates were described by the TraPPE force field. [49][50][51] Here, CO 2 is modelled as a tri-atomic molecule with partial atomic charges of 0.70e and −0.35e for C and O respectively to represent the quadrupole moment. 52 Framework charges 53 were computed using the REPEAT/ChelpG method 54,55 by fitting atomic charges on the ab initio electrostatic potential. The charge-charge interactions were computed using the Ewald summation 56 with a precision of 10 −6 and the van der Waals cut-off radius was set to 12.0 Å. Since MIL-47(V) does not show significant breathing upon adsorption, the framework atoms were kept rigid.
The heat of adsorption is computed using two methods: (i) the Widom insertion method 36,43 and (ii) the fluctuation theorem. 57,58 The heat of adsorption obtained from the Widom insertion method 36,59 is computed according to with 〈U hg 〉 the host-guest (MOF-adsorbate) energy, 〈U g 〉 the energy of the guest and 〈U h 〉 the energy of the host. In the gas phase, we assume ideal behaviour and replace the PV-term by RT with R being the gas constant and T the temperature. Since the framework and guest are considered rigid, 〈U g 〉 = 〈U h 〉 = 0. The (isosteric) heat of adsorption within the energy-fluctuation theorem in the grand canonical ensemble 60 is often applied to non-zero loading and is computed according to 57,58 ΔH ¼ hU ÁNi μ À hUi μ hNi μ hN 2 i μ À hNi μ hNi μ À hU g i À β À1 : with 〈··〉 being the average of either N, the amount of adsorbates or U, the internal energy, μ denotes the chemical potential and β the inverse temperature.

Lattice parameters and physical properties
The optimized lattice parameters (Table 1) show that the amino functionalized MIL-47(V) supercells have angles that are not exactly 90°. The MIL-47 family members have orthorhombic space groups with angles of 90 degrees. 11,61 This implies that various orientations of the unit cells and various linker-pair combinations give on average, in the long range, cell parameters of α = β = γ = 90°. The lattice axis a and lattice angle α are similar in their structures which can be attributed to the rigid metal-oxide bonds along the a-axis. Functionalizing MIL-47(V) reduces the pore volume and the effective surface area by >150 A 3 and 200 m 2 g −1 respectively. Two effects can be attributed to the pore volume and surface area reduction: the finite size of the substituent and the induced torsion due to the rotation of the BDC-linker described by the dihedral angles ϕ A and ϕ B ( Table 2). Although A2B3 has the smallest torsions, the cell parameters deviate more than those of A2B6 which has significantly larger dihedral angles (−14.20°and −10.88°vs. −2.21°a nd 2.5°).
Upon full substitution, configurations A2B2 and A2B5 show a decrease in the dihedral angles with respect to the partial substitution. This might imply that the full substitution reduces strain from the supercell for these configurations. Large dihedral angles of A2B6 are due to repulsive steric interactions caused by the amino groups and therefore result in significantly lower pore volume and effective surface area. For A2B3, the sign of the dihedral angles differ, indicating that the BDC-linkers are twisted in opposite directions. Fig. S1-A2B3 † shows indeed that the geometry between V atoms is convex.

Relative ab initio energies and elastic constants
Given the small energy difference (≤2.00 kJ mol −1 ) at 0 K, A2B2 and A2B5 are the most stable amino functionalized MIL-47(V) supercells for both partial and full substitution. At finite temperature, it is expected that from an energetic point of view, supercells A2B2 and A2B5 are equally likely. Despite the small dihedral angle of A2B3, the high relative ab initio energy indicates that this structure is energetically unfavourable. It is expected that this structure is strained due to the opposite positioning of the BDC-linkers. For A2B6 there is an unfavourable cooperative interaction upon substitution, since there is a factor 4.4 gain in ab initio energy going from partial to full substitution (Table 3). The amino groups of the adjacent linkers are in proximity to each other for A2B5, see Fig. 2 and S3. † Note that the linkers have a torsion and thus the nitrogens have different positions along the b-axis. Straining along the b-axis reduces the nitrogen-nitrogen distance resulting in a repelling interaction which gives a more rigid structure. This repelling nitrogen-nitrogen behaviour is not present in A2B2.
An inverse trend between the relative ab initio energies and the elastic constants can be observed. NH 2 -MIL-47(V) supercells that have high ab initio energies (A2B3 and A2B6) have lower elastic constants, whereas supercells that have low

Single component adsorption isotherm
The amino functionalized MIL-47(V) structures have lower surface areas than MIL-47(V) ( Table 1) and this can also be observed from the saturation capacity for methane and CO 2 as shown in Fig. 4 and for benzene in Fig. S2. † In particular for A2B6 this effect is more profound. At lower pressures the CO 2 isotherms differ, whereas for methane the isotherms are identical. For CO 2 this difference is expected due to the dominating sorbent-adsorbate interaction between the quadrupole moment of CO 2 and the amino group. Comparing heat of adsorptions at infinite dilution Table 4 shows that the simulated trend agrees with the experimental trend: larger heat of adsorption for the amino functionalized MIL-47(V). A similar behaviour has been observed between MIL-53(Al) and NH 2 -MIL-53(Al). 64 At finite loading, it is also observed, see Fig. S4, † that the heat of adsorption is more favourable for NH 2 -MIL-47(V) than for MIL-47(V) for methane and CO 2 adsorption. At higher loadings, the heat of adsorption changes due to adsorbate-adsorbate interactions (Fig. S5 †). The heat of adsorption for CO 2 is higher than for methane as expected, although there is a slight difference between experimental and computational results. This difference is attributed to the rigidity of the force field framework; upon CO 2 adsorption, the amino group is no longer in the same plane as the BDC-linker. Optimized VASP structures from Leus et al. 61 show a dihedral angle of −15.20°with the hydrogens pointing away. The heat of adsorption difference among the various supercells is negligible and is not surprising since CO 2 will adsorb most-likely on the nitrogen of the amino group which is irrelevant of the substitution pattern of the linkers. The combined effect of the substituted linker can alter the adsorption positions at higher loadings.
At saturation loadings the isotherms change for both adsorbates. Leus et al. 61 measured experimentally the single component isotherms of methane and CO 2 at T = 303 K for MIL-47   (V) and NH 2 -MIL-47(V) with a saturation capacity for CO 2 and CH 4 in MIL-47(V) of 7.7 and 4.1 mol kg −1 respectively. The simulated saturation capacities are higher than the experimental ones, which is a general feature since for simulated isotherms, defects and impurities are not taken into account. For larger adsorbates, the influence of the substitution patterns is more significant. Fig. S2 † shows the computed single component adsorption isotherms for benzene at 328 K. Steeper uptake of the benzene adsorption isotherms for NH 2 -MIL-47(V) is observed and indicates that benzene adsorbs stronger in NH 2 -MIL-47(V) than in MIL-47(V). This is not a confinement effect, which has been observed for benzene adsorption in zeolite MFI, 65 since benzene is approximately 5.5 Å in size and the pore diameter of MIL-47(V) is 10.5 and 11.0 Å. 66 However, at high pressures, the saturation capacity of MIL-47(V) for benzene is higher than that for NH 2 -MIL-47(V) due to the better packing of the benzene adsorbates.

Conclusions
We demonstrated that various linker-pair combinations of amino functionalized MIL-47(V) affect structural and adsorption properties of the sorbent. Functionalizing MIL-47(V) with amino groups reduces the accessible surface area, pore volume and the saturation capacities for methane, CO 2 and benzene. Heat of adsorption at infinite dilution of various linker-pair combinations is equal, indicating that different substitution patterns do not affect adsorption in the Henry regime. As a result of linker substitution, torsions of the BDC-linker induce strain in the material leading to both energetically and mechanically unstable structures. An inverse relationship between the relative ab initio energy and elastic constants was observed. It is therefore crucial to consider an appropriate unit cell in simulation studies. For amino functionalized MIL-47(V), A2B2 and A2B5 are the most stable structures.