Dynamic effects on ligand field from rapid hydride motion in an iron(ii) dimer with an S = 3 ground state

Hydride complexes are important in catalysis and in iron–sulfur enzymes like nitrogenase, but the impact of hydride mobility on local iron spin states has been underexplored. We describe studies of a dimeric diiron(ii) hydride complex using X-ray and neutron crystallography, Mössbauer spectroscopy, magnetism, DFT, and ab initio calculations, which give insight into the dynamics and the electronic structure brought about by the hydrides. The two iron sites in the dimer have differing square-planar (intermediate-spin) and tetrahedral (high-spin) iron geometries, which are distinguished only by the hydride positions. These are strongly coupled to give an Stotal = 3 ground state with substantial magnetic anisotropy, and the merits of both localized and delocalized spin models are discussed. The dynamic nature of the sites is dependent on crystal packing, as shown by changes during a phase transformation that occurs near 160 K. The change in dynamics of the hydride motion leads to insight into its influence on the electronic structure. The accumulated data indicate that the two sites can trade geometries by rotating the hydrides, at a rate that is rapid above the phase transition temperature but slow below it. This small movement of the hydrides causes large changes in the ligand field because they are strong-field ligands. This suggests that hydrides could be useful in catalysis not only due to their reactivity, but also due to their ability to rapidly modulate the local electronic structure and spin states at metal sites.

.1005(4) 11.9085(4) 11.8026 (3) 11.8009 (5) b (Å) 11.3249(4) 11.0443(4) 10.9628 (3) 10.9639 (4) c (Å) 18.9146 (4) 17.4990 (7) 17.4178 (5) 17.4087 (7) a (º) 103 Figure S2. Rotation photographs as a function of temperature during the phase change of 1 to 1-LT. Note that the split diffraction spots return to the original unsplit spots upon warming.   Attempted Boltzmann fit of temperature-dependent Mössbauer data. The experimental zero field Mössbauer spectra of toluene-free crystals of the hydride dimer shown in Figure 2 are simulated in Figure S6 with a relaxation model for the electric quadrupole interaction in order to explain the collapse of quadrupole splitting from two distinct iron site with increasing temperature. For this model we assume correlated 'flipping' of the electronic environment of both 57 Fe-nuclei between two situations, 'T' and 'P', with two different electric field gradient (efg) tensors; the Mössbauer spectra for 'P' and 'T' at the two iron sites are indistinguishable.
For the sake of simplicity, co-axial efgs of the same sign are assumed for 'T' and 'P', and crosswise correlated 'hopping' between 'T' and 'P' for both iron sites, with a common rate that increases with temperature. In the spectra, the relaxation occurs then between pairs of (adjacent) At 80 K, and practically also at 140 K, the system is in the static limit where two quadrupole doublets of equal absorption depth and Lorentzian line width are observed. At 140 K the isomer shifts and quadrupole splittings for both iron sites are slightly reduced as compared to 80 K (see Table S2 below). This is due to the second-order Doppler shift, and the presence of close lying excited spin-orbit states for the 3d 6 configuration of iron(II).
In the range 140 -170 K the coalescence processes of the two quadrupole doublets was simulated by the onset of quadrupole relaxation with temperature-dependent relaxation rates 1/τ, given in units of the Mössbauer line width of 0.2 mm/s (close to the natural line width).
The same set of the other Mössbauer parameters is used throughout the simulations.   Table S2. Note that the model was discarded because the rates (1/τ) yield an unrealistic Arrhenius plot ( Figure S7) due to an unreasonably steep temperature dependence. The dotted blue line shows the simulation for 140 K superimposed on the 80 K spectrum, to enable comparison.
However, when the rates are placed on an Arrhenius plot ( Figure S7), the best-fit line implies Ea = 10.4 kcal/mol and A = 2 x 10 19 s -1 . The latter value is too high to be physically meaningful, and we conclude that the interconversion involves a cooperative process instead (see below).

Arrhenius plot
Cooperative model for temperature-dependent Mössbauer data. The phase change observed for 1 at ca. 160 K is accompanied by a marked decoalescence of the Mössbauer spectrum from a single quadrupole doublet into two distinct subspectra ( Figure 2 of main text). Table S3 summarizes the result of the fits for the spectral series with two Lorentzian quadrupole doublets of equal intensity; the temperature dependence of the quadrupole splitting for both subspectra is depicted in Figure S8 by symbols 'D' and '+'.  Coalescence of Mössbauer spectra arises from dynamic processes and has been observed in the disorder of O2 bound to iron porphyrin complexes, 6 or the Jahn-Teller distortion of the cubic iron(II) site in cesium iron phosphate hexahydrate. In the latter case, time averaging of the distortion over the three canonical axes gives pseudo-cubic symmetry and collapses the otherwise large quadrupole splitting of the high-spin iron(II) ions to a single doublet above 100 K without significant line broadening. 7 In the following model, the coalescence of the Mössbauer spectra of 1 arises from a phase transition from a static low-temperature to a dynamic high-temperature phase, the latter having Moreover, the swaps of the dimer configurations cannot be stationary but must 'hop around' the ensemble, also on a fast time scale, since no low-temperature fraction of subspectra coexists with the coalescence spectra. Due to the time and ensemble average, the weight factors are given by the temperature-dependent fraction x(T) of molecules in the high-temperature phase Here, the factor 0.5 takes into account the 50:50 average of both electronic environments for both iron sites in the high temperature limit.
Since the coalescence of the quadrupole spectra occurs within a temperature range of less than 40 K, a Boltzmann-driven process of the Arrhenius-type is excluded (see above for our unsuccessful attempt to fit this model). Thus this second model involves cooperativity and an entropy-driven phase transition occurring within domains. As shown in Figure S8, a good fit of the quadrupole splittings was possible with the Sorai-Seki domain model 8,9 describing the fraction of molecules x(T) in the high-temperature phase where Tc is the transition temperature, R is the gas constant and nΔH is an enthalpy factor that depends on the size n of the domains, and which determines the width of the transition. Leastsquares optimization yielded the transition temperature Tc = 150 K, and an enthalpy factor, nΔH = 1525 cm -1 .
In accordance with the magnetic Mössbauer spectra described in the main text, the quadrupole splitting due to EFG(I) was taken to be negative, ΔEQ1 = -0.74 mm/s whereas ΔEQ2 was taken to be positive. Moreover, the asymmetry parameters η of the electric field gradients were taken to be zero. With this choice, the average quadrupole splitting derived from the 80 K values, (ΔEQ1 + ΔEQ2)/2 = 0.74 mm/s (if ΔEQ1 is negative) is far from the coalescence value ΔEQ1,2 = 1.25 mm/s at high temperature. This suggests that we cannot simply average the quadrupole splittings, because the EFG(I) and EFG(I) tensors are not collinear. This agrees with the fits to the magnetic Mössbauer spectra, the SQUID data, and the computational results.
Note that quadrupole splittings in general are not additive, but electric field gradients are.
This holds also here for the superposition of field gradients in a hopping or other coalescence process. Before addition of EFGs, the components Vij of both tensors have to be expressed in a common reference frame, and the result has to be diagonalized to determine the new quadrupole splitting. Since field gradients are traceless, Vxx + Vyy + Vzz = 0, the Vxx and Vyy components have signs opposite to the main component Vzz, which determines the quadrupole splitting.
Therefore an Euler rotation of the electric field gradient tensor EFG(II) with respect to the principal axis system of EFG(I) was essential for correct superposition of EFGs at high and intermediate temperatures and the proper simulations shown in Figure S8 (green and blue lines).
The best fit for a rotation around the y-axis yielded an Euler angle β = 110°. We note that this solution may be not unique since consideration of the asymmetry parameters η would lead to a variety of possible solutions that cannot be disentangled from the limited number of data available. For the sake of simplicity we kept η at zero for both sites and limited the Euler rotations to the angle β. Interestingly, the angle β = 110° is close to the value found independently from the magnetic Mössbauer spectra (β = 95°, see below).
In summary, the (de)coalescence of the quadrupole spectra of 1 at 150 K can be well  The data were integrated with the APEX2 software package and absorption corrections were applied with SADABS. 10 The structure was solved with SHELXT and was refined against F 2 on all data by full-matrix least squares with SHELXL. 11 All non-hydrogen atoms were refined anisotropically. Most hydrogen atoms were included in the model at calculated positions and refined using a riding model. The isotropic displacement parameters of all hydrogen atoms were fixed to 1.2 times the U value of the atoms to which they are linked (1.5 times for methyl groups). The only exceptions are H1 and H2, which were found in the difference map and freely refined. One reflection was recorded improperly due to instrument artifacts and subsequently omitted. The absolute structure parameter was 0.028 (6) Figure S13.  Figure S13. Plot of the isomer shifts of the two peaks as a function of temperature, from Mössbauer spectra collected with closely spaced temperatures. This shows that the coalescence temperature of 163 ± 3 K is not significantly different between H and D isotopologues.
where the first term is the Heisenberg-Dirac-van Vleck Hamiltonian with the exchange coupling constant J, the g symbols denote the local electronic g matrices, and D and E/D are the single-ion axial and rhombic zero-field splitting (zfs) parameters for the pseudo tetrahedral (Td) and square planar (SP) sites. The hyperfine interaction for 57 Fe on Td and SP sites was calculated with the usual nuclear Hamiltonian (7) where A is the hyperfine coupling matrix and S denotes the local spin of Td and SP sites. HQ is the usual Hamiltonian for electric quadrupole interaction. Isomer shifts were taken to be additive.
have to be considered as described in the following.

Nuclear Hamiltonian and Hyperfine Coupling
For the evaluation of magnetically split Mössbauer spectra within the spin-Hamiltonian formalism the purely S-dependent Hamiltonian must be extended by an appropriate nuclear Hamiltonian for the nuclear spin I: The leading term in is usually the magnetic hyperfine coupling which connects the electron spin S and the nuclear spin I. It is parameterized by the hyperfine coupling tensor . The I-dependent nuclear Zeeman interaction and the electric quadrupole interaction are included as 2 nd and 3 rd terms. Their detailed description for 57 Fe is provided in section 4.3 and 4.4. The total spin Hamiltonian for electronic and nuclear spin variables is then: . (

4.75)
To calculate Mössbauer spectra, which consist of a finite number of discrete lines, the nuclear Hamiltonian, and thus also , has to be set up and solved independently for the nuclear ground and excited states. The electric monopole interaction, i.e. the isomer shift, can be omitted here since it is additive and independent of M I . It can subsequently be added as an increment δ to the transition energies of each of the obtained Mössbauer lines.

Separation of I-and S-Dependent Contributions
The obvious basis for the diagonalization of are the product functions |S,M S > ⊗ |I,M I >, which yield (2S+1)(2I+1) eigenstates and eigenvalues. The interpretation Details on Magnetic Susceptibility. We were able to fit the data with different spin Hamiltonian (SH) models (eq. S6 above), one with two S = 2 subsites (S2S2 model) and one with an S = 2 subsite and an S = 1 subsite (S2S1 model, shown in Figure 7 of main text).
The S2S2 model with two equivalent S = 2 subsites required g = 2.0 -2.2, large negative zero-field splitting (zfs) with D of -52 cm -1 , and exchange coupling of J = +3.9 cm -1 ( Figure   S14). Importantly, this fit required us to postulate the presence of a diamagnetic impurity comprising 15% of the sample (to match in particular the high-temperature level cMT(T)    Figure 7, we also tested a S2S1 model with one S = 2 subsite and one S = 1 subsite. The S2S1 model did not require impurities to fit the data, though it requires anisotropic g-values that are significantly shifted from 2 at both the pseudo-tetrahedral (FeTd) and square-planar (FeSP) iron sites, with distinctly different zfs parameters. We were able to model the cMT(T) data as well as the multi-field M(µBB/kT) data well with the following SH parameters: STd = 2, gTd = 2.47, 2.47, 2.61, and DTd = -26 cm -1 , and SSP = 1, gSP = 3.65, 3.65, 1.0, and DSP = +40 cm -1 . The exchange coupling constant was J = +63 cm -1 (using the convention -2J STd×SSP). We restrained the local D and g matrices to be axial to avoid overparametrization of the model. However, for fitting the DSP tensor and gSP matrix of the FeSP site had to be rotated relative to the magnetic axes of the FeTd site by an Euler angle β = 90(±10)° around the y axis. The rotation makes the local magnetic quantization axes, zSP and zTd, of SP and Td sites perpendicular to each other, and xSP to be along the zTd direction (Figure 8, top).
These SH axes could be associated with the molecular structure of 1 as indicated in Figure 8  Interestingly, in spite of the low projection factor (2/30) for the SP site, the rotated DSP tensor was also important for a satisfying global fit of cMT(T) as well as the multi-field M(µBB/kT) data. In the sense of a perturbation contribution, the negative xx and yy components of DSP being aligned by the Euler rotation along the negative main component of DTd (shown as blue double-arrows in Figure 8) helped to sufficiently isolate the low-lying pair of highest ms levels and therefore to foil nesting of the simulated magnetization curves at different fields.
The large zfs of the Stotal = 3 ground state, being dominated by axial DTd = -26 cm -1 , causes the Zeeman splitting of the low-lying "ms,total = ±3" levels to be strong for fields in the zdirection but weak in x/y directions. This feature renders an "easy axis of magnetization." For 1, the magnetic anisotropy is increased by g anisotropy of the ground state septet ( Figure S16   hyperfine coupling values are inferred from the field dependence of the spectra ( Figure S16), which exhibit a significant increase of the overall magnetic splitting for both subspectra due to a positive sign of the internal field Bint = -<S>A /gnβn, where <S> is the spin expectation value (negative for the lowest ms level in the field). The trend of the spectra is reproduced by the corresponding spin Hamiltonian parameters, but we cannot reproduce the narrow experimental linewidths of these spectra. This may arise from torqueing of sample crystallites by the stronger fields. Therefore, we refrained from further refinement of the spin Hamiltonian description of the hydride dimer from the Mössbauer field dependence.
Considering the results of CASSCF studies (see below), we also modeled an isolated Stotal = 3 system. This model can also provide a good global fit of magnetic susceptibility and Mössbauer data. In this model the nesting of the M(µBB/kT) curves did not fit quite as well ( Figure S17); however note that E/D was fixed to 0 in the fit, whereas the CASSCF calculations predicted a large E/D of 0.29, which might explain the deviation. In this simulation, axial symmetry of the g matrix was assumed, because the magnetic data are not sensitive to gt,y, and its value cannot be accurately determined. As shown in Figure S18, the low-lying doublet hardly splits when the applied field aligns along y.  Magnetic Mössbauer spectra were recorded on a conventional spectrometer with alternating constant acceleration of the g-source (MPI-CEC, own construction). The sample temperature was maintained constant in a cryogen-free closed-cycle Mössbauer cryostat from Cryogenic Ltd, equipped with top-loading variable-temperature insert (VTI) and a split-pair superconducting magnet for fields up to 7 T, which were oriented perpendicular to the γ rays here.
The 57 Co/Rh source (1.8 GBq) was kept at room temperature and was positioned inside the gap of the magnet system by using a re-entrant bore tube. The source was adjusted horizontally to a zero-field position in ca. 9 cm distance from the sample. The sample (KCM-II-121) was kept The Mössbauer spectroscopic parameters were computed using the same density functional as for the geometry optimization step. We used the CP(PPP) 28 basis set for Fe, the TZVP 29 basis set for N, and the bridging H atoms and the SV(P) basis set 30 for remaining atoms.
Isomer shifts d were calculated from the electron densities r0 at the Fe nuclei by employing the linear regression: (8) Here, C is a prefixed value, and a and b are the fit parameters. Their values for different combinations of the density functionals and basis sets can be found in our earlier work. 31 Quadrupole splittings DEQ were obtained from electric field gradients Vij (i = x, y, z; Vii are the eigenvalues of the electric field gradient tensor) by using a nuclear quadrupole moment Q( 57 Fe) = 0.16 barn: 32 (9) Here, is the asymmetry parameter.
For complete active space self-consistent field (CASSCF) computations, 33 the active space was chosen to distributed 14 electrons into 11 orbitals including all Fe 3d orbitals of the Td and SP sites, as well as the -bonding counterpart of the FeSp dyz orbital. Quasi-restricted orbitals derived from DFT calculations were used as the initial orbitals. To capture dynamic correction, on top of CASSCF wavefunctions, the second-order N-electron valence perturbation theory (NEVPT2) 34 was employed.
The geometry-optimized computational models at the BP86/def2-TZVP, B3LYP/def2-TZVP and TPSSh/def2-TZVP levels of theory unanimously demonstrated the lowest-energy nonet state of 1 to be higher in energy than the corresponding septet (by 29.3, 10.5 and 16.8 kcal/mol, respectively). Note that these differences exceed the ±7 kcal/mol uncertainty level that is estimated for spin-state energetics by DFT; 35 hence, these strongly support the S2S1 model of complex 1 with an S = 3 ground state. For both the S2S1 and S2S2 models, the dihedral angles between the two Fe-nacnac planes are roughly 90°, likely arising from steric repulsion between the bulky aryl groups. As elaborated below, the potential energy surfaces with respect to the dihedral angle ( ) between the equatorial plane of the FeSP site and the Fe2H2 core are rather flat, and the location of the energy minimum in the calculations is not robust. Because of this ambiguity, we also calculated the Mössbauer parameters for these geometries. As summarized in Table S4, the S2S1 model at the experimental geometry (from neutron crystallography) predicted Mössbauer parameters most accurately. Therefore, the geometry obtained from neutron diffraction was used for the more detailed information on the electronic structure.