First-Principles Calculation of 1H NMR Chemical Shifts of Complex Metal Polyhydrides: The Essential Inclusion of Relativity and Dynamics

1H NMR spectroscopy has become an important technique for the characterization of transition-metal hydride complexes, whose metal-bound hydrides are often difficult to locate by X-ray diffraction. In this regard, the accurate prediction of 1H NMR chemical shifts provides a useful, but challenging, strategy to help in the interpretation of the experimental spectra. In this work, we establish a density-functional-theory protocol that includes relativistic, solvent, and dynamic effects at a high level of theory, allowing us to report an accurate and reliable interpretation of 1H NMR hydride chemical shifts of iridium polyhydride complexes. In particular, we have studied in detail the hydride chemical shifts of the [Ir6(IMe)8(CO)2H14]2+ complex in order to validate previous assignments. The computed 1H NMR chemical shifts are strongly dependent on the relativistic treatment, the choice of the DFT exchange–correlation functional, and the conformational dynamics. By combining a fully relativistic four-component electronic-structure treatment with ab initio molecular dynamics, we were able to reliably model both the terminal and bridging hydride chemical shifts and to show that two NMR hydride signals were inversely assigned in the experiment.


INTRODUCTION
Transition-metal hydride complexes play an important role in catalytic transformations, including transfer hydrogenation. 1−4 The study of their molecular and electronic structure is crucial for a better understanding of the reaction mechanisms and for the design of more efficient catalysts. Because metallic hydrides are difficult to locate by X-ray diffraction, their prime characterization is by NMR spectroscopy. Transition-metal hydride complexes with partially filled d shells often occupy extreme positions in the 1 H NMR shift range, with lowfrequency shifts as low as about −50 ppm. 5−7 These unusual 1 H NMR shift values are primarily attributed to relativistic spin−orbit (SO) coupling, usually denominated by the SO-HALA (heavy atom on the light atom) effect, 8 which induces relativistic shielding at the hydrogen nuclei, leading to the characteristic negative 1 H NMR chemical shifts. 6,8−13 A detailed discussion of the mechanisms that dictate the size and sign of the SO-HALA effect has recently been provided by Vícha et al. 7 In brief, partially occupied heavy-atom valence shells induce relativistic shielding at the light-atom nuclei, while empty heavy-atom valence shells induce relativistic deshielding. In particular, the light-atom nuclei are relativistically shielded by 5d 2 −5d 8 and 6p 4 metals. Interestingly, the maximum shielding in 5d transition-metal complexes is observed in 5d 6 (Ir III ) complexes. 7 Various novel iridium deactivation products formed in the catalytic conversion of glycerol into lactic acid 14 have been studied and characterized by combining experimental and computational approaches. 15−18 These iridium species have unique bow-tie structures generated by a central tetra-or hexairidium core bound to multiple N-heterocyclic carbene (NHC) and hydride ligands. For instance, a surprisingly high hydride content was found in the iridium hexamer [Ir 6 (IMe) 8 (CO) 2 H 14 ] 2+ complex (1; IMe = 1,3-dimethylimidazol-2-ylidene; Figure 1a). 15 The hydride ligands were not located by single-crystal X-ray diffraction studies, and a computational approach was instead adopted to estimate the hydrogen positions, converging into the complex structure shown in Figure 1, with 10 bridging and 4 terminal hydride ligands. 15 A complete assignment of all 1 H and 13 C NMR resonances, including the 1 H NMR of hydride ligands, was achieved by combining the information gained from both density-functional-theory (DFT) geometry optimization and 1D/2D NMR experiments. Although this strategy was consistent with the experimental spectra and DFT data, the question then became whether the assignment of the hydride signals made as described above can be supported by more advanced methods. In this regard, the use of computational NMR spectroscopy based on relativistic electronic-structure theory provides an alternative route to reproducing the experimental spectra and support the characterization of these species. Therefore, we turned to the challenging task of computing the hydride 1 H NMR chemical shifts of complex 1 by accurately taking into consideration relativistic, solvent, and dynamic effects. Although several individual effects can influence the quality of the calculated NMR parameters, various studies have demonstrated that the inclusion of both relativistic and solvent effects is essential for proper prediction of the NMR chemical shifts in transition-metal complexes. 19−22 Calculation of NMR parameters in transition-metal complexes requires both an accurate representation of the system (e.g., geometry) and environment (solution) and a high-level electronic-structure method for the NMR calculation itself. 21,23,24 Recent theoretical and implementational advances have made it possible to carry out all-electron quantumchemical NMR calculations using a quasi-relativistic (twocomponent) or a fully relativistic (four-component) model with Hamiltonians including both scalar relativistic (SR) and SO interactions. 25,26 Likewise, ab initio molecular dynamics (AIMD) appears to be a particularly useful simulation technique to investigate the conformational flexibility of a system, including solvation effects and dynamical averaging of the calculated NMR chemical shifts. 27−30 Such first-principles calculations inherently account for anharmonicities, and the temperature can be set to study systems under more realistic conditions than static quantum-mechanical calculations at 0 K. 20

RESULTS AND DISCUSSION
2.1. Computational Relativistic Approach. The reported complex 1 is based on a polynuclear Ir 6 H 14 core bound to eight NHC ligands. Each iridium center, NHC ligand, and metal hydride has a symmetrically equivalent partner due to the inversion center at the core of the complex. Thus, the spectroscopic data of 1 showed 7 inequivalent hydrogen signals (in a narrow range between −16 and −20 ppm) that were assigned to 14 classical hydride ligands. 15 Figure 1b displays the iridium atomic numbering, hydride labeling (H a −H g ), and color code used in this study.
As a first approximation, the hydride 1 H NMR chemical shifts δ( 1 H) of complex 1 were calculated based on a static (fully optimized) structure in the gas phase, reported by Campos et al. 15 at the ωB97xd/LANL2TZ*(Ir),6-311G**(N, C, H) level. The ωB97xd functional 31 produced the best results when benchmarked against other functionals in a similar iridium complex. 15,32 For the static δ( 1 H) shift calculations, the performance of the PBE 33,34 and KT2 35 functionals was assessed and the role of relativity was analyzed by comparing the nonrelativistic (NR) approach, the scalar relativistic zeroth-order regular approximation (SR-ZORA), 36−38 the two-component SO relativistic zeroth-order regular approximation (2c-ZORA), 36−40 and the four-component (4c) fully relativistic approach based on the Dirac− Coulomb Hamiltonian. 41,42 The 2c relativistic corrections with the ZORA Hamiltonian were performed using the ADF program. 43,44 To calculate the 4c relativistic corrections, we used the ReSpect program 45 with a four-component Dirac− Kohn−Sham method (4c-DKS); see the Computational Methods section for more details. For a direct comparison with the experiment, all calculated 1 H shielding constants were converted to chemical shifts δ( 1 H) (in ppm) relative to the shielding of dichloromethane (CH 2 Cl 2 ), computed at the same level of theory.
The resulting hydride 1 H NMR chemical shifts, plotted as deviations from the experimental values (Δδ), are shown in Figure 2. The δ( 1 H) values of the hydrides calculated at the NR and SR-ZORA levels using the KT2 and PBE functionals deviate significantly from the experimental values with Δδ ≈ 3−8 ppm; see the blue and green columns in Figure 2. By contrast, inclusion of the SO contribution at the 2c-ZORA level improves the results significantly; the Δδ values at the 2c-ZORA level range between +1.4 and −1.1 ppm with the KT2 functional and between +1.8 and −0.3 ppm with the PBE functional; see Tables S1 and S2. The improved accuracy of the hydride shifts with inclusion of the SO interaction is to be expected, bearing in mind that the HALA effect is mostly due to SO coupling. 46−48 We conclude that the observed lowfrequency signals of complex 1 hydrides are mostly due to strong SO effects. Surprisingly, the 4c-DKS method produced slightly larger deviations than the 2c-ZORA method ( Figure  2). The Δδ values at the 4c-DKS level range between −1.3 and −2.9 ppm for KT2 and between −0.3 and −2.0 ppm for PBE. According to these results obtained using a static optimized structure, the improvement in the relativistic treatment from   Tables S1 and S2 for numerical data and standard deviations. The labeling of the hydrides is the same as that in Figure 1.
Inorganic Chemistry pubs.acs.org/IC Article 2c-ZORA to the fully relativistic 4c-DKS approach does not, in this particular case, lead to an improvement of the accuracy. However, this behavior may reflect an error cancellation at the 2c-ZORA level of theory, leading to a fortuitously good agreement with experiment. Certainly, this striking difference between the 2c-ZORA and 4c-DKS levels will be corrected after incorporating dynamical averaging in the 1 H NMR shift calculations, as shown in section 2.2.
For the purpose of validation, special attention was paid to the results obtained at the 2c-ZORA and 4c-DKS levels of relativistic treatment with the KT2 and PBE functionals; see the Computational Methods section. The overall performance of these approaches was evaluated on the basis of the total root-mean-square deviations (RMSDs) between the calculated and experimental values. The 2c-PBE and 2c-KT2 levels show the best performance, with nearly identical RMSD values of 0.8 and 0.9 ppm, respectively (Table 1). However, the deviations vary significantly among the hydrides: the 2c-PBE method yields large deviations (Δδ > 1.0 ppm) for H b and H d , while the 2c-KT2 method yields large deviations for H b , H f , and H g . At the four-component level, the 4c-PBE method can be considered to be an acceptable approximation with RMSD = 1.2 ppm, whereas the 4c-KT2 method shows the largest deviation from experimental values with RMSD = 2.4 ppm. In short, the δ( 1 H) values calculated at the 4c-DKS level depend strongly on the choice of the DFT functional.
In this study, we are particularly interested in a reasonable estimation of the 1 H NMR chemical shifts to reproduce and confirm the assignment of experimental spectra. We have therefore analyzed the experimental and calculated relative δ( 1 H) values by using H g as the reference; see Table 2. Note that the experimental Δδ values of H a −H f hydrides are in a very narrow range of only 3.1 ppm. This range is well reproduced by the calculations, with maximum Δδ values of ∼4.0 and ∼5.0 ppm at the 4c-DKS and 2c-ZORA levels, respectively. Likewise, the relative δ( 1 H) values calculated at the 4c-DKS level are close to the 2c-ZORA results, differing by 1.3 ppm or less (see Δδ in Table 2). Compared with the relative experimental values, the 2c-ZORA method showed large deviations for the terminal H a and H b and bridging H d hydrides, while the 4c-DKS method showed deviations for H b and H d (see ΔΔδ in Table 2). According to the reported 2D nuclear Overhauser effect spectroscopy (NOESY) spectrum, which indicates the interactions between the hydrides and NHC ligands, 15 the H a , H b , and H d hydrides exhibited the largest number of noncovalent interactions with the methyl groups of the ligand. Note that these interactions can also alter the NMR resonances of the hydrides and may thus be important for the calculated δ( 1 H) values.
The suitability of the selected methods was also analyzed by comparing the absolute calculated and experimental δ( 1 H) trends from the H a to H f hydrides (Figure 3). At the 2c-ZORA level, we find H b > H d with a large difference (Δδ > 1.0 ppm), while the 4c-DKS level leads to H b ∼ H d . Furthermore, these levels of theory do not fully reproduce the experimental trend; the observed H a > H b and H c > H d trends, for instance, are not reproduced by the calculations. Hence, although the calculated H e −H g shifts are in agreement with the experiment (decreasing in the order H e > H f > H g ), it is possible that the calculated 1 H NMR chemical shifts are not sufficiently well described using a static optimized structure. In addition to the level of theory and relativistic effects, an important source of error in the NMR  Our first approach to investigating the solvent effects was to compute the static 1 H NMR chemical shifts by including an implicit continuum model approach for describing the bulk solvation in CH 2 Cl 2 . Solvation corrections were calculated at the 2c-ZORA level using the implicit conductor-like screening solvent model (COSMO) 49,50 and at the 4c-DKS level using the polarizable continuum solver model (PCM). 51 However, inclusion of the COSMO and PCM models leads only to minor changes in the shifts, up to 0.08 and 0.26 ppm, respectively (see Δδ solv in Table S4). The solute−solvent interactions that affect the hydride δ( 1 H) chemical shifts of complex 1 therefore cannot be properly modeled with a continuum model, even though such models have been successfully used in several cases to describe the bulk solvent effects on NMR calculations. 52−54 Among the selected methods, the 2c-ZORA-COSMO method (estimated using the PBE exchange−correlation functional) shows the best performance with RMSD = 0.8 ppm (Table S3).
Additionally, to determine the importance of dynamics and the effect of solvation, we performed AIMD simulations using the CP2K program package 55 where complex 1 was surrounded by solvent molecules (CH 2 Cl 2 ) and followed over time (see the Computational Methods section for details). From this trajectory, a total of 40 snapshots were used to yield dynamic (averaged) 1 H NMR chemical shifts for the seven different hydrides (H a −H g ) of complex 1. The molecular dynamics of complex 1 reveals drastic conformational changes on the Ir 6 H 14 core region (Table S5 and Figures S1−S5). As a result, a large range of variation in the calculated 1 H NMR shift values is observed along the full trajectory (40 ps), with changes of up to 9.2 ppm (estimated at the 4c-DKS level using the KT2 functional) for each of the hydrides (Figures S6−S9).
The dynamic 1 H NMR hydride chemical shifts were calculated at the 2c-ZORA and 4c-DKS levels with the KT2 and PBE functionals. Interestingly, inclusion of molecular dynamics caused significant changes in the δ( 1 H) chemical shift values and trends. In particular, the dynamical averaging significantly increased the δ( 1 H) values of hydrides H b −H g and slightly decreased the δ( 1 H) value of H a (Figures S10 and S11). Analysis of the dynamic relative δ( 1 H) results by using H g as the reference and a comparison with the experimental values are shown in Table 3. Notably, some deviations between the calculated and experimental results still persist in the dynamic approach; namely, the 2c-ZORA method shows large deviations in the relative δ( 1 H) values for H b and H d hydrides, while the 4c-DKS method shows deviations for H a and H d (see ΔΔδ in Table 3).
Because the dynamic δ( 1 H) chemical shift calculations are not in full agreement with the experimental trend, we systematically examined all of the results obtained for each of the hydrides. Note first that the largest difference between the 2c-ZORA and 4c-DKS relativistic levels is mostly due to the changes in the two terminal H a and H b hydrides ( Figure  4a). Raising the level from the quasi-relativistic 2c-ZORA to the fully relativistic 4c-DKS approximation causes a large shielding on these two hydrides, decreasing the chemical shifts by approximately −2.0 ppm (using PBE). Certainly, the   (Figure 4a). Among the selected methods, the dynamic 4c-DKS results using the PBE functional must be considered to be the most reliable because they represent the highest level of theory employed and feature the smallest deviation (RMSD = 0.6 ppm). On the basis of the above analysis, our results indicate that two signals, the terminal H a and bridging H d hydrides, were assigned inversely in the experimental study. To support this proposal, we analyzed the dynamic δ( 1 H) chemical shifts by swapping the experimental H a by H d values (Figure 4b). The reassignment of H a to H d significantly improved the correlation with the experimental data. For instance, the dynamic δ( 1 H) chemical shift values calculated at the 4c-DKS level showed an excellent agreement between theory and experiment for all hydrides (RMSD = 0.2 ppm); see Figure S12. At the 2c-ZORA level, all of the bridging hydrides (H c −H f ) were obtained in good agreement with the experiment, but substantial deviations (>1.0 ppm) were observed for the terminal H a and H b hydrides. Although some minor differences between the 4c-DKS and 2c-ZORA approaches can be expected, for example, from the different types of basis sets (Gaussian vs Slater basis sets in ReSpect and ADF, respectively), our results demonstrate the need to include the relativistic effects at a four-component level.

Compatibility of Reassignment with the Experimental Data.
To check if the proposed reassignment is compatible with the original experimental spectra, we have investigated the reported information from both DFT calculations and 1D/2D NMR experiments. In the experimental study, NOE cross-peaks were traced for the seven hydride resonances. The interactions between the hydrides and N−Me groups gave 18 intense and 7 weaker NOE signals. 15 The DFT-optimized structure was used to rationalize the NOESY spectrum on the basis of the calculated H···H distances between the hydrides and N−Me groups, with only one of the possible assignments being consistent with all of the 2D NMR and DFT data.  Inorganic Chemistry pubs.acs.org/IC Article Important for the assignment, only three of the resonances due to the N−Me groups present a single strong NOE peak with a metal hydride. In agreement with the NOESY spectrum, the optimized structure showed only three methyl groups whose protons exhibited a single H···H distance below 3.0 Å, namely, C(1), C(2), and C(3). Nevertheless, it should be noted that the NHC ligands are very flexible and present conformational changes that are not considered in the static optimized structure. Thus, we have analyzed the H···H distances (Å) between the metal hydrides and the N−CH 3 wingtips. The interactions of the C(1), C(2), and C(3) methyl groups with H a , H b , and H d are shown in Figure 5. Both static and dynamic H···H bond distances show that C(1), C(2), and C(3) exhibit a strong interaction with hydrides H a , H b , and H d , respectively. Nevertheless, analysis of the dynamic H a −C(3) distances along the trajectory reveals also considerable strong interactions (below 3.0 Å) between C(3) and H a ( Figure 5). Thus, while the static optimized structure shows that the C(3) methyl group interacts with a single metal hydride (H d ), the dynamic structure indicates that C(3) can interact with two metal hydrides (H a and H d ).
2.4. Validation of the Optimized Protocol. Additionally, we studied in detail the case of the [Ir 4 (IMe) 8 H 10 ] 2+ complex 2, for which the exact positions of the hydrides were known from neutron diffraction experiments. 18 This species is based on a polynuclear Ir 4 H 10 core similar to complex 1, with four terminal hydrides, four equivalent hydrides bridging the shorter Ir···Ir edges, and two equivalent hydrides bridging the longer Ir···Ir edges (Figure 6a). In agreement with the neutron diffraction structure, the 1 H NMR spectrum showed three inequivalent H signals integrating to 10 classical hydride ligands in a 4:4:2 ratio. 18 The hydride chemical shifts of complex 2 were first examined with a static protocol, based on a fully optimized structure at the ωB97xd/LANL2TZ*(Ir),6-311G**(N, C, H) level and subsequent 1 H NMR calculations at the 2c-ZORA and 4c-DKS levels. The static approach performs well, reproducing the experimental trend H a > H b > H c ( Figure  S15). The calculated shifts yield small RMSD values that range between 0.61 ppm (2c-KT2) and 1.05 ppm (4c-PBE), although the deviations seen for each hydride differ appreciably among the selected methods (Tables S8 and S9). The dynamic protocol was then tested at the 4c-DKS level using the KT2 and PBE exchange−correlation functionals (Tables S10 and S11). As shown in Figure 6b, the hydride 1 H NMR chemical shifts δ( 1 H) calculated at the 4c-PBE method yielded the trends closest to experiment and showed the smallest deviation (RMSD = 0.41 ppm). Notably, the optimized dynamic protocol at the 4c-PBE approach was proven to be the most accurate in the prediction of the 1 H NMR hydride chemical shifts of complex 2, which is consistent with the results obtained for complex 1. Hence, this is a reliable protocol to calculate the NMR chemical shifts of the terminal and bridging hydrides and has potential utility for the assignment of hydride signals in other challenging metal polyhydrides. 57−59

CONCLUSIONS
In this study, we calculated the 1 H NMR hydride chemical shifts of the iridium polyhydride complex 1 by accurately taking into consideration relativistic, solvent, and dynamic effects. The reliability of the computed 1 H NMR hydride values was assessed by comparing them with the experimental signals reported by Campos et al., providing a new strategy to predict and support the characterization of these species. The calculated 1 H NMR hydride chemical shifts are strongly dependent on the relativistic treatment; the SO contribution is the most important for the accurate reproduction of the NMR shift ranges, in particular, for the terminal H a and H b hydrides. The static 1 H NMR chemical shifts calculated using both the 2c-ZORA and 4c-DKS approaches do not fully follow the experimental trend, showing large deviations for the H a , H b , and H d hydrides. Furthermore, the 4c-DKS level was shown to be highly dependent on the choice of the DFT exchange− correlation functional, where the PBE functional performs better than KT2.
The role of an implicit solvent model was found to be negligible for the hydride chemical shifts, with minor changes of up to 0.2 ppm (estimated at the 4c-DKS level using the PBE functional). By contrast, the effect of dynamical averaging using AIMD simulations resulted in significant changes in the 1 H NMR chemical shift values and trends. The dynamic 1 H NMR chemical shifts at the 4c-DKS level showed large deviations for the H a and H d hydrides, which were inversely assigned in the experiment. The reassignment of H a and H d gave significantly improved correlations with the experimental data. The 4c-DKS level showed an excellent agreement between theory and experiment for all hydrides (RMSD = 0.2 ppm). In contrast, the 2c-ZORA level gave good results for the bridging H c −H g hydrides but failed for the terminal H a and H b hydrides, showing the importance of including a fourcomponent relativistic methodology.
Among the selected methods, the 4c-DKS level using the PBE functional and including dynamical averaging when calculating the 1 H chemical shifts showed the best perform-   (2) complexes were calculated using the optimized structure reported at the ωB97xd/LANL2TZ*(Ir),6-311G**(N, C, H) level, which proved best in simulating the X-ray structure. 15 Additionally, dynamic δ( 1 H) chemical shifts were obtained from quantum-chemical calculations based on an ensemble of structures from AIMD simulations. AIMD simulations of complex 1 were run in an explicit dichloromethane (CH 2 Cl 2 ) solvent according to the Born−Oppenheimer approximation using the CP2K program package. 55 The initial model system, created using the PACKMOL package, 60 consists of complex 1 (optimized at the ωB97xd/LANL2TZ*(Ir),6-311G**(N, C, H) level) surrounded by 233 CH 2 Cl 2 molecules in a cubic box of 30.0 Å 3 edge to reproduce the appropriate density of 1.325 g/mL. The simulation cell was treated under periodic boundary conditions for all AIMD calculations performed using the Kohn−Sham DFT with the PBE exchange−correlation functional, 33,34 in a mixed DZVP Gaussian 61 and auxiliary plane-wave (200 Ry cutoff) basis set. Core electrons were described using pseudopotentials of the Goedecker− Teter−Hutter type. 62 Dispersion forces were taken into account using Grimme's D3 model. 63 AIMD simulations of complex 2 were performed in similar conditions, using a smaller cubic box of 25.0 Å 3 edge containing 136 CH 2 Cl 2 solvent molecules. The initial configuration was relaxed by AIMD simulation using a microcanonical (NVE) ensemble, until an average temperature of 298 K was reached. After equilibration, the simulation was then run by using a canonical (NVT) ensemble with a temperature of 298 K maintained with the CSVR algorithm. 64 The trajectory was extended up to 40 ps (complex 1) and 30 ps (complex 2), with a time step of 0.25 fs. From the long 40 ps simulation, a total of 40 snapshots taken at regular 1 ps intervals were used for obtaining the dynamic average of the 1 H NMR shielding constants. Because we are interested in the chemical shifts of the thermodynamic ensemble of structures, the molecular coordinates were taken as provided by AIMD simulation without further geometry optimizations.
4.1. Relativistic 1 H NMR Chemical Shift Calculations. The calculations of the hydride 1 H NMR chemical shifts were performed using the 4c-DKS method in combination with the Dirac−Coulomb Hamiltonian, 41,42 as implemented in the ReSpect program. 45 The PBE 33,34 and KT2 35 functionals were tested because both functionals have been proven to give reliable data for NMR chemical shifts; the KT2 functional is specifically optimized to provide high-quality shielding constants for light main-group nuclei. Hybrid exchange− correlation functionals require computationally more demanding calculations and were not considered in this study. Mixed uncontracted basis sets were employed to preserve high accuracy at lower computational cost, combining Dyall's VDZ 65−68 basis set for iridium and the IGLO-II 69 basis set for the rest of the light atoms; Dyall's basis sets are designed for relativistic calculations and perform well for calculation of the NMR parameters, 70 while IGLO basis sets are designed for computation of the magnetic properties. For acceleration of the calculations, the relativistic electron repulsion integrals and related two-electron Fock contributions were calculated using the resolution-of-identity for the two-electron Coulomb term, 71 which has demonstrated good accuracy for large systems while reducing the computational cost. All of the hydride δ( 1 H) chemical shifts were performed in the gas phase, and the solvent effects were assessed in the static δ( 1 H) calculations using the PCM 51 for CH 2 Cl 2 .
The four-component results were compared with those obtained by a quasi-relativistic 2c-ZORA(SO) method, 36−40 as implemented in the Amsterdam Density Functional (ADF) program, 44 using either the PBE or KT2 functional in conjunction with all-electron Slater-type orbital basis sets designed for relativistic ZORA calculations, combining the triple-ζ quality plus two sets of polarization functions (TZ2P) 73 basis set for iridium and the even-tempered polarized valence quadruple-ζ (ET-pVQZ) 72 basis set for the rest of the light atoms. Solvation effects were taken into account through the implicit COSMO 49,50 model for simulating bulk solvation in CH 2 Cl 2 . The gauge-origin dependence was handled using gauge-including atomic orbitals (GIAO) approach. 42,74 A dataset collection of the computational results is available in the ioChem-BD repository 75