Prediction of Multiple Hydrogen Ligation at a Vanadium(II) Site in a Metal-Organic Framework

Densifying hydrogen in a metal-organic framework (MOF) at moderate pressures can circumvent challenges associated with high pressure compression. The highly tunable structural and chemical composition in MOFs affords vast possibilities to optimize binding interactions. At the heart of this search are the nanoscale characteristics of molecular adsorption at the binding site(s). Using density functional theory (DFT) to model binding interactions of hydrogen to the exposed metal site of cation-exchanged MFU-4 l , we predict multiple hydrogen ligation of H 2 at the first coordination sphere of V(II) in V(II)-exchanged MFU-4 l . We find that the strength of this binding between the metal site and H 2 molecules can be tuned by altering the halide counterion adjacent to the metal site and that the fluoride containing node affords the most favourable

interactions for high density H 2 storage.Using energy decomposition analysis, we delineate electronic contributions that enable multiple hydrogen ligation and demonstrate its benefits for hydrogen adsorption and release at modest pressures.
A major challenge for a future hydrogen-based economy is the need for energy intensive compression and liquefaction of H 2 to counteract its low volumetric energy density. 1,2Current machinery for hydrogen distribution for vehicular storage involves compression above 350 bar, 3 where the majority of the operating cost at the refuelling station is due to the compressor. 4In addition to environmental and energy costs of such high-pressure treatment, storage of a highly compressed flammable gas is a safety concern for vehicle operators.Metalorganic frameworks (MOFs) offer a potential solution to this problem as a broad class of new materials that can, if appropriately designed, adsorb and release H 2 at more moderate pressures. 5 a materials design problem, optimization of MOFs for hydrogen storage is an enormous challenge, with a complex interplay between constraints imposed by the desired functional characteristics, and non-trivial synthetic considerations.The high surface to volume ratio in a MOF combined with its structural and chemical tunability result in combinatorial numbers of possibilities for the design of an optimal storage material. 5,6Nonetheless, the syntheses and measurements of new MOF combinations are time-consuming and resource-intensive, prompting the need for computationally aided approaches.Reliable quantum chemical calculations can screen candidate materials to identify promising leads, at a minimum.Design principles for high-density hydrogen storage can also be derived from computational prediction of the electronic structure of hydrogen interacting with the MOF binding site and the associated pore. 7In-silico characterizations of the feature space in MOFs can identify the leading contributors to interactions and thereby offer a pathway to release H 2 below operating pressures of 100 bar.
In the list of technical targets for vehicular storage, the US Department of Energy (DOE) has earmarked volumetric storage capacity at 40 g/L by the year 2025 at operating condi-tions in the range of [−40, 60] • C and 5 to 12 bar. 8The current record for physisorptive storage is held by Ni 2 -m-dobdc (dobdc=1,4-dioxido-1,3-benzenedicarboxylate) with a volumetric capacity of 11 g/L at 25 • C for a pressure swing between 100 and 5 bar. 9Evidently, physisorption-driven storage alone does not suffice for room temperature storage applications, and stronger binding interactions are desirable. 5,6,10,11Given the demand for high density storage of hydrogen in the solid state, 12,13 it is highly desirable to bind multiple hydrogens at an open metal site in a MOF with binding enthalpies in the 15-25 kJ/mol range though this has not yet been demonstrated.[16][17] Strong H 2 -binding to a five-coordinate V(II) center in V 2 Cl 2.8 (btdd), which displays an isosteric heat of 21 kJ/mol, 18 provides promise for ambient temperature storage of hydrogen.
Guided by this insight, we wanted to explore analogous H 2 binding to V(II) in a fourcoordinate scaffold with two available coordination sites, which could potentially be accessed within MFU-4l.Toward this objective, we predict feasible ligation of two hydrogen molecules to V(II) in V(II) exchanged MFU-4l which has been previously prepared and characterized in a hydrated form. 19Furthermore the parent material, MFU-4l has the benefit of being able to readily tune the peripheral anion, which can modulate the binding energetics of H 2 .
The parent MFU-4l framework is highly robust and demonstrates exceptional thermal stability. 20,21The defining feature of MFU-4l is scorpionate-type coordination with BTDD 2− ligands (( ) at these tetrahedral metal sites, where a tridentate ligand speared with nitrogens binds the metal in a fac manner.The term scorpionate 22 derives from two donor sites grasping at the metal like the pincers of a scorpion followed by a third and final donor, the stinger, reaching over the plane formed by the metal site and the two donors.4][25][26] This material piqued our interest for two reasons: i) it features 4 outer metal sites per node that have tetrahedral coordination, and ii) metal-exchange for V(II), which is known to favour octahe- dral coordination, has been reported 19 for these sites.Figure 1a shows the cubic unit cell of Zn(II)MFU-4l.Here we consider each node to behave independently.The pentanuclear cluster model shown in Fig. 1b was built from the parent scaffold.Binding environment around the V(II) site was modeled by truncation of the node at its benzotriazolate extremities and a triazolate truncation of the linker was employed elsewhere giving a molecular formula of V(II) 1 Zn(II) 4 Cl 4 (btz) 3 (tz) 3 (btz = benzotriazolate, tz = triazolate) to the cluster.An octahedral Zn(II) at the center of the node is surrounded by 4 peripheral tetrahedral sites, one of which is substituted with V(II).Following experimental characterization of H 2 binding to V(II) in a similar ligand-field environment in V 2 Cl 2.8 (btdd) 18 we model V(II) in its high-spin (S=3/2) ground state.
The MOF crystal is periodic and three dimensional, and upon approach of a guest molecule to a metal site, there may be strong binding contributions reflecting local chemical interactions at the site, supplemented by weaker long-range dispersion interactions.To ensure accurate description of the critical local chemistry, we employ the ωB97M-V functional, 27 which is an accurate range-separated hybrid metaGGA that includes non-local VV10 dispersion, 28 designed by the combinatorial "survival of the most transferable" protocol.ωB97M-V was the top-performing hybrid density functional in several large assessments, including the MGC84 database, 29 the large and diverse GMTKN55 benchmark dataset, [30][31][32] and the TMC151 transition metal database. 33These findings are buttressed by a recent extensive benchmark 34,35 of density functionals for hydrogen storage.Employing ωB97M-V limits calculations to a cluster node, as exact exchange is not feasible for periodic calculations.We have chosen our cluster model to be very faithful to the full MOF for capturing metal-ligand interactions at short distances, and adequate for capturing the leading contributions from comparatively long-range dispersion forces that prevail in MOF pores.Geometries were optimized in the def2-SVP basis that contains f-polarization at the transition metal site, and p polarization at the hydrogens.The def2-ECP effective core potential was used for elements heavier than Rb (Z=37) such as iodine. 36Single point calculations were performed with this basis augmented with the fairly large Karlsruhe def2-TZVPPD basis 37,38 in the binding region with 2f and 1g polarization functions at the metal site and 3p and 1d polarization functions for hydrogen.Geometries were converged to 3 × 10  [40][41][42] Here the MOF node and its ligated hydrogen(s) are partitioned into two interacting fragments.Separate single-point calculations were performed on them followed by a block diagonal reassembly of MO coefficients to give the frozen density matrix and energy (∆ Frz).Polarization is described using a basis of fragment electric-field response functions (FERFs) which allows a fragment orbital to adjust to the weak electric field of its counterpart resulting in dipolar and quadrupolar responses. 43Finally, the charge transfer term ∆ CT is computed from the difference between the completely relaxed and polarized electron density.Preparation energy (∆ Prep) is the energy required to distort the unbound MOF node and free hydrogen into their bound configuration.∆ E is the sum of all electronic contributions to binding and ∆ H evaluated at 77 K includes the effect of vibrations of nuclear modes at finite temperature.Columns I and II in  to engage in backbonding where the effective range is shorter than that of the σ bond. 44,45e pentacoordinate V(II) in MFU-4l is distorted trigonal bipyramidal (TBP) and transitions to a distorted octahedron on loading of the second hydrogen.To explore this further we used a perturbative approach 40 to separate forward from back-donation.Donation from the hydrogen σ orbital to the empty metal d z 2 is the dominant charge transfer contributor offering -23.4 and -26.1 kJ/mol of stabilization for the first and second hydrogen loading.At shorter distances, favorable π overlap is afforded by the metal d π backdonation to hydrogen σ * resulting in sizeable contributions of -6.6 and 6.7 kJ/mol for binding the first and the second H 2 .Strikingly, the presence of the first bound H 2 causes forward donation from the second H 2 to be strengthened with the back-donation remaining nearly identical.
Complementary Occupied-Virtual Pairs (COVPs) constructed from unitary transformations within the occupied and virtual subspaces of the two fragments (the MOF node and hydrogen) provide a physical picture for orbitals engaged in forward-and back-bonding.
While the set of all COVPs exactly describes charge transfer between the fragments, the major energetic (and charge) contributions to bonding arises from one or two significant COVPs. 40Charge transfer orbitals for V(II) to H 2 backbonding in Fig. 2b plotted with an isodensity contour of 0.06 shows d π orbitals of V(II) engaging with the hydrogen σ * orbital.
A closer inspection of the back-donating orbital at the MOF node reveals contribution from the p z orbital on the chloride ligand bound to V(II) (see Fig. 2b).This motivated an investigation of the effect of halide substitution on the binding strength of hydrogen with the MOF node.Table 1 lists the binding energies for the sequential loading of two hydrogens onto the V(II)MFU-4l node represented as a V(II)L 3 X complex containing an X -halide bound to the metal.Among different halides, an F -was found to promote the most favorable H 2 binding energies.Heavier halide ligands, especially in the case of Br -, show strong comparative stabilization on binding the second H 2 .This effect is driven by a stronger charge transfer term dominated by forward donation from H 2 to V(II).Thus, the halide ligand can be used to finely tune the strength of multiple H 2 binding with V(II) in this material.
Energy lowering upon H 2 binding is countered by the entropic tendency of the bound hydrogen to gain translational and rotational freedom.The free energy change on binding describes the balance between these two terms and determines the adsorption isotherm.We estimate the enthalpy and Gibbs free energy change of binding muliple hydrogens to the node with a standard frequency analysis at the B3LYP-D2/def2-SVP/RRHO level of theory and employ anharmonic 1D corrections to zero point vibrational energies (ZPVEs) for the frustrated modes of H 2 .Low frequency (≤ 100 cm −1 ) normal modes were raised to 100 cm −1 for evaluation of the enthalpic and entropic terms to correct for the well-known breakdown of the harmonic oscillator model for soft vibrations. 46The binding of each hydrogen adds 5 modes to the system due to the frustrated translations and rotation of bound H 2 on the MOF surface.Anharmonic 1D corrections to the ZPVEs of these modes were computed by employing the Fourier Grid Hamiltonian 47 to solve for the vibrational spectrum by sampling with 1D PES with finite displacements scans around the equilibrium binding configuration with the framework held rigid.Details pertaining to anharmomic corrections to the ZPVE are included in Section I of the SI.A scaling factor of s = 0.9893 obtained from quartic force fields was added to the B3LYP-D2 ZPVEs to account for anharmonic coupling between vibrational modes. 48ble 2 lists changes in binding energy, enthalpy, and free energies for the sequential binding of two hydrogens to V(II)MFU-4l at the temperature of liquid nitrogen (77 K) and 1 atm pressure.Thermal excitation of soft framework modes (< 100 cm −1 ) are dampened at such low temperatures and the entropic penalty of binding is small compared to the Table 1: Binding energy and its decomposition for ligation of two hydrogens to the V(II)MFU-4l node represented as V(II)L 3 X cluster model with L standing for the chelating nitrogens of the scorpionate SBU and X -as the placeholder for the halide counterion.All contributions are in kJ/mol.2: Binding Energy, Enthalpy, and Free Energy at 77 K and 1 atm for H 2 binding to V(II)XMFU-4l, X = F -, Cl -, Br -, I -.∆∆ values correspond to changes on second hydrogen loading (∆ 2 − ∆ 1 ).All values are in kJ/mol.
The binding curve, determined by the average number of gas molecules v = ⟨i⟩ bound as a function of partial pressure p at temperature T , is given by which is the binding polynomial that sums over the ligation states V, V−H 2 and V−2 H 2 .
Equilibrium constants K i (T ) were determined by the free energy of adsorption ∆G  which contains the same type of nodes studied here, results in denser packing of metal sites and an increased maximum allowed hydrogen uptake of 21.18 g/L.1 The secondary axis of binding curves in Fig. 3 shows predictions of total hydrogen uptake in g/L in V(II)MFU-4l.Since these predictions correspond to H 2 bound at the V sites and do not account for additional H 2 that will physisorb within the pore, they are estimates of the lower bounds for uptake.Advanced Grand Canonical Monte Carlo (GCMC) models for surface coverage such as the ones formulated by Snurr and coworkers 10,51 can be used predict surface coverage due to dispersion interactions of H 2 with the pore that also takes into account other structural parameters of the MOF such as surface area and free volume.
Strong uptake is seen at liquid nitrogen temperature at low pressures (0-1 bar) with almost all sites saturated by 0.5 bar for the F -and Cl -ligated node.The isotherm for the Br -containing nodes shows sigmoidal shape with comparatively slower uptake in the low pressure regime.The strongly negative ∆∆E for the Br -containing node is attenuated when one considers the effect of ZPVEs and vibrational contributions to the enthalpy of and free energy binding.Consequently, the binding isotherm at 77 K does not show a clear step on the ligation of the second H 2 (see fig.C) and 1 atm for H 2 binding to V(II)XMFU-4l, X = F -, Cl -, Br -, I -.Table .3. V(II)(F)MFU-4l emerges with the steepest uptake for the 100-5 bar pressure swing.
We estimate the usable capacity of the F -containing V(II)MFU-4l node at 4.1 g/L at 228 K.For the MFU-4 scaffold, which has a higher density of metal sites whilst retaining the same nodes, we predict a usable capacity of 11.9 g/L under these conditions.Interestingly, the predicted entropic penalty for binding the second hydrogen is greater when compared to that for the first for V(II)MFU-4l nodes containing heavier halides (Cl -, Br -, and I -) which offsets the stronger electronic binding energy for ligation of the second H 2 .
A key element in the materials design challenge for high volumetric density hydrogen storage in MOFs is the ability to reversibly bind multiple hydrogens at the first coordination sphere of a metal site.We predict that the V(II)MFU-4l node containing fluoride ligands can bind two hydrogen molecules per metal with binding enthalpies in the range of 15-25 kJ/mol.When compared to conventional open site metals in MOFs (e.g Cu(I)MFU-4l ), Figure 1: a) Cubic cell unit of Zn(II)MFU-4l, the parent material for the transmetallated V(II)MFU-4l studied here.b) Cluster model for V(II)MFU-4l with Cl -ions.Figures c) and d) are closeup views of the equilibrium geometries of singly and doubly loaded V(II)MFU-4l.Vanadium sites are rendered in grey, nitrogen in blue, chloride in green and Zn in violet grey.Distances are labeled in Å.

Figure 2 :
Figure 2: a) Comparison of Energy Decomposition of the first and second H 2 binding to V(II)(Cl)MFU-4l.Cooperative binding of H 2 to V(II) in MFU-4l is predicted by increase in binding energy on adsorption of the second hydrogen.b) Complimentary Occupied-Virtual pair (COVP) of charge transfer orbitals for backdonation from V(II) to the second ligated H 2 .

L 3 = 7 .
25 g/L.The smaller unit cell dimension (L = 21.63Å) in MFU-4, 50 3a).Weak binding in V(II)(I)MFU-4l results in a small percentage of sites being occupied.Open metal sites in MOFs that saturate after binding one H 2 would give overwhelmingly positive ∆∆H and ∆∆G.Only the fluoride containing node shows a small positive ∆∆H(T)= 4.8 kJ/mol (∆∆G = 4.6 kJ/mol).This node accomodates two hydrogens with strong binding enthalpies of ∆H 1 (T ) = -19.3and ∆H 2 (T ) = -14.5 kJ/mol at liquid nitrogen temperatures.Consequently, the binding curves show a very steep uptake in the low pressure (0-0.1 bar) range at 77 K. Cooling to cryogenic temperatures (-90 • C/183 K or lower) is expensive and for practical applications of a nanoporous material one would want steep uptake at cool temperatures that are relatively inexpensive to access.Hence we look at the binding curve at -45 • C (228 K) in Fig. 3b.The corresponding energy, enthalpy and free energy changes are reported in

Figure 3 :
Figure 3: Predicted binding curve for halide variants of V(II)MFU-4l at the temperature of liquid nitrogen (a) and at 228 K (-45 °C) (b).

two H 2
s binding to V(II) in V(II) exchanged MFU-4l doubles the volumetric capacity of the material at liquid nitrogen temperatures, and offers strong uptake in the 100-5 bar pressure swing at -45 • C.Consistent with the completion of octahedral coordination at the V(II) site, binding of the second hydrogen in the presence of the first is accompanied by an increase in charge transfer between the metal and hydrogen.Forward donation from H 2 to V(II) increases with the ligation of the second hydrogen to the node with back-donation from the metal site remaining largely unchanged.Nodes containing heavier halides, especially Br -, show stronger binding for ligation of the second hydrogen at the V(II), but this energy lowering is counteracted by a commensurate increase in the entropic penalty.Change in binding behaviour with different halide ions in the node suggests that the anion coordinated to the binding site in this material can be used to finely tune the binding strength with H 2 or other guests.The MFU-4l parent scaffold, which has been shown to allow a broad range of metal substitutions at its tetrahedral metal sites, is critical in accommodating two hydrogens at the V(II) site.Provided the synthetic challenge of realizing a stable material with suitable density of V(II) sites is met, this strategy can be used to enable high density H 2 storage at close to ambient temperatures of -45 • C. Thus, V(II)-MFU-4l and V(II)-MFU-4 represent promising candidate materials for the coordination of multiple H 2 molecules.
Figures 1c. and 1d.show close-ups of equilibrium geometries for single and doublehydrogen loaded V(II)MFU-4l nodes, with Cl − .The second hydrogen binds closer to the metal site with a center of mass distance of 2.09 Å (the first H 2 binds at 2.16 Å).The singly loaded node assumes a distorted trigonal bipyramidal configuration around the V(II) site (see Fig.1c) with the equatorial plane formed by the chloride ion and nitrogens I and II, and the axial sites occupied by hydrogen and nitrogen III.The axial H 2 and N III make angles of 90.6°a nd 91.31°with the equatorial N II .The doubly loaded node forms a distorted octahedron −4kJ/mol in energy, and 3 × 10 −5 a.u in the maximum gradient component.Counterpoise corrections proposed by Boys and Bernardi 39 were employed to reduce basis set superposition error (BSSE) in energy evaluations.aroundthe V(II) site with the bonding axes of two hydrogens aligned nearly perpendicular to one-another.Completion of the octahedron at the metal site brings both hydrogens inward in the coordination sphere which is indicative of stronger binding interactions.In going from zero to one, to two hydrogens bound, the V(II)-Cl distance changes from 2.30, to 2.33, to 2.39 Å.To decompose metal-H 2 binding into physically meaningful interactions, we use Energy Decomposition Analysis (EDA), a tool based on Absolutely Localized Molecular Orbitals (ALMOs).
Cl-H 2 V(II)L 3 ClH 2 -H 2 ∆∆ Table 2a list energy contributions ∆ for single (V(II)-H 2 ) and double (V(II)H 2 -H 2 ) loading, and Column III shows changes E(kJ/mol) V(II)L 3 Table 49thalpy of binding.We employ a stoichiometric two-site model suggested by Dill and Bromberg49to predict binding isotherms.Hydrogen binds the V(II) center in the sequence