Density Functional Study of 17O NMR Chemical Shift and Nuclear

5200

J. Phys. Chem. B 2000, 104, 5200-5208

Density Functional Study of 17O NMR Chemical Shift and Nuclear Quadrupole Coupling Tensors in Oxyheme Model Complexes Martin Kaupp,*,† Carme Rovira,‡ and Michele Parrinello† Max-Planck-Institut fu¨ r Festko¨ rperforschung, Heisenbergstr. 1, D-70569 Stuttgart, Germany, and UniVersitat de Barcelona, Departament de Quı´mica Fı´sica, Facultat de Quı´mica, Martı´ i Franque` s 1, ES-08028 Barcelona, Spain ReceiVed: December 1, 1999; In Final Form: February 24, 2000

The 17O NMR chemical shift tensors for a dioxygen ligand bound to iron porphyrin model complexes, with or without inclusion of picket-fence-type substituents, have been calculated by density functional IGLO approaches, based on structures optimized in Car-Parrinello molecular dynamics simulations. The calculations confirm the experimentally found, extremely large 17O shifts, and the assignment of the higher-frequency signal to the terminal oxygen position. Metal-ligand π-back-bonding and the presence of low-lying excited states influence the chemical shift tensors characteristically. The magnitude and orientation of the shift tensors also reflect the close analogy to the bonding in the ozone molecule. Possible explanations for the temperature dependence of the solid-state NMR spectra are discussed. The computed 17O nuclear quadrupole coupling constants of ca. 11 MHz and ca. 17 MHz for bridging and terminal oxygen positions are close to the results for ozone, where theory and experiment agree well. This contradicts earlier assumptions of very small 17O field gradients for picket-fence oxyheme models and raises some questions about the analysis of the solidstate 17O NMR spectra. Possible reasons for the temperature dependence of the line widths in the solution spectra are suggested. We also discuss the implications of the 17O NMR spectra for the nature of the ground state of oxyheme complexes.

1. Introduction The binding of dioxygen to hemoproteins and synthetic model compounds remains a topic of great ongoing interest.1-3 NMR spectroscopy is a valuable tool for studying molecular and electronic structure, and it has also been used widely on iron porphyrin complexes.4 17O NMR spectroscopy of the bent endon coordinated dioxygen ligand itself in heme complexes has only become feasible during the past decade.5 Gerothanassis et al. studied the solution 17O NMR spectra of a number of synthetic oxyheme model systems.6,7 Oldfield et al.8 reported solid-state NMR results for the picket-fence porphyrin system, as well as directly for the myoglobin and hemoglobin proteins. 17O and 15N NMR spectra for Co-NO and Fe-RNO analogue systems have also been reported.9 The oxyheme complexes have been found to exhibit 17O resonances at particularly high frequencies (ca. 1200-2500 ppm relative to liquid water),6-8 consistent with the presence of lowlying excited states. The oxygen chemical shifts are closely related to the intriguing electronic structure of these important end-on dioxygen complexes. In the continuing discussion on the nature of the ground state of oxyheme systems, the 17O NMR spectra have been considered6-8 to support Pauling‘s closedshell singlet formulation,10 rather than the open-shell singlet model suggested by Weiss.11 In the solid-state 17O spectra of the picket-fence porphyrin system, both the isotropic chemical shifts and the shift tensors * Corresponding author. Present address: Institut fu¨r Anorganische Chemie, Universita¨t Wu¨rzburg, Am Hubland, D-97074 Wu¨rzburg, Germany. E-mail: [email protected]. † Max-Planck-Institut fu ¨ r Festko¨rperforschung. ‡ Universitat de Barcelona, Departament de Quı´mica Fı´sica.

show considerable dependence on temperature, in particular for the bridging oxygen positions.8 Furthermore, significant differences were found between terminal oxygen shifts in the solidstate and solution spectra.7,8 Unusually small 17O nuclear quadrupole coupling constants (NQCC) for both atoms of the dioxygen ligand had to be assumed to fit the static solid-state spectra.8 On the other hand, the temperature dependence of the line widths in the solution spectra was explained by a dominant relaxation via the electric field gradients, i.e., the NQCCs were found to be important.7 These are just a few of the apparent contradictions and open questions that the experimental NMR results pose. Quantum-chemical calculations should be able to provide further insight into the relations between the NMR parameters, molecular and electronic structure, and molecular dynamics. It has been demonstrated during the past few years that state-ofthe-art density functional theory (DFT) methods do now allow quantitative calculations of NMR chemical shifts even for transition metal complexes.12-14 Indeed, Oldfield and co-workers have recently applied DFT approaches to calculate 13C and 17O chemical shifts, as well as 57Fe NMR and Mo¨ssbauer parameters, of iron porphyrin carbonyl complexes.15,16 Calculations of 15N and 17O shifts for analogue complexes with NO and RNO ligands have also been performed, again accompanied by calculations of Mo¨ssbauer parameters (the latter were also computed for O2 complexes).9 However, to our knowledge studies of the particularly sensitive 17O NMR chemical shift or electric field gradient tensors for the O2 complexes by quantitative computational methods have not yet been reported. Recent large-scale ab initio molecular dynamics simulations of the Car-Parrinello type (CPMD) on various iron porphyrin complexes17-19 have provided a large amount of structural and

10.1021/jp994234k CCC: $19.00 © 2000 American Chemical Society Published on Web 05/10/2000

17O

NMR Studies of Oxyheme Model Complexes

J. Phys. Chem. B, Vol. 104, No. 21, 2000 5201

Figure 2. DFT-optimized structures used to calculate the NMR parameters for (a) O3 and (b) O3(H2O). Figure 1. View of the DFT-optimized structure of Fe(TpivPP)(2-MeIm)(O2) (cf. ref 18). The most important optimized structure parameters around the metal are: d(Fe-O) ) 1.778 Å, d(O-O) ) 1.295 Å, d(Fe-N, porphyrin) ) 2.008 Å (average), d(Fe-N, imidazole) ) 2.125 Å, ∠(Fe-O-O) ) 122.5°.

dynamical information. In this work, we have employed these data as a basis for quantum chemical studies of the 17O NMR chemical shift and nuclear quadrupole coupling tensors for both bridging and terminal oxygen positions of the end-on dioxygen ligand. We study the influence of the pivaloyl-amido substituents (Tpiv) of picket-fence iron porphyrin on the shifts and field gradients. Using results of the CPMD simulations,19 we examine the origin of the temperature dependence of the shifts in the solid state, and of the line widths in the solution spectra. The implications of the NMR data for the nature of the electronic ground state of oxyheme complexes will be discussed. 2. Computational Methods Our chemical shift calculations are based on the recently reported CPMD optimized20 structures18,19 for the picket-fence porphyrin complex Fe(TpivPP)(2-MeIm)(O2) (Figure 1), and for a simpler model system with an unsubstituted porphyrin ligand and an unsubstituted axial imidazole, FeP(Im)(O2).18 The actual models for which we studied the 17O NMR chemical shifts were (i) the optimized structure of the FeP(Im)(O2) model; (ii) an FeP(2-MeIm)(O2) model cut out of the optimized structure of Fe(TpivPP)(2-MeIm)(O2) (cf. Figure 1), with the TpivP substituents replaced by hydrogen atoms (in contrast to model i, this structure accounts for the ruffling of the porphyrin ring and the tilting of the imidazole ligand, due to the steric effects of TpivPP and 2-MeIm substituents); (iii) model ii augmented by partial point charges to simulate the electrostatic effects of the pivaloylamido substituents (charges were taken from the Mulliken population analyses of ref 18); and (iv) a more realistic model with the TpivP substituents replaced explicitly by smaller -CHd CH-NH-C(O)H amido groups (Tam), again based on the optimized structure of the larger picket-fence system (the methyl substituent was removed from the imidazole, i.e., the model consists of Fe(TamP)(Im)(O2)). Calculations at partly optimized structures, with the dioxygen ligand forced to eclipse one of the Fe-N vectors within the porphyrin plane, were carried out for models i and iv to estimate the influence of the O2 rotation above the plane on the chemical shifts. For model i, we also varied the Fe-O-O angle stepwise, starting from the fully optimized structure, and tilted the Fe-O bond vector from the normal of the porphyrin plane to study the dependence of the

shifts on these structural parameters. For details of the optimizations and CPMD simulations, the reader is referred to refs 1719. The ozone molecule and its hydrogen-bonded complex with one water molecule, O3(H2O), were also included for comparison. The structures of these smaller systems were also optimized at the DFT level, using the deMon program,21 the BP86 functional,22 effective-core potentials (ECPs) and DZP valence basis sets for oxygen,23,24 as well as a DZVP basis for hydrogen.25 The optimized structures, which may be considered to be of comparable quality as the CPMD results for the larger systems, are shown in Figure 2.26 The 17O nuclear shieldings for the dioxygen ligand (and for ozone) were calculated using the SOS-DFPT(IGLO) approach with the Loc1 approximation of the “Malkin correction”,27,28 as implemented in the deMon-NMR program,21,28 and for comparison also without the correction, i.e., at the uncoupled DFT (UDFT) level.28,29 The correction term has previously been found to improve the agreement with experimental shifts for main group compounds with low-lying excited states,27,28,30 and it was also found to be important for the ozone molecule,27,28,31 the electronic structure of which is closely related to the present systems of interest. We used the PW91 exchange-correlation functional,32 as well as FINE angular integration grids with 32 radial shells.28,33 For the two oxygen atoms of the dioxygen ligand, and for all atoms of O3 and O3(H2O), we employed the all-electron basis sets BII of Kutzelnigg et al.34 (often referred to as IGLO-II). The iron atom was treated with a scalar relativistic small-core ECP and an (8s7p6d)/[6s5p3d] valence basis set.35 The other nonhydrogen atoms were also treated with ECPs, and with DZ valence basis sets,23 augmented by polarization d-functions24 for those nitrogen atoms directly coordinated to the metal. A DZV basis was used for hydrogen.25 Auxiliary basis sets to fit the charge density and exchange-correlation potential were of the sizes (5,2) for O with all-electron basis, (3,4) for Fe, (3,3) for C, N, and O atoms with ECPs, and (5,2) or (4,0) for hydrogen with BII and DZV orbital bases, respectively (n,m denotes n s-functions and m spd-shells21,28). An extra iteration without fit of the potential and with larger grids was appended after SCF convergence was reached, to obtain numerically accurate Kohn-Sham orbitals (cf. refs 21 and 28 for further details). Computed absolute 17O nuclear shieldings σ have been converted to relative chemical shifts δ using a shielding value of σ ) +291 ppm for liquid water, derived from the absolute shielding scale of ref 36. 17O electric field gradients (nuclear quadrupole coupling tensors, NQC, given in MHz) were computed from the same

5202 J. Phys. Chem. B, Vol. 104, No. 21, 2000

Kaupp et al.

TABLE 1: Computed Isotropic 17O NMR Chemical Shifts (ppm vs H2Oliq) for Different Oxyheme Models Compared to Experimental Data for the Picket-Fence Porphyrin Complexa model complexb

Obridge

Calculated (i) FeP(Im)(O2)//opt 1805(2130) (i) eclipsed 1993(2400) (ii) FeP(2-MeIm)(O2)// 1767(2076) Fe(TpivPP)(2-MeIm)(O2) (iii) (ii) + point charges 1631(1837) 1716(1964) (iv) Fe(TamP)(Im)(O2)// Fe(TpivPP)(2-MeIm)(O2) (iv) eclipsed 1968(2262) Experimental solid-state 77Kc 1190 solid-state 298 Kc 1600 CH2Cl2 solution, 298 Kd 1739

Oterminal 2420(2915) 2723(3356) 2431(2925) 2111(2430) 2290(2675) 2704(3149) 1967 2017 2507

a Computational data are SOS-DFPT results with UDFT results in parentheses. b See computational details section. c See ref 8. d See ref 7.

TABLE 2: Computed Isotropic 17O NMR Chemical Shifts (ppm vs H2Oliq) for O3 and for the O3(H2O) Complexa O3 calcd O3 exptl O3(H2O) calcd

Obridge

Oterminal

1124(1264) 1032 1033(1128)

1548(1771) 1598 1449(1600)b 1521(1691)c

a Computational data are SOS-DFPT results with UDFT results in parentheses. b Hydrogen-bonded terminal oxygen. c Non-hydrogenbonded terminal oxygen (cf. Figure 2). The 17O shift for the coordinated H2O molecule is computed to be +41.0 (+40.5) ppm.

Kohn-Sham wave functions as the chemical shifts, using the electric-field gradient (EFG) module37 of the deMon-NMR code. Our previous experience gives us confidence that this computational level should provide accurate oxygen NQC tensors.38 NQC tensors and spin densities of the lowest triplet state for model i were computed at the optimized structure for the closedshell singlet (which, however, is similar to that of the triplet, except for a larger Fe-O-O angle18). An unrestricted KohnSham wave function was used for the triplet. 3.Results and Discussion 3.1. Isotropic 17O NMR Chemical Shifts. Table 1 summarizes the isotropic 17O shifts computed for the different model systems discussed above, for both bridging (i.e., directly adjacent to the metal) and terminal positions. Experimental solid-state NMR results of Oldfield et al. at 77 K and at 298 K, and solution NMR data of Gerothanassis et al. at 298 K, for the picket-fence porphyrin system are also included. We first note that, irrespective of the particular model employed, the calculations confirm that the shifts of the dioxygen ligand are at the highest 17O frequencies known for a diamagnetic compound. They also support the previous6,8,39 assignment of the larger shift value to the terminal oxygen atom. The UDFT shifts in parentheses are generally considerably larger than the SOS-DFPT data, which include the Malkin correction27 terms. This is also the typical behavior for main group compounds when a few individual low-energy electronic excitations dominate the paramagnetic shielding contributions.27,28 The magnitude of the correction terms is even larger than found previously31 for the ozone molecule (cf. also Table 2), with its closely related40-42 electronic structure. In view of the good agreement with experiment for the latter case, we believe that the corrected SOS-DFPT results are preferable over

the UDFT data also for the present oxyheme complexes. We will therefore concentrate in the following on the SOS-DFPT data, keeping in mind the inherently large uncertainties associated with such unusually large paramagnetic (deshielding) contributions to the shifts. Note also that the error margins estimated for the isotropic shifts in the solid-state NMR experiments are (200 ppm,8 but they might be even larger due to possible problems in the simulation of the spectra (cf. section 3.3). In going from model i to model ii (i.e., taking the structure from the optimization of the full picket-fence system), the shift of the bridging position is reduced moderately (by ca. 40 ppm), and that of the terminal position increases slightly (by ca. 10 ppm). This indicates that the ruffling of the porphyrin ring, and the slight tilting of the imidazole ligand, due to the steric influences of TpivP substituents and of the extra methyl group on the imidazole ligand influences the oxygen chemical shifts only relatively little (Table 1). In contrast, inclusion of point charges in model iii, to simulate the electrostatic effect of the TpivP substituents, reduces the shifts significantly. As expected, this electrostatic effect is most pronounced for the terminal oxygen atom, as this is in closer vicinity of the TpivP amido groups. However, the point charge model probably overestimates the substituent effect. Upon explicit inclusion of amido groups (Tam) in model iv, the low-frequency shift relative to model ii is still significant (and still largest for Oterminal) but less pronounced than with model iii. We believe model iv to be the most realistic of the four models studied here, and it will form the basis of our comparison with experiment. The situation on the experimental side is far from clear-cut (Table 1): our computed isotropic shifts (SOS-DFPT results for model (iv) might be considered in reasonable agreement with the 298 K solid-state data, overestimating them by ca. 116 ppm and ca. 173 ppm for Obridging and Oterminal, respectively, i.e., within the stated experimental error bars. However, the 17O bridging resonance was reported to exhibit a dramatic lowfrequency shift upon lowering the temperature to 77 K.8 We will consider below whether this large change truly reflects a dynamical effect, or whether it is partly exaggerated by artifacts of the spectra simulations. On the other hand, the 298 K solution data are relatively close to the 298 K solid-state results (and thus to our calculations) for the bridging position but deviate significantly to high frequencies for the terminal oxygen atom. Thus, while our SOS-DFPT data may be considered to be within the range of the experimental data, the large variation of the latter makes a more detailed assessment of the accuracy of the computational results difficult at this point. This indicates a further need for both experimental and theoretical study. Table 2 shows the agreement with experiment that is obtained at the same computational level for the ozone molecule. We expect the accuracy of our computational method to be comparable for the oxyheme models. We also included results for the hydrogen-bonded complex of ozone with one water molecule to probe the effect that hydrogen bonds may have, also for the dioxygen complexes. The hydrogen bonding enhances nuclear shielding, both for the hydrogen-bonded terminal and for the bridging oxygen position (the effect for the remote terminal oxygen is small). A shielding effect due to hydrogen bonding with distal moieties in various modified porphyrin model systems has been inferred by Gerothanassis et al.7 from their solution NMR data. In those cases, the effect of hydrogen bonding was considered shielding by ca. -35 to -40 ppm for the terminal oxygen atoms, i.e., significantly less than computed here for the O3(H2O) complex. The effect was

17O

NMR Studies of Oxyheme Model Complexes

even a slightly deshielding one for Obridging (ca. +10 ppm). These trends appear consistent with weaker and longer hydrogen bonds in these systems. Note that similar influences of hydrogen bonding might be expected for the actual protein environment, due to weak hydrogen bonds with the distal histidine.5,7 Oldfield et al.8 associated the large temperature dependence of δ17O(Obridging) in the solid-state experiments with the rotation of the O2 ligand above the porphyrin plane (similar arguments were previously used to rationalize the temperature dependence of the 57Fe Mo¨ssbauer data43). A ”freezing in” of a conformational substate was considered at 77 K, whereas essentially free rotation should take place at 298 K. We may obtain further insight from recent CPMD simulations on FeP(Im)(O2).19 The room-temperature simulations showed that the O2 ligand prefers a conformation that bisects a N-Fe-N angle within the porphyrin plane. It spends several ps above one quadrant of the plane before ”jumping” over one of the relatively low barriers (ca. 5-6 kJ/mol for FeP(Im)(O2),19 somewhat larger for the picket-fence system18) associated with an ecplipsed position of the O2 ligand above one of the porphyrin Fe-N bonds. Indeed, as shown in Table 1, we compute considerably larger 17O shifts for the eclipsed conformation. This increase is ca. 190 ppm and ca. 300 ppm for bridging and terminal oxygen nuclei, respectively, with model i. It is even ca. 240 ppm and ca. 415 ppm, respectively, with the amido substituents present in model iv. The electrostatic effect of the amido groups enhances the conformational dependence of the shifts, by reducing the shifts additionally at the bisecting conformation, in particular for Oterminal. The larger shifts at the ecplipsed conformations might contribute to the temperature dependence of the isotropic shifts. However, in the simulations these eclipsed positions are populated only sparsely, as the O2 ligand has its largest population density in bisecting conformations.19 Moreover, contributions from the eclipsed conformations should enhance particularly the shift for Oterminal, whereas that of Obridging would be affected less. This contradicts the much larger observed temperature dependence for the latter position (cf. Table 1). We have also computed shifts and field gradients for different relatiVe orientations of the axial O2 and imidazole ligands (data not shown), but we found only minor changes. The Fe-O-O angle has to be considered as another motional degree of freedom. It has been found previously that the competitiveness of open-shell states increases for larger angles.18 Excitations into low-lying excited states should thus enhance the paramagnetic part of the shielding even more upon opening the angle. While the average Fe-O-O angle in the roomtemperature simulation of FeP(Im)(O2) was 123°, only slightly larger than the optimized equilibrium angle of the closed shell singlet at the same computational level (122°), a significant population of configurations with angles above 125° was found.19 Figure 3 shows the computed oxygen chemical shifts as a function of the Fe-O-O angle. Most interestingly, there is a particularly dramatic increase for Obridging above 125°, while the changes for Oterminal are much less pronounced. This would indeed be consistent with a large temperature effect on δ17O(Obridging) but not on δ17O(Oterminal). Of course, other motional effects may have to be considered as well. For example, we have tested the tilt of the Fe-O bond vector relative to the normal of the porphyrin plane. However, within the ca. (5° range of tilt angles sampled in the room-temperature molecular dynamics simulations,19 the bridging and terminal oxygen shifts only deviated by ca. (50 ppm from those at the equilibrium position. Positive tilting of the Fe-O vector, i.e., away from

J. Phys. Chem. B, Vol. 104, No. 21, 2000 5203

Figure 3. Isotropic oxygen chemical shifts (ppm relative to H2Oliq) in model i as a function of Fe-O-O bond angle (all other structural parameters have been kept fixed at their fully optimized values).

TABLE 3: Computed 17O NMR Chemical Shift Tensor Components (ppm vs H2Oliq) for Different Oxyheme Models Compared to Experimental Solid-State NMR Data for the Picket-Fence Porphyrin Complexa Obridging

Oterminal

model

δ11

δ22

δ33

δ11

δ22

δ33

(i) eclipsed (ii) (iii) (iv) (iv) eclipsed exptl 77 Kb exptl 298 Kb O3 O3(H2O)

4241 4749 4138 3815 4021 4706 2300 2650 2247 1963

955 1019 942 865 904 978 1170 2150 829 830

221 211 220 212 223 219 160 0 295 304

5743 6571 5829 4990 5449 6554 4200 2850 2936 2483c 2808d

1292 1365 1257 1158 1213 1331 850 1600 1355 1459c 1393d

223 232 208 184 207 226 850 1600 354 406c 262d

(i)

a Computational data are SOS-DFPT results. b Reference 8. c Hydrogen-bonded terminal oxygen. d Non-hydrogen-bonded terminal oxygen.

the quadrant sampled by the FeO2 moiety, decreased both shifts, whereas negative tilting increased them slightly. 3.2. 17O Chemical Shift Tensors. Table 3 summarizes the computed principal components of the shift tensors for both oxygen positions in the model complexes, compared to both the 77 and 298 K solid-state data for the picket-fence system. Changes between the different models are seen mainly in the dominant, most deshielded component, δ11, whereas δ22 is already less sensitive, and the most shielded δ33 varies very little. In particular, the electrostatic effect of the picket-fence substituents (modeled by point charges in model iii and by amido groups in model iv) affects mainly δ11 of Oterminal. Comparison of our model iv to the 77 K solid-state data, which have been interpreted essentially as rigid-lattice tensors,8 indicates that our δ11 values are considerably larger than the experimental ones (Table 3). This holds in particular for the bridging oxygen, where the computed value is more than 1700 ppm larger than the experimental one, whereas differences in the isotropic shifts are much less pronounced. Computed and experimental (77 K) δ22 and δ33 values for Obridging may be considered to be in reasonable agreement. In contrast, there are striking differences for Oterminal. While experimentally an axially symmetrical tensor has been assumed (i.e., δ22 ) δ33 within the experimental uncertainties8), we find a significantly nonaxial tensor, resembling the situation for ozone (see below). Given the reasonable magnitude of the computed isotropic shifts, and in view of our previous experience, it appears unlikely that the calculations would provide such a poor qualitative

5204 J. Phys. Chem. B, Vol. 104, No. 21, 2000

Kaupp et al.

Figure 4. Computed orientation of the 17O chemical shift tensors relative to the molecular framework. (a) Ozone. (b) Oxyheme model i.

description of the shift tensors. We believe that this may indicate either that (i) even at 77 K the shift tensor for the terminal oxygen is still averaged by significant motion of the O2 ligand, that (ii) the simulations of the spectra were inaccurate, or that (iii) contributions from an open-shell singlet may have to be considered (cf. discussion in section 3.4). Comparing the experimental tensors at 77 and 298 K, we find it difficult to understand why a motional averaging process should actually increase δ11 for Obridging, whereas δ33 is reduced to 0 ppm. For comparison, a temperature increase from 298 K through 393 K has been found to decrease the span of the 17O shift tensor (via a decrease of δ11 and an increase of δ33) for the nitrosyl ligand in Co(OEP)(NO) (OEP ) octaethylporphyrin).9 These are further indications that there are problems in the interpretation of the static solid-state spectra for the oxy picket-fence system. Turning to our reference system ozone, the general pattern of the shift tensors (and their orientation, cf. below) is very similar to the oxyheme case, except for the significantly smaller δ11. The reduced shifts upon binding of a water molecule, for the hydrogen-bonded terminal as well as for the bridging oxygen, are largely due to a reduction in δ11. Figure 4 shows the orientations of the oxygen shift tensors for the oxyheme complexes (Figure 4b; the picture is based on calculations for model i, but the other models give virtually identical orientations), compared to those for the closely related ozone molecule (Figure 4a). In both systems, the dominant principal components δ11 are roughly oriented parallel to the O-O bond. The δ22 component may be termed in-plane perpendicular, whereas δ33 is out-of-plane. Upon closer inspection, differences between the two systems are apparent: while δ11 is aligned almost exactly with the O-O bond (with δ22 perpendicular to it) for both oxygen nuclei in the oxyheme models, in ozone it deviates from the line of the nuclei by ca. 21° for Oterminal and ca. 28° for Obridging. The latter angle is determined by the fact that the O-O-O angle in ozone is exactly bisected by δ22, due to symmetry. Interestingly, the computed tensor orientation for the oxyheme models is essentially identical to those computed for both 15N and 17O shift tensors in Fe-RNO and Co-NO analogue systems.9 The DFT-IGLO approach used here allows a breakdown of the shielding tensors into contributions from individual localized occupied and canonical virtual MOs. The analysis indicates clearly that the orientations of the tensors are sensitive measures of the character of the out-of-plane π* LUMO of the Fe-O-O and O-O-O moieties, respectively (Figure 5). The same type of LUMO has already been implied in previous discussions of bonding and optical spectra for oxyheme complexes.40-42 The paramagnetic contributions to δ11 of the central oxygen nucleus in ozone are dominated by the magnetic-field induced coupling between the oxygen lone pair localized molecular orbital (LMO)

Figure 5. Lowest unoccupied Kohn-Sham MOs (LUMOs) for the closed-shell singlet ground states of (a) ozone and (b) oxyheme model i. HOMO and LUMO energies are -7.2 and -5.69 eV, respectively, for ozone, and -4.27 and -3.40 eV, respectively, for model i.

on this atom and the LUMO. The fact that the LUMO in ozone has identical contributions from both terminal oxygen atoms (Figure 5a) determines the orientation of δ11. Similarly, δ11 of terminal oxygen is dominated by the coupling of the terminal oxygen lone pair LMOs (in this case two of them) with the LUMO. The larger shift for the terminal atoms is due to the fact that two lone-pair LMOs may couple to the LUMO. O-O bonding LMOs provide the main contributions to δ22 (they also contribute moderately to δ11). The results of our analysis are qualitatively similar to but of course quantitatively different from those of the early CHF-IGLO study of the 17O shieldings in ozone by Schindler and Kutzelnigg.44 The situation in the oxyheme models is very similar to that in ozone but modified by the inequality of the iron porphyrin and terminal oxygen sides of the bent Fe-O-O unit. The LUMO again corresponds to an out-of-plane π* MO of the Fe-O-O moiety (Figure 5b). The paramagnetic contributions to δ11 for the bridging oxygen nucleus are dominated by couplings of two occupied LMOs to this LUMO, (a) the inplane oxygen lone pair LMO, and (b) an LMO with in-plane Fe-O bonding character. The computed orientation of δ11 along the O-O bond is the result. The larger shift compared to the bridging oxygen of ozone is largely due to the smaller HOMOLUMO gap: the Kohn-Sham-MO energy difference is 0.81 eV compared to 1.67 eV in ozone (cf. footnotes to Figure 5). The δ11 component for the terminal oxygen atom is again dominated by the coupling of two in-plane oxygen lone pairs with the out-of-plane LUMO. The orientation of δ11 along the O-O bond arises from the asymmetrical nature of this LUMO (Figure 5b). Thus, the orientation of the oxygen shift tensors in such bent end-on dioxygen complexes depends on the composition of the out-of-plane π*-type LUMO. This in turn reflects the magnitude of the Fe-O out-of-plane π-back-bonding, and thus also the weight of different resonance structures in a valence-bond description.42,47 We note that the presence of such orbitals, which are largely centered on the Fe-O-O unit, is also thought to be responsible for the characteristic near-infrared features in the optical spectra of oxyhemes.41,42 3.3. 17O Nuclear Quadrupole Coupling Tensors. Table 4 gives the principal components of the NQC tensors computed for both bridging and terminal oxygen nuclei in the oxyheme model complexes. Again, results for ozone and for the hydrogenbonded ozone-water complex are given for comparison. In the case of ozone, experimental gas-phase data are available,45 and the calculated absolute principal values agree well with these

17O

NMR Studies of Oxyheme Model Complexes

TABLE 4:

17O

J. Phys. Chem. B, Vol. 104, No. 21, 2000 5205

Nuclear Quadrupole Coupling Tensor Components (MHz) Obridge

Oterminal

modela

q11 (ip-⊥)b

q22 (op)c

q33 (ip-|)d

q11 (ip-⊥)b

q22 (op)c

q33 (ip-|)d

(i) eclipsed (ii) (iii) (iv) (iv) eclipsed (i) triplete

+11.47 +11.72 +11.18 +10.71 +10.99 +11.04 -10.09 (op)c

-8.91 -8.63 -8.97 -9.75 -9.32 -8.83 +6.78 (ip-⊥)b

-2.56 -3.09 -2.21 -0.97 -1.67 -2.21 +3.32 (ip-|)d

+17.52 +17.56 +17.67 +16.45 +17.03 +17.30 -13.16 (op)c

-14.36 -14.10 -14.12 -14.80 -14.44 -14.10 +7.19 (ip-⊥)b

-3.16 -3.46 -3.54 -1.65 -2.58 -3.19 +5.97 (ip-|)d

q11 (ip-⊥)b

q22 (op)c

q33 (ip-|)d

q11 (ip-⊥)b

q22 (op)c

q33 (ip-|)d

(i)

O3 calcd O3 exptlf O3(H2O) calcd

+9.78 (-)8.52f +8.47

-8.75 (+)7.46f -8.24

-1.03 (+)1.06f -0.23

+19.01 +19.96 +17.91g +19.08h

-15.00 -15.52 -15.14g -14.78h

-4.01 -4.44 -2.78g -4.30h

See computational details section. b In-plane perpendicular component, except for the triplet state. c Out-of-plane component, except for the triplet state. d In-plane parallel component. e Triplet state computed at the optimized singlet structure (which is close to triplet structure, cf. refs 17 and 18). The orientations of the NQC tensor components differ from those of the singlet state. f Cf. ref 45. It is not clear how the sign of the principal components was determined. gHydrogen-bonded position. h Non-hydrogen bonded position. a

(the signs derived experimentally for Obridging are just opposite to our computed data, possibly due to the use of a different system of axes). The computed tensor components for the ironporphyrin models are very similar to those for ozone. Closer examination shows slightly increased q11 for Obridging but slightly reduced q11 for Oterminal compared to ozone. Thus, not only the shielding tensors but also the electric field gradients reflect the close similarity in the bonding of our model oxyheme complexes and ozone. The main effect of hydrogen bonding in O3(H2O) is a reduction of q11 for the hydrogen-bonded terminal oxygen atom, with only minor changes for the other parameters. The good agreement with experiment for ozone gives us confidence in the NQC tensors for the oxyheme systems computed at the same level. The ratio between the NQCCs for bridging and terminal oxygen nuclei is also consistent with the ratio of the line widths in the solution spectra6,7 (cf. below). Our present results contradict the conclusions of Oldfield et al.,8 who assumed NQCCs of less than 2 MHz in the simulations of the static solid-state 17O NMR spectra. They argued that the use of larger oxygen NQCCs would reduce the quality of the fit. Thus, there is obviously a discrepancy between our quantum chemical results and the previous interpretation of the solidstate spectra. As the quality of the simulations will also affect the reliability of the 17O NMR chemical shift tensors deduced, a reconsideration of the solid-state NMR spectra, or even a redetermination at higher field strengths would appear desirable. The orientations of the NQC tensor components (qii in Table 4) are the same for the oxyheme complexes and for O3, but they differ from those of the shielding tensors (δii in Figure 4): for the central oxygen in ozone, q11 is in the molecular plane, bisecting the O-O-O angle. q22 is out-of-plane, whereas q33 is in the plane (”tangential”). The q33 component for Oterminal is oriented almost exactly along the O-O bond, with q22 again out-of-plane and q11 in-plane perpendicular. For Obridging in the oxyheme complexes, q11 is still in the plane but tilted considerably toward the iron atom, making an angle of ca. 37° with the Fe-O bond. Correspondingly, q33 is almost oriented along the O-O bond. q22 is again out-of-plane. For Oterminal, q11 is in the plane, essentially perpendicular to the O-O bond. The q22 component is out-of-plane, whereas q33 is along the bond. Gerothanassis et al.7 have taken the significant temperature dependence of the 17O solution-NMR line widths as an indication of a conformational excitation mechanism, in which

TABLE 5: Dependence of 17ONQCCs on Fe-O-O Anglea ∠Fe-O-O 110° 115° 120.8° 125° 130° 140°

Obridging q11 η 10.74 11.02 11.47 11.73 12.19 13.18

0.63 0.60 0.55 0.55 0.50 0.44

Oterminal q11 η 17.47 17.44 17.52 17.44 17.57 17.85

0.67 0.66 0.64 0.66 0.64 0.64

q11(term)/ q11(bridge)b

χterm/χbridge.c

1.63 1.58 1.53 1.49 1.44 1.35

2.68 2.56 2.40 2.30 2.18 1.96

a The angle was varied for model i, starting from the fully optimized structure. b Ratio between NQCCs for terminal (term) and bridging oxygen nuclei. c Ratio between χ ) q211(1 + η2/3) for terminal and bridging oxygen nuclei.

two conformers with very different electric field gradients at oxygen coexist at higher temperatures, whereas one preferred conformer freezes out at low temperature. A similar explanation had been proposed previously for the temperature dependence of the 57Fe electric field gradients.43 However, in the latter case more recent experimental and computational studies of analogue systems9 have provided no indication of a conformational excitation mechanism. Instead, fast rotation of the Fe-O2 unit was suggested to account for the temperature dependence of the EFGs at iron.9 We have computed the 17O NQC tensors for different orientations of the O2 ligand in models i and iv relative to the porphyrin plane (entries in Table 4) or relative to the axial imidazole (data not shown). Changes in the computed field gradients are minor for all of these calculations. The computations thus provide no support for the suggested conformational excitation mechanism.7 Table 5 shows the influence of a variation in the Fe-O-O angle on the computed NQCCs (q11) and asymmetry parameters η ) (q33 - q22)/q11 for model i. While the variations in both q11 and η for Oterminal are minor and nonmonotonic, there is a notable increase of q11 and a decrease of η for Obridging with increasing angle. Gerothanassis et al. assumed nuclear quadrupole relaxation to provide the dominant relaxation mechanism.7 If this is true, and given that the rotational correlation times for Oterminal and Obridging are the same, the ratio between the line widths for the two positions should be proportional to the ratio between the quantities χ ) q211(1 + η2/3) for both positions.7 The ratio χterm./χbridg. decreases significantly with increasing Fe-O-O angle (Table 5). The average angle increases with increasing temperature, due

5206 J. Phys. Chem. B, Vol. 104, No. 21, 2000 to the anharmonicity of the bending mode. This motion would therefore be expected to lead to a monotonic decrease of the line width ratio upon increasing the temperature. This is indeed what has been observed.7 However, experimentally the individual line widths for both positions show a nonmonotonic behavior (a decrease followed by an increase at still higher temperatures). Such a behavior could not be deduced from the numbers in Table 5, which would suggest a monotonic increase for Obridging and very little change for Oterminal. If we exclude a significant influence of the other motional degrees of freedom on the NQCCs (see above) we have to conclude that either (a) relaxation mechanisms other than nuclear quadrupole relaxation may contribute to the line widths (e.g., the unusually large shift anisotropy46) or (b) contributions from an open-shell singlet state at lower temperatures may not be excluded (see section 3.4). Table 4 includes also the NQCCs for the triplet state of model i, computed at the optimized structure for the closed-shell singlet. Compared to the closed-shell state, q11 and q22 are interchanged, both for Oterminal and for Obridging. The in-plane perpendicular component has decreased for both positions (moderately so for Obridging, dramatically for Oterminal), whereas the out-of-plane component (q11 for the triplet) is slightly more negative for Obridging but slightly less negative for Oterminal. The in-plane parallel component (q33 in both states) has switched from a small negative value to a somewhat larger positive one. As a consequence, the asymmetry of the NQC tensor is considerably smaller for the triplet state. All of these differences are consistent with the removal of one electron from the inplane π-system, and its addition to the out-of-plane π*-LUMO (cf. Figure 5). As the Weiss-type open-shell singlet has the same orbital configuration as the lowest triplet state, we expect the overall electron distribution and thus the electric field gradients to be similar for the two open-shell states. The open-shell singlet should thus exhibit an NQC tensor very different from the closed-shell singlet. Therefore, even a small admixture of openshell singlet character into the ground-state wave function is expected to have a pronounced influence on the line widths. This should be kept in mind in the discussion of the nature of the ground state of oxyheme systems, to which we now turn. 3.4. Implications of the 17O NMR Spectra for the Nature of the Ground State of Oxyheme Complexes. The 17O NMR spectra of oxyheme complexes, in particular the absence of any apparent Fermi-contact shifts,5-8 have been taken as clear indications of a closed-shell singlet ground state, as suggested originally by Pauling10 (note that the “ozone-like” bonding model suggested by Goddard et al.40 differs from the Pauling closed-shell model only in its preferred decomposition into fragment states47). Unfortunately, things are not that straightforward. The absence of contact shifts48 does of course effectively rule out a triplet ground state or a thermally populated low-lying excited triplet state, consistent with early magnetic measurements.10a However, a Weiss-type11 antiferromagnetically coupled open-shell singlet state would not exhibit contact shifts. Recent unrestricted gradient-corrected DFT calculations17,18 could not discriminate between a closed-shell and an open-shell singlet ground state. More precisely, the spin-unrestricted calculations pertained to a strongly spin-polarized brokensymmetry state corresponding to the singlet contaminated with triplet character. The energy of the open-shell singlet after spin projection should in fact be even lower than that of the brokensymmetry state.49 The very small energy differences involved (the broken-symmetry state was computed to be ca. 6 kJ/mol below the closed-shell singlet, whereas the triplet was found to

Kaupp et al. be ca. 13 kJ/mol above the latter17,18) and the difficulties in describing accurately the open-shell singlet wave function for a molecule of the size of an oxyheme, make it currently impossible to pinpoint the ground state of the system clearly by such energy calculations alone. Our calculations of 17O chemical shift tensors have been based on a closed-shell singlet wave function. The reasonable agreement of the computed isotropic shifts with experiment might thus be taken as confirmation of a closed-shell ground state (the solid-state 17O NMR spectra of hemoglobin and myoglobin have a similar appearance but a less favorable signal-to-noise ratio compared to those of synthetic model complexes; they have also been interpreted in terms of a closed-shell ground state8). On the other hand, we do not know what the shifts should be for the open-shell singlet, as we presently have no reliable way of estimating the shifts in such a case. There are various inconsistencies between different experiments (e.g., the unusual temperature dependence of the solid-state shifts for Obridging or the difference between solution and solid-state shifts for Oterminal), and between experiment and computation (e.g., the much larger computed shift anisotropies compared to the supposedly rigid-lattice tensors measured at 77 K, the very low nuclear quadrupole coupling constants estimated from the solidstate spectra, and the nonmonotonic temperature dependence of the line widths in the solution spectra), which still appear to leave room for alternative interpretations. One of them could be a temperature-dependent interconversion or mixing between open and closed-shell singlet states, as suggested earlier in the context of Mo¨ssbauer spectra.50 This possibility would be consistent with the very small energy differences computed17 previously between these two states. We must therefore conclude that even the present, successful calculations of 17O chemical shifts and NQCCs with a closed-shell wave function do not rule out completely that an open-shell singlet state may contribute to some extent. Conclusions In the present work we have used density functional calculations to highlight the interdependence of structure, dynamics, and oxygen NMR parameters of oxyheme systems. The unusually large 17O shifts found experimentally have been confirmed in calculations based on a closed-shell state. The analogy of the electronic structure of oxyheme complexes with that of the ozone molecule has been strengthened by the closely similar orientations of the 17O chemical shift and nuclear quadrupole coupling tensors. The out-of-plane π*-type lowest unoccupied MO of the Fe-O-O unit of the oxyheme models dominates the appearance of the 17O shift tensors in a way very similar to that found for the ozone analogue. Slight differences in the orientation, and particularly in the magnitude of the tensor components, are mainly related to the asymmetry of the Fe-O-O moiety, and to the smaller HOMO-LUMO gap in the oxyheme complexes. By combining knowledge drawn from recent ab initio molecular dynamics simulations on oxyhemes with calculations of NMR parameters, we have also gained insight into the effects of molecular dynamics on the spectra. The comparison of computed and experimental shift tensor components suggests that, in contrast to previous interpretations,8 even at 77 K the rotation of the O2 ligand above the porphyrin plane is not frozen in, consistent with the very small computed barriers for this rotation.18,19 Our predictions of 17O nuclear quadrupole coupling constants (ca. +11.5 MHz and ca. +17.5 MHz for bridging and terminal oxygen, respectively) agree well with the line width

17O

NMR Studies of Oxyheme Model Complexes

ratio of the solution spectra,6,7 and with experiment as well as computation for the related ozone molecule. This contradicts the very small (