Search tips
Search criteria 


Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
Theor Chem Acc. Author manuscript; available in PMC 2011 March 1.
Published in final edited form as:
Theor Chem Acc. 2010 March; 125(3-6): 305–317.
doi:  10.1007/s00214-009-0566-4
PMCID: PMC2863024

Quantum cluster size and solvent polarity effects on the geometries and Mössbauer properties of the active site model for ribonucleotide reductase intermediate X: a density functional theory study


In studying the properties of metalloproteins using ab initio quantum mechanical methods, one has to focus on the calculations on the active site. The bulk protein and solvent environment is often neglected, or is treated as a continuum dielectric medium with a certain dielectric constant. The size of the quantum cluster of the active site chosen for calculations can vary by including only the first-shell ligands which are directly bound to the metal centers, or including also the second-shell residues which are adjacent to and normally have H-bonding interactions with the first-shell ligands, or by including also further hydrogen bonding residues. It is not well understood how the size of the quantum cluster and the value of the dielectric constant chosen for the calculations will influence the calculated properties. In this paper, we have studied three models (A, B, and C) of different sizes for the active site of the ribonucleotide reductase intermediate X, using density functional theory (DFT) OPBE functional with broken-symmetry methodology. Each model is studied in gas-phase and in the conductor-like screening (COSMO) solvation model with different dielectric constants ε = 4, 10, 20, and 80, respectively. All the calculated Fe-ligand geometries, Heisenberg J coupling constants, and the Mössbauer isomer shifts, quadrupole splittings, and the 57Fe, 1H, and 17O hyperfine tensors are compared. We find that the calculated isomer shifts are very stable. They are virtually unchanged with respect to the size of the cluster and the dielectric constant of the environment. On the other hand, certain Fe-ligand distances are sensitive to both the size of the cluster and the value of ε. ε = 4, which is normally used for the protein environment, appears too small when studying the diiron active site geometry with only the first-shell ligands as seen by comparisons with larger models.

Keywords: Ribonucleotide reductase intermediate X, Active site, COSMO, Dielectric constant, DFT, OPBE, Mössbauer properties

1 Introduction

Quantum mechanical calculations have long been used to study protein active site properties. The initial geometry of the protein is normally taken or constructed from the X-ray crystal structure. Irrespective of the kind of ab initio quantum mechanical method (e.g., Hartree Fock (HF), configuration interaction (CI), many-body perturbation theory, and density functional theory (DFT), etc.) used, in practice it is not possible to treat the whole protein structure quantum-mechanically. The main focus is on the active site region for the quantum calculations. For metalloproteins, the “quantum region” or the “quantum cluster” includes at least the first coordination sphere (or shell) around the metal site(s). There are different ways, like QM/MM [13] and QM/Electrostatics (Poisson–Boltzmann) [4, 5], to include part of the protein and solvent effects during the electronic structure and property calculations [6, 7]. However, for simplicity, very often the bulk protein and solvent environment is either neglected or considered as a continuum dielectric medium around the quantum cluster using a certain solvation model with a fixed dielectric constant (ε) [8, 9].

In recent years, our group has successfully applied the COSMO [1013] (conductor-like screening model) solvation model implemented in the Amsterdam Density Functional (ADF) package [1416] together with the DFT broken-symmetry and spin-projection methodology [1719] to study the geometries, relative energies, oxidation states, redox potentials, net-spin populations, pKa's, as well as Mössbauer isomer shifts, quadrupole splittings, and various magnetic hyperfine properties of the Fe–oxo, Fe–S, and Fe–Cu active sites of ribonucleotide reductase (RNR) [2024], methane monooxygenase (MMO) [25, 26], nitrogenase [6, 27], cytochrome c oxidase [28], and other iron–sulfur proteins [29].

Many active site complexes in electron-transfer metalloproteins are multiply charged. Therefore, the interaction energy between the active site and the surrounding protein and solvent environment is very large. Including solvation effects during quantum mechanical calculations (including geometry optimization) is expected to produce more realistic predictions than the pure “gas-phase” calculations. However, most current continuum solvation models allow only one fixed dielectric constant representing the polarity of the solvent. How to choose a proper dielectric constant for the protein plus solvent environment is still an open question. Although the dielectric value ε = 4 is commonly used for the protein interior, since this is the value of the dielectric constants of crystalline and polymeric amides[30] and dry protein and peptide powders [3134], many studies show that higher effective dielectric constant values (4–30) for protein interiors are needed in reproducing the pKa values of certain internal ionizable groups [3441]. The necessary high protein dielectric constant in calculations has been rationalized from the water penetration [34, 37], the reaction field of bulk solvent [42, 43], solvation by permanent dipoles in the protein [44], the structural reorganization [40], fluctuations of surface charged groups [39], and the influence of nonlinear electrostatics [45]. How the value of the dielectric constant applied in the continuum solvation model will influence the calculated metalloprotein active site properties, especially the geometries and Mössbauer properties, is not fully investigated. On the other hand, the first-shell ligands in the metalloproteins we studied also interact with the second (and then the third) shell residues through H-bonding effects. It is also not fully understood to what extent the calculated active site properties vary with the size of the quantum cluster, which increases by inclusion of the 1st, 2nd, and additional shells of H-bonding residues around the metal centers.

In the current paper, we will take the candidate active site model of the class-I ribonucleotide reductase intermediate state X as an example to see how the size of the quantum cluster and the different values of the dielectric constant in COSMO solvation model will influence the calculated geometries, Heisenberg J coupling constants, Mössbauer and ENDOR (electron–nuclear double resonance) properties of the Fe–oxo center. By comparing these results, we will examine how to choose the dielectric constant of the continuum environment for the quantum active site models of metalloproteins with different sizes.

2 Ribonucleotide reductase intermediate X

Ribonucleotide reductases (RNRs) catalyze the reduction of ribonucleotides to deoxyribonucleotides providing the required building blocks for DNA replication and repair. Class-I RNR, found in eukaryotes as well as in some prokaryotes and viruses, consists of a homodimer of two dissimilar protein subunits: R1 and R2 in an α2β2 architecture. The ribonucleotide-to-deoxyribonucleotide reactions occur by a long range radical (or proton-coupled-electron-transfer) propagation mechanism initiated by a fairly stable tyrosine (Tyr122 in E. coli) radical (“the pilot light”). It is the subunit R2 that contains a dinuclear iron cluster which initially generates and stabilizes this tyrosine radical. The subunit R1 contains the substrate binding site, and catalyzes the dehydroxylation of the 2′-hydroxyl group of the ribose ring. Once the tyrosine radical (Tyr122·) is formed by the diiron center of R2, the catalytic reaction is started by a long range radical transfer (or proton-coupled-electron-transfer) to generate a thiyl radical on cysteine 439 of subunit R1, which then performs the nucleotide reduction [46, 47].

Presently, the detailed structure of the active diferric form of R2 which contains the tyrosine radical is not known. This neutral tyrosine radical has been identified in the oxidized metalloprotein form and is stable for days at room temperature [48]. Once this radical is lost, the active form can be regenerated by a complicated sequence of steps involving changes in oxidation state and structural rearrangement with coupled electron and proton transfers. First the resting oxidized diferric met form of R2 (R2ox, see Fig. 1) is reduced by two electrons from a reductase protein to the diferrous form, R2red (see Fig. 1). Next, a molecular oxygen (O2) binds to the diiron center of R2red and an electron is transferred from Trp48 to one of the iron sites. Afterwards, a transient high-oxidation R2 intermediate state, named X, is kinetically and spectroscopically observed. X is evidently the species which regenerates the tyrosine radical. RNR-X has captured the attention of many researchers over the past 10 years to elucidate its chemical and structural nature [2024, 4969]. A combination of Q band ENDOR and Mössbauer data on Y122F-R2 indicate the iron centers of X are high spin Fe(III) (S = 5/2) and high spin Fe(IV) (S = 2) sites that antiferromagnetically couple to give an Stotal = 1/2 ground state [56].

Fig. 1
The X-ray active site structures of the reduced diferrous (R2red) and oxidized(met) diferric (R2ox(met)) RNR-R2 from E. coli [83, 85]

In the past 5 years, we have been studying the properties of a set of active site model clusters for RNR-X using broken-symmetry DFT methods [2024], and have compared them with the available experimental data, including Mössbauer [56], 57Fe, 1H, 17O2, and H217O ENDOR [59, 61], EXAFS (extended X-ray absorption fine structure) [60], and MCD (magnetic circular dichroism) [62]. Based on the detailed analysis and comparisons, we have proposed that the experimentally observed, in particular the ENDOR observed RNR-X active site contains two μ–oxo bridges, plus one terminal water which binds to Fe1(III) (where Fe1 is the iron site closer to Tyr122) and also H-bonds to both side chains of Asp84 and Glu238, and one bidentate carboxylate group from the side chain of Glu115 (see the colored picture in Fig. 2) [21]. However, the active site models to which we applied Mössbauer and hyperfine calculations together with COSMO solvation have usually contained only the first shell residues plus the Tyr122 (or mutant Phe122) side chains [21, 23]. (Some gas-phase calculations were performed on larger models [23, 24].) Considering the charged groups in the active site and the H-bonding environment, the dielectric constant for COSMO calculations was set to 80 (as for water) [21], which is much higher than the commonly used ε = 4 for protein interior. Now we will study how the calculated results change if we vary the dielectric constant and increase the size of the active site model.

Fig. 2
Small RNR-X model-A with labels and its optimized central Fe-ligand distances (Å) and the calculated properties [Mössbauer isomer shift (δ), quadrupole splittings (ΔEQ), J coupling constant, and spin projected energy ( ...

3 Computational methodology

3.1 DFT calculations

All density functional spin-unrestricted calculations have been performed using the ADF package [1416]. The parametrization of Vosko, Wilk, and Nusair (VWN) [70] is used for the local density approximation term. In our previous RNR-X studies, the corrections of Perdew and Wang (PW91) [71] were used for the nonlocal exchange and correlation terms. However, our recent studies on MMO intermediate Q show that PW91 potential overestimates the Fe-ligand covalencies for some of the Fe(IV)–Fe(IV) active site models [25]. We therefore also applied OPBE [7274] functional in our MMO intermediates Q and P studies [25, 26]. OPBE is the combination of Handy's optimized exchange (OPTX) [74] and PBE correlation (PBEc) functionals [72, 73]. Swart et al. [75, 76] have reported systematic studies on the performance of several exchange-correlation functionals for various properties of different systems. They found that the OPBE potential correctly predicted the spin states for all the iron complexes they tested, and in general performs well in predicting other properties. In the current study, we will apply only OPBE for the nonlocal exchange and correlation terms, and will compare with our previous corresponding PW91 calculations.

The geometries of the quantum clusters are optimized in gas-phase and in the COSMO solvation model with dielectric constant ε = 4.0, 10.0, 20.0, and 80.0, respectively. Optimizations are stopped when the maximum Cartesian gradient is less than 0.01 Hartree/Å, the estimated uncertainty in the Cartesian coordinates is less than 0.01 Å, and the change in energy is less than 0.001 Hartree. In COSMO, the quantum cluster is embedded in a molecular shaped cavity surrounded by a continuum dielectric medium. The van der Waals radii for atoms Fe, C, O, N, and H were taken as 1.5, 1.7, 1.4, 1.55, and 1.2 Å, respectively. The probe radius for the contact surface between the cluster and solvent was set to 2.0 Å. The triple-ζ polarization (TZP) Slater-type basis sets with frozen cores [C(1 s), N(1 s), O(1 s) and Fe(1s,2s,2p) are frozen] are applied for geometry optimizations.

Since our previous calculations show that the site Fe1, which is closer to Tyr122 should be the Fe(III) site in X [21], we will only study the Fe1(III)Fe2(IV) high-spin antiferromagnetically (AF) coupled {S1 = 5/2, S2 = 2, Stotal = 1/2} state for all calculations on the model clusters. Usually the AF spin-coupled state cannot be obtained directly from the normal DFT calculations in ADF. As in previous work, we represent the AF spin-coupled state in DFT by a “broken-symmetry” (BS) state [1719], where a spin-unrestricted determinant is constructed in which one of the Fe site has spin-up electrons as majority spin and the other site has spin-down electrons. To obtain this broken-symmetry solution, first we construct a ferromagnetically (F) spin-coupled (Stotal = 9/2) determinant, where the spins on both irons are aligned in a parallel fashion. Then we rotate the spin vector located on either atom Fe1 or atom Fe2 by interchanging the α and β fit density blocks on the site Fe1 or Fe2 from the output file TAPE21 created by this F-coupled calculation in ADF. Using the modified TAPE21 as a restart file and reading the starting spin density from there, we then obtain the expected BS state through single-point energy calculation or geometry optimization.

After geometry optimization, a high-spin F-coupled single-point energy calculation (in COSMO) with all-electron TZP Slater-type basis sets (i.e., without frozen core approximation) is performed at the BS optimized geometry, and the energy EF is obtained. The corresponding BS state calculation then follows to obtain the electron density [ρ(0)] and the electric field gradient (EFG) at the Fe nucleus, the A-tensors (electron-nuclear magnetic hyperfine interaction), and the BS state energy EBS.

3.2 Mössbauer isomer shift and quadrupole splitting calculations

The Mössbauer isomer shifts δ are calculated based on ρ(0):


In our previous studies [25, 77], the parameters α and C have been fitted separately for the Fe2+,2.5+ and Fe2.5+,3+,3.5+,4+ complexes for PW91, OPBE, and OLYP, with all-electron TZP Slater type basis sets. For the Fe2.5+,3+,3.5+,4+ complexes, we have obtained A = 11,877.0, α = −0.312, and C = 0.373 mm s−1 for OPBE method.

For calculating the Mössbauer quadrupole splittings (ΔEQ), the EFG tensors V are diagonalized and the eigenvalues are reordered so that |Vzz| ≥ |Vxx| ≥ |Vyy|. The asymmetry parameter η is defined as


Then the ΔEQ for 57Fe of the nuclear excited state (I = 3/2) can be calculated as


where e is the electrical charge of a positive electron, Q is the nuclear quadrupole moment (0.15 barns) of Fe [78].

3.3 Pure-spin ground state energy estimation

The BS state of the X models obtained from DFT calculations is a mixture of pure spin states. When the following Heisenberg Hamiltonian H (with Heisenberg coupling J) is applicable,

H=2JS1[center dot]S2

the energy difference between the F-coupling (Stotal = Smax = S1 + S2 = 5/2 +2 = 9/2) and the BS (Stotal = Smin = |S1S2| = 5/2 − 2 = 1/2) states can be described by

EFEBS=4JS1S2=20J(for RNRX)

(Note that for a completely delocalized-mixed valence dimer, a more general spin Hamiltonian is H = −2J0S1·S2 ± B(Stotal + 1/2), where B is the resonance delocalization parameter [79, 80]. If the total spin Stotal is small and B is not large, the resonance stabilization energy −B(Stotal + 1/2) is often neglected. In this case, this term is largely quenched by vibronic effects or other environmental effects, including solvation.) J is then obtained from Eq. 5, and the pure-spin ground state energy E0 for the particular spin state (S1, S2) coupled to Smin according to the broken-symmetry geometry is estimated as:

E0=EF+JSmax(Smax+1)JSmin(Smin+1)=EF+24J=EBS+4J(for RNRX)

Note that for more accurate calculations, both BS and F-coupled high-spin state geometries need to be optimized. The structure with the minimum E0 can be obtained by extrapolating the geometries between the optimized BS and F-coupled geometries. For current X models, since the J coupling constant is small, the diiron centers are weakly coupled, and the BS and F-coupled states do not differ much. For simplicity, the model geometries are only optimized at the BS state, and an F-coupled high-spin single-point energy calculation is performed at the BS optimized geometry to get the EF energy. The J and E0 values are then calculated from Eqs. 5 and 6.

3.4 Hyperfine A-tensor calculations

The ligand 57Fe, 1H, and 17O hyperfine coupling constants are predicted based on the A-tensor calculations in ADF, which assumes that there is only one unpaired electron in the system when the total spin Stotal = 1/2 and z component Ms = 1/2. The spin-orbit coupling contributions to the A-tensors are neglected. For the present systems with high spin AF coupled sites, we therefore need to rescale the ADF-obtained A-tensors by the spin projection coupling factors |KA/2SA| for Fe(III) (KA = 7/3, SA = 5/2) and |KB/2SB| for Fe(IV) (KB = −4/3, SB = 2) [24, 81, 82]. Absolute values of the coupling factors are used here, since the broken symmetry state carries the proper A-tensor sign. For the two bridging oxo atoms, the coupling factor is taken as the average of |KA/2SA| and |KB/2SB|.

4 Models and results

Three models with different sizes for our proposed RNR-X active site structure are studied here. They are named Model-A, Model-B, and Model-C in order of increasing cluster size. The details of the models and the calculated geometries, Heisenberg J coupling constants, energies, and Mössbauer isomer shift and quadrupole splitting properties are given in the Sects. 4.1–4.3. The predicted 1H, 17O, and 57Fe hyperfine coupling constants for these models are compared in Sect. 4.4.

4.1 Small Model-A

Model-A is shown in Fig. 2. This model contains only the first-shell ligands plus the Tyr122 sidechain. The initial positions of the first shell ligand sidechains in this model are taken from chain A of the oxidized RNR-R2 (met) X-ray crystal structure (PDB code: 1RIB) [83], by breaking the Cβ−Cα bonds and adding a linking hydrogen atom along the Cβ−Cα direction to fill the open valence of the terminal carbon atom [84]. No preliminary alteration or optimization procedure is applied to the X-ray crystal structure. As shown in Figs. 1 and and2,2, these first-shell ligands are Asp84, His118, Glu115, His241, Glu204, and Glu238. The μ–oxo bridge (O1) which lies between His118 and His241 and the terminal water which binds with Fe1 in the active site of R2ox(met) are maintained and the additional oxygen (O2) which lies between Asp84 and Glu204 is added to construct the initial structure of Model-A. The orientation of the Asp84 sidechain is modified so that one of the oxygen atoms in the carboxylate group H-bonds to both the terminal water and Tyr122.

Model-A is then optimized in gas-phase and in COSMO solvation model with ε = 4, 10, 20, and 80, respectively. The linking H atom on Tyr122 is fixed during geometry optimizations. Mössbauer and hyperfine calculations are then performed at these optimized geometries. The calculated Fe-ligand distances, J coupling constants, spin-projected energies (E0), Mössbauer isomer shift (δ) and quadrupole splitting (ΔEQ) properties are given in Fig. 2 and compared with available Mössbauer experimental data [56].

From gas-phase to ε = 80, the Fe1–Fe2 distance gradually increases by 0.05 Å. The largest change for the Fe-ligand distances is from gas-phase to ε = 4. When ε > 10, most of the distances are relatively stable, with changes less than 0.02 Å. However, the Fe1–N1(His118) distance is an exception. The interaction between Fe1 and N1(His118) is extremely sensitive to the polarity of the environment. This distance is decreased by 0.07 Å from gas-phase to ε = 4, and is further decreased by 0.03 Å when ε increased to 10. From ε = 10 to ε = 80, it continues to decrease by 0.17 Å. The Fe2–N2(His241) distance is also sensitive to the polarity of the solvent, but not as much as Fe1–N1. In general, from gas-phase to solvent environment, the distances of Fe1–Fe2, Fe1–O1, Fe1–O2, Fe1–O3, Fe2–O7 increase, and all other Fe-ligand distances decrease (see the labels in the colored structural model in Fig. 2). However, with a further decrease of the Fe1–N1 and Fe2–N2 distances (from ε = 20 to ε = 80), the Fe–O-carboxylate distances of Fe1–O4, Fe1–O5, and Fe2–O6 may increase a little.

In gas-phase, the electron contributions to the Fe1(III)–HOMO orbital mainly come from Fe1(28.46%), O1(25.24), O2(21.41), O4(7.24%), O5(2.13%), and Fe2(1.0%). In ε = 80, these contributions change to Fe1(24.44%), O1(24.46), O2(17.59), O4(8.02%), Fe2(1.06%), and N1(2.58%). Therefore, the decreasing of the Fe1–N1 distance from gas-phase to ε = 80 correlates with the slight increasing (from <1.0 to 2.58%) of the electron contribution of N1 to the Fe1(III)–HOMO molecular orbital. The potential energy surface along the Fe1–N1 direction should be very shallow.

The absolute value of the J coupling constant decreases with increasing solvent polarity. The maximum change of J is 11 cm−1. From gas-phase to ε = 4, the electronic energy of Model-A drops sharply by 2.75 eV, consistent with the net cluster charge of −1. With increasing solvent polarity, the energy change is larger in the low dielectric range (0.74 eV from ε = 4 to ε = 10) than in the high dielectric range (only 0.27 eV from ε = 20 to ε = 80). This trend is roughly consistent with a Born solvation model.

It is valuable to observe that the calculated Mössbauer isomer shift values are very stable, virtually unchanged with variations in the polarity of the environment. Our calculations predict δ(Fe1) = 0.55 (or 0.54) mm s−1 and δ(Fe2) = 0.23 (or 0.22) mm s−1 for Model-A in all dielectric conditions. These results are very consistent with the experimental observations (0.56 and 0.26 mm s−1). Our previous PW91 calculations with ε = 80 also predicted very similar isomer shifts (0.57 and 0.22 mm s−1) for the two iron sites in a similar small X model [21]. By contrast, the calculated quadrupole splittings do change with the polarity of the environment although not very significantly. The asymmetry parameter η of ΔEQ(Fe1) gets larger in polar media, and ΔEQ(Fe1) changes sign when ε = 80. The observed quadrupole splitting for Fe(III) (corresponding to Fe1 here) is −0.9 mm s−1 with η = 0.5, and the ΔEQ for Fe(IV) (Fe2 here) is −0.6 mm s−1 with η = 2.7 (η > 1). If we trace back and reorder the eigenvalues to have |Vzz| ≥ |Vxx| ≥ |Vyy| for ΔEQ{Fe(IV)}exp, we will have ΔEQ{Fe(IV)}exp = −0.6 mm s−1 with η = 0.08. Considering both the signs and the absolute values, none of the calculated quadrupole splittings for Model-A agree with the observed ones. Previous PW91 calculations on different RNR-X models drew the same conclusion [21]. We mainly focus on the absolute values of the quadrupole splittings when comparing our calculations with experiment, since our other studies also show that if η is close to 1, the sign of the calculated quadrupole splitting can also change with the chosen atomic basis sets and computational methods. Later we will see that the sign may also vary by increasing the quantum cluster size. For RNR-X Model-A, the calculated ΔEQ absolute values are smaller than the observed ones, however, the relative orders of |ΔEQ(Fe1)| > |ΔEQ(Fe2)| are reproduced, which is consistent with our previous PW91 calculations [21].

Now we move onto a larger quantum cluster Model-B, and see how the calculated results of Model-A in different dielectric constants will be related to different forms of Model-B.

4.2 Model-B

Compared with Model-A, additional second and third shell H-bonding residue sidechains are included in Model-B (see Fig. 3). These residues are: Asp237, Trp48, Gln43, and Trp111. The linking H atoms of these residues and Tyr122 are fixed during geometry optimizations. Several water molecules are found around the carboxylate groups of Asp237 and Glu204. Since there are already three H-bonding interactions with Asp237 (from His118, Trp48, and Gln43, respectively), we therefore only include one water molecule (Wat621) which H-bonds to Glu204 in Model-B.

Fig. 3
Model-B of RNR-X with more second and third shell H-bonding residues (as labeled in the structural picture) and its central optimized Fe-ligand distances (Å) and the calculated properties [Mössbauer isomer shift (δ), quadrupole ...

Calculated results for Model-B are shown in Fig. 3. The geometry of Model-B also changes with increasing solvent polarity. Although the tendencies of the variations are similar to Model-A, the magnitudes are much smaller. From gas-phase to ε = 80, the Fe1–Fe2 and Fe1–N1(His118) distances change by 0.02 and 0.11 Å, respectively, in Model-B, comparing with the corresponding changes of 0.05 and 0.27 Å in Model-A. Using these two distances as indicators, the geometry of the diiron center of Model-B in gas-phase is similar to the geometry of Model-A in the environment with 20 < ε < 80. Therefore, the influence of the second and third shell H-bonding residues on the geometry of the diiron center cannot be neglected. To reasonably reflect such influence when studying the small active site model of X with only first coordination shell, the dielectric constant in the solvation model is suggested to be 20 or larger.

In gas-phase, the calculated J coupling constants are almost the same (−136 and −135 cm−1) for Models A and B. However, the J value of Model-B obtained with ε = 4 (−128 cm−1) lies between the J values of Model-A in ε = 20 (−131 cm−1) and ε = 80 (−125 cm−1). This reflects that calculations of Model-A with high dielectric constant (>20) are equivalent to the results of Model-B in low dielectric environment.

Just as observed in Model-A calculations, the polarity of the bulk solvent has very little influence (up to 0.01 mm s−1) on the calculated Mössbauer isomer shifts for Model-B. The influence of the extra 2nd and 3rd shell residues in Model-B on the calculated isomer shifts are also small (no more than 0.01 mm s−1). Model-A and Model-B, therefore, yield almost the same isomer shifts in all dielectric environment.

For Model-A, the absolute values of the calculated quadrupole splittings decrease with increasing solvent polarity. From gas-phase to ε = 80, the calculated |ΔEQ(Fe1)| and |ΔEQ(Fe2)| change from 0.66 and 0.65 mm s−1 to 0.38 and 0.34 mm s−1. A similar trend but with smaller changes was also observed for Model-B ΔEQ calculations from gas-phase to ε = 10. However, the ΔEQ values become nearly constant for Model-B when ε ≥ 10.

4.3 Model-C

As mentioned above, several water molecules are found around the second shell H-bonding residues in the diferric met RNR-R2 X-ray crystal structure. In addition to Wat621 which was included in Model-B (Fig. 3), there is another water molecule, Wat520, also close to the carboxylate group of Glu204. Three water molecules labeled as Wat612, Wat727, and Wat734 are found in the vicinity of the carboxylate group of Asp237. Two of them, Wat612 and Wat727, H-bond directly to the two oxygen atoms of Asp237, and the third one, Wat734, is in H-bonding distance with Wat612 (see Fig. 4). Next we will increase the size of Model-B by including these four extra water molecules in the quantum cluster and see how these third shell explicit H-bonding water molecules will influence the calculated active site properties. The structure of Model-C and its calculated results are given in Fig. 4.

Fig. 4
Model-C of RNR-X with four more water molecules as labeled comparing with Model-B, and its central optimized Fe-ligand distances (Å) and the calculated properties [Mössbauer isomer shift (δ), quadrupole splittings (ΔEQ ...

In gas-phase, these water molecules have little effect on the central active site structure. The corresponding Fe–Fe and Fe-ligand bond lengths are similar to each other in Models B and C. However, these water molecules further strengthen the polarization when the clusters are put in the solvent model. From gas-phase to ε = 4, the distance of Fe1–N1(His118) in Model-B is only decreased by 0.03 Å from 2.491 to 2.462 Å. However, in Model-C, this distance is shortened by 0.13 Å from 2.494 to 2.369 Å. Since the Fe1–N1 distance in the X models is very sensitive to both the size and the polarity of the bulk solvent, in Fig. 5 we compare how this distance changes in Models A–C with increasing solvent polarity from gas-phase to ε = 80.

Fig. 5
Comparison of the Fe1–N1(His118) (see Fig. 2) distance in Models A, B, and C, with increasing solvent polarity from gas-phase to ε = 80

For all three models, a sharp decrease of the distance Fe1–N1 is observed from gas-phase to ε = 10. When ε > 10, the Fe1–N1 distance in small Model-A continues decreasing and meets with the line of Model-B around ε = 80. Previously, we have used ε = 80 in studying the small X models with the similar size as Model-A [21]. Now it seems that it was a good choice since the geometry for Model-A obtained with ε = 80 is very similar to that of Model-B obtained with ε = 10. For Models B and C, the Fe1–N1 distance changes very slightly (around 0.01 Å) from ε = 10 to ε = 80, suggesting that for larger active site models including 2nd and 3rd shell H-bonding residues, ε = 10 is near the convergence point. Note that there is always a gap (about 0.05 Å) between the two Fe1–N1 distances in Models B and C in the region of (10 < ε < 80), therefore the structural influence of the 3rd-shell explicit water molecules on the diiron centers, especially on this Fe1–N1 distance cannot be replaced by the continuum solvation model with any high dielectric constant. For studying the geometrical properties of the active sites of metalloproteins, whenever possible, it is important to include more outer shell H-bonding residues in the quantum cluster model.

Again as seen in Models A and B, Mössbauer isomer shift calculations for Model-C yield essentially the same results in gas-phase and in COSMO models, and also the same as those for Models A and B. The RNR-X active site has a high-spin AF-coupled Fe(III)Fe(IV) center. Presently, no literature iron model compounds have such oxidation and spin states simultaneously. High-spin Fe(IV) is also rare in model compounds. Some Fe–oxo complexes show the isomer shift of Fe(III) site around 0.55–0.56 mm s−1, [25, 77] but most are lower, 0.40–0.50 mm s−1 [25, 77]. Therefore, it is not sure that experimentally the isomer shifts are environmentally sensitive or insensitive, since this influence has to be separated from that of the ligand environment and of the Fe–O–Fe bridge.

The predicted J coupling constant for the larger Model-C is almost converged when ε > 10, and the quadrupole splittings are almost unchanged in the region of 4 < ε < 80. These show that the influence of the continuum solvation effect on the active site properties decreases with increasing size of the quantum cluster.

4.4 Hyperfine A-tensor calculations

All 1H, 17O, and 57Fe A-tensors are calculated and rescaled based on spin projection coefficients (see Sect. 3.4) at all the optimized geometries of the three RNR-X models obtained in gas-phase and in COSMO solvation model with ε = 4, 10, 20 and 80, respectively. It turns out that the three X models in all solvation conditions yield similar hyperfine A-tensors, regardless of the cluster size and the dielectric constant. There are some obvious changes for the 17O A-tensors of each model going from gas-phase to ε = 4 solvation environment. We therefore only present in Table 1 the hyperfine results for Model-A in gas-phase, and for all three models in the solvation medium with ε = 10. The calculated A-tensors are also compared with the available experimental data.

Table 1
Calculated hyperfine coupling constants (MHz) for E. coli RNR-R2-X Models A, B, and C in Fe1(III)Fe2(IV){S1 = 5/2, S2 = 2, Stotal = 1/2} state, and compared with experimental results

The exchangeable proton signals observed from the CW and pulsed Q-band [1, 2] H ENDOR measurements were assigned to two protons which were associated with a terminal water or a twofold disordered terminal hydroxo bound to Fe(III) [59]. The observed [1, 2] H ENDOR hyperfine signals are nearly axial, consistent only with terminal HxO, not with bridging hydroxo. Our predicted 1H A-tensors (see Ht1 and Ht2) are almost the same for the three X models in different solvation conditions, and are in agreement with the experiment. The calculated isotropic values (Aiso) are near to zero, and the A3 values (~19 MHz for Ht1 and ~17 MHz for Ht2) are very close to the observed data (20.5 and 17.6 MHz). Very similar results are also obtained from previous PW91 calculations [21, 23].

From gas-phase to solvation environment (see results for Model-A in Table 1), about 2.5–4.9 MHz changes are observed for the A-tensor components A1 of the oxo bridges O1 and O2, and for A3 of the bridge O2 (atom labels are given in Fig. 2). Then, within solvation model with ε [set membership] [4, 80], all calculations for Models A–C yield similar corresponding A-tensors for atoms O1, O2, and O3. The calculated A-tensors for O3 in the terminal water are in good agreement with the observed 17O hyperfine constants for the terminal oxygen (Ot) [61]. Our previous PW91 calculations (with ε = 80) yield (−15.68, −17.87, −32.91 MHz) for O3 [21, 23], which are even closer to the ENDOR experimental data (−17.0, −20.5, −34.0 MHz) [61]. Note that only the relative signs of the three principal values were determined in the 17O2 and H217O ENDOR experiments. We have set the signs according to our calculations, respecting the relative signs set by experiment. We also reordered the 17O A-tensor components for convenience so that |A3| is the largest.

Experimentally only one bridging (br) oxo 17O signal was found [61]. There are two bridging oxo atoms O1 and O2 in our RNR-X models. However, for each model in different solvation environments, both the calculated A-tensor principal values and the principal axes of the two bridging oxo atoms are very close to each other. This may explain why only one Obr signal has been observed from the 17O2 ENDOR experiments.

The observed 57Fe hyperfine spectra from ENDOR and magnetic Mössbauer spectroscopies show a highly isotropic hyperfine tensor for the high-spin Fe1(III) site and moderate anisotropy for the high-spin Fe2(IV) site. Theoretically it is very difficult to predict the isotropic hyperfine coupling constants of metal centers. Our previous PW91 calculated 57Fe isotropic hyperfine coupling constants (Aiso) for small model X are −32.40 MHz for Fe1 and 11.17 MHz for Fe2 [21, 23], which are less than half of the observed ones (−73.2(Fe1) and 33.7(Fe2) MHz) [56]. Current OPBE calculated 57Fe isotropic hyperfine coupling constants are closer to the experimental data. The values of Aiso for Fe1 are about 60% of the experimental value. All current X models in different solvation conditions yield very similar Aiso values for both Fe1 (around −44 MHz) and Fe2 (~14 MHz). Normally the 57Fe anisotropic components (Aaniso) can be more accurately predicted by the DFT calculations. Current OPBE calculations also yield very similar 57Fe Aaniso components for Fe1 and Fe2, respectively, for all X models with different solvent polarities. The predicted 57Fe Aaniso components are in very good agreement with the experimental data.

5 Conclusions

Three active site models (A, B, and C) for class-I RNR intermediate X have been studied using broken-symmetry density functional theory OPBE potential incorporated with the COSMO solvation model. These models have the same kind of diiron core structure, but with different sizes. Model-A contains the first Fe-ligand coordination shell plus the Tyr122 phenol ring. Model-B includes more second and third shell H-bonding residue sidechains. Model-C is larger than Model-B by including more explicit H-bonding water molecules around the third shell, which were found in the X-ray crystal structure of the oxidized diferric RNR. The purpose of this paper is to study how the size of the diiron active site model and the polarity of the environment represented by the solvent continuum model will influence the calculated geometries, Heisenberg J coupling constants, Mössbauer and hyperfine properties of the Fe–oxo active site center. We optimized the geometries of these three RNR-X active site models in gas-phase, and in COSMO solvation model with dielectric constant ε = 4, 10, 20, and 80, respectively, and computed the Heisenberg J coupling constant, Mössbauer isomer shift, quadrupole splitting, and various hyperfine A-tensor properties at the optimized geometries.

It is found that the calculated properties of Mössbauer isomer shift, 57Fe, and the 1H and 17O (for the terminal water ligand O3) hyperfine coupling constants are very stable. Almost the same values of these properties are obtained for different sizes of the quantum clusters in different solvation conditions. Therefore, the small first-shell active site model is adequate for studying these properties.

The calculated Fe-ligand distances, especially the Fe1–N1(His118) distance, Mössbauer quadrupole splittings, Heisenberg J coupling constant are sensitive to both the quantum cluster size and the polarity of the environment. In larger Models B and C, these properties are almost converged when ε > 10. If the quantum cluster for the diiron active site includes both the first and second shell residues (or residue sidechains), ε = 10 can be a good compromise representing both the outer-shell H-bonding effect and the bulk solvation environment. However, the explicit water molecules around the 3rd coordination shell still have strong influence to the diiron central geometry, which cannot be reproduced by the continuum solvent model with any dielectric constant (see Fig. 5).

In order to obtain equivalent geometry and J coupling constant for Model-A which are similar to those obtained from Models B and C in lower dielectric solvents, the dielectric constant for Model-A calculations should not be less than 20. ε = 4 for the bulk solvation environment seems too small when only the first shell residues are included in the quantum cluster of the diiron active site model.


We thank NIH for financial support (GM43278 to L.N.). The support of computer resources of the Scripps Research Institute is also gratefully acknowledged.


Dedicated to Professor Sandor Suhai on the occasion of his 65th birthday and published as part of the Suhai Festschrift Issue.


1. Senn HM, Thiel W. Atomistic approaches in modern biology: from quantum chemistry to molecular simulations. Vol. 268. Springer; Berlin: 2007. pp. 173–290.
2. Klahn M, Braun-Sand S, Rosta E, et al. J Phys Chem B. 2005;109:15645–15650. doi: 10.1021/jp0521757. [PMC free article] [PubMed] [Cross Ref]
3. Lill SON, Siegbahn PEM. Biochemistry. 2009;48:1056–1066. doi: 10.1021/bi801218n. [PubMed] [Cross Ref]
4. Asthagiri D, Dillet V, Liu TQ, et al. J Am Chem Soc. 2002;124:10225–10235. doi: 10.1021/ja020046n. [PubMed] [Cross Ref]
5. Li J, Nelson MR, Peng CY, et al. J Phys Chem A. 1998;102:6311–6324. doi: 10.1021/jp980753w. [Cross Ref]
6. Noodleman L, Lovell T, Han WG, et al. Chem Rev. 2004;104:459–508. doi: 10.1021/cr020625a. [PubMed] [Cross Ref]
7. Orozco M, Luque FJ. Chem Rev. 2000;100:4187–4225. doi: 10.1021/cr990052a. [PubMed] [Cross Ref]
8. Tomasi J, Persico M. Chem Rev. 1994;94:2027–2094. doi: 10.1021/cr00031a013. [Cross Ref]
9. Cramer CJ, Truhlar DG. Chem Rev. 1999;99:2161–2200. doi: 10.1021/cr960149m. [PubMed] [Cross Ref]
10. Klamt A, Schüürmann G. J Chem Soc Perkin Trans. 1993;2:799–805. doi: 10.1039/p29930000799. [Cross Ref]
11. Klamt A. J Phys Chem. 1995;99:2224–2235. doi: 10.1021/j100007a062. [Cross Ref]
12. Klamt A, Jonas V. J Chem Phys. 1996;105:9972–9981. doi: 10.1063/1.472829. [Cross Ref]
13. Pye CC, Ziegler T. Theor Chem Acc. 1999;101:396–408. doi: 10.1007/s002140050457. [Cross Ref]
14. Vrije Universiteit; Amsterdam, The Netherlands: ADF 2006 01; SCM, Theoretical chemistry.
15. te Velde G, Bickelhaupt FM, Baerends EJ, et al. J Comput Chem. 2001;22:931–967. doi: 10.1002/jcc.1056. [Cross Ref]
16. Guerra CF, Snijders JG, te Velde G, et al. Theor Chem Acc. 1998;99:391–403.
17. Noodleman L. J Chem Phys. 1981;74:5737–5743. doi: 10.1063/1.440939. [Cross Ref]
18. Noodleman L, Case DA. Adv Inorg Chem. 1992;38:423–470. doi: 10.1016/S0898-8838(08)60070-7. [Cross Ref]
19. Noodleman L, Lovell T, Han WG, et al. In: Comprehensive coordination chemistry II, from biology to nanotechnology. Lever AB, editor. Elsevier; Amsterdam: 2003.
20. Han WG, Lovell T, Liu T, et al. Inorg Chem. 2003;42:2751–2758. doi: 10.1021/ic020465l. [PubMed] [Cross Ref]
21. Han WG, Liu T, Lovell T, et al. J Am Chem Soc. 2005;127:15778–15790. doi: 10.1021/ja050904q. [PubMed] [Cross Ref]
22. Han WG, Liu TQ, Lovell T, et al. Inorg Chem. 2006;45:8533–8542. doi: 10.1021/ic060566+. [PubMed] [Cross Ref]
23. Han WG, Liu T, Lovell T, et al. J Inorg Biochem. 2006;100:771–779. doi: 10.1016/j.jinorgbio.2006.01.032. [PubMed] [Cross Ref]
24. Han WG, Lovell T, Liu T, et al. Inorg Chem. 2004;43:613–621. doi: 10.1021/ic0206443. [PubMed] [Cross Ref]
25. Han WG, Noodleman L. Inorg Chim Acta. 2008;361:973–986. doi: 10.1016/j.ica.2007.06.007. [PMC free article] [PubMed] [Cross Ref]
26. Han WG, Noodleman L. Inorg Chem. 2008;47:2975–2986. doi: 10.1021/ic701194b. [PubMed] [Cross Ref]
27. Pelmenschikov V, Case DA, Noodleman L. Inorg Chem. 2008;47:6162–6172. doi: 10.1021/ic7022743. [PMC free article] [PubMed] [Cross Ref]
28. Fee JA, Case DA, Noodleman L. J Am Chem Soc. 2008;130:15002–15021. doi: 10.1021/ja803112w. [PMC free article] [PubMed] [Cross Ref]
29. Torres RA, Lovell T, Noodleman L, et al. J Am Chem Soc. 2003;125:1923–1936. doi: 10.1021/ja0211104. [PubMed] [Cross Ref]
30. Gregg EC. Handbook of chemistry and physics. Chemical Rubber; Cleveland: 1976.
31. Harvey SC, Hoekstra P. J Phys Chem. 1972;76:2987–2994. doi: 10.1021/j100665a011. [PubMed] [Cross Ref]
32. Bone S, Pethig R. J Mol Biol. 1982;157:571–575. doi: 10.1016/0022-2836(82)90477-6. [PubMed] [Cross Ref]
33. Bone S, Pethig R. J Mol Biol. 1985;181:323–326. doi: 10.1016/0022-2836(85)90096-8. [PubMed] [Cross Ref]
34. Dwyer JJ, Gittis AG, Karp DA, et al. Biophys J. 2000;79:1610–1620. doi: 10.1016/S0006-3495(00)76411-3. [PubMed] [Cross Ref]
35. Sham YY, Muegge I, Warshel A. Biophys J. 1998;74:1744–1753. doi: 10.1016/S0006-3495(98)77885-3. [PubMed] [Cross Ref]
36. Bashford D, Karplus M. Biochemistry. 1990;29:10219–10225. doi: 10.1021/bi00496a010. [PubMed] [Cross Ref]
37. Fitch CA, Karp DA, Lee KK, et al. Biophys J. 2002;82:3289–3304. doi: 10.1016/S0006-3495(02)75670-1. [PubMed] [Cross Ref]
38. Antosiewicz J, McCammon JA, Gilson MK. Biochemistry. 1996;35:7819–7833. doi: 10.1021/bi9601565. [PubMed] [Cross Ref]
39. Simonson T, Brooks CL. J Am Chem Soc. 1996;118:8452–8458. doi: 10.1021/ja960884f. [Cross Ref]
40. Karp DA, Gittis AG, Stahley MR, et al. Biophys J. 2007;92:2041–2053. doi: 10.1529/biophysj.106.090266. [PubMed] [Cross Ref]
41. Harms MJ, Schlessman JL, Chimenti MS, et al. Protein Sci. 2008;17:833–845. doi: 10.1110/ps.073397708. [PubMed] [Cross Ref]
42. King G, Lee FS, Warshel A. J Chem Phys. 1991;95:4366–4377. doi: 10.1063/1.461760. [Cross Ref]
43. Loffler G, Schreiber H, Steinhauser O. J Mol Biol. 1997;270:520–534. doi: 10.1006/jmbi.1997.1130. [PubMed] [Cross Ref]
44. Warshel A, Aqvist J, Creighton S. Proc Natl Acad Sci USA. 1989;86:5820–5824. doi: 10.1073/pnas.86.15.5820. [PubMed] [Cross Ref]
45. Gong HP, Hocky G, Freed KF. Proc Natl Acad Sci USA. 2008;105:11146–11151. doi: 10.1073/pnas.0804506105. [PubMed] [Cross Ref]
46. Fontecave M. Cell Mol Life Sci. 1998;54:684–695. doi: 10.1007/s000180050195. [PubMed] [Cross Ref]
47. Stubbe JA, vanderDonk WA. Chem Biol. 1995;2:793–801. doi: 10.1016/1074-5521(95)90084-5. [PubMed] [Cross Ref]
48. Wallar BJ, Lipscomb JD. Chem Rev. 1996;96:2625–2658. doi: 10.1021/cr9500489. [PubMed] [Cross Ref]
49. Bollinger JM, Tong WH, Ravi N, et al. J Am Chem Soc. 1994;116:8015–8023. doi: 10.1021/ja00097a008. [Cross Ref]
50. Bollinger JM, Tong WH, Ravi N, et al. J Am Chem Soc. 1994;116:8024–8032. doi: 10.1021/ja00097a009. [Cross Ref]
51. Bollinger JM, Jr, Edmondson DE, Huynh BH, et al. Science. 1991;253:292–298. doi: 10.1126/science.1650033. [PubMed] [Cross Ref]
52. Bollinger JM, Jr, Tong WH, Ravi N, et al. In: Methods in enzymology. Klinman JP, editor. Academic Press; New York: 1995.
53. Bollinger JM, Jr, Stubbe J, Huynh BH, et al. J Am Chem Soc. 1991;113:6289–6291. doi: 10.1021/ja00016a066. [Cross Ref]
54. Ravi N, Bollinger JM, Jr, Huynh BH, et al. J Am Chem Soc. 1994;116:8007–8014. doi: 10.1021/ja00097a007. [Cross Ref]
55. Ravi N, Bominaar EL. Inorg Chem. 1995;34:1040–1043. doi: 10.1021/ic00109a008. [Cross Ref]
56. Sturgeon BE, Burdi D, Chen SX, et al. J Am Chem Soc. 1996;118:7551–7557. doi: 10.1021/ja960399k. [Cross Ref]
57. Burdi D, Sturgeon BE, Tong WH, et al. J Am Chem Soc. 1996;118:281–282. doi: 10.1021/ja952651e. [Cross Ref]
58. Veselov A, Scholes CP. Inorg Chem. 1996;35:3702–3705. doi: 10.1021/ic951544i. [Cross Ref]
59. Willems JP, Lee HI, Burdi D, et al. J Am Chem Soc. 1997;119:9816–9824. doi: 10.1021/ja9709942. [Cross Ref]
60. Riggs-Gelasco PJ, Shu LJ, Chen SX, et al. J Am Chem Soc. 1998;120:849–860. doi: 10.1021/ja9718230. [Cross Ref]
61. Burdi D, Willems JP, Riggs-Gelasco P, et al. J Am Chem Soc. 1998;120:12910–12919. doi: 10.1021/ja9824270. [Cross Ref]
62. Mitić N, Saleh L, Schenk G, et al. J Am Chem Soc. 2003;125:11200–11201. doi: 10.1021/ja036556e. [PubMed] [Cross Ref]
63. Baldwin J, Krebs C, Ley BA, et al. J Am Chem Soc. 2000;122:12195–12206. doi: 10.1021/ja001278u. [Cross Ref]
64. Krebs C, Chen SX, Baldwin J, et al. J Am Chem Soc. 2000;122:12207–12219. doi: 10.1021/ja001279m. [Cross Ref]
65. Siegbahn PEM. Inorg Chem. 1999;38:2880–2889. doi: 10.1021/ic981332w. [PubMed] [Cross Ref]
66. Siegbahn EM. Q Rev Biophys. 2003;36:91–145. doi: 10.1017/S0033583502003827. [PubMed] [Cross Ref]
67. Younker JM, Krest CM, Jiang W, et al. J Am Chem Soc. 2008;130:15022–15027. doi: 10.1021/ja804365e. [PubMed] [Cross Ref]
68. Jiang W, Yun D, Saleh L, et al. Biochemistry. 2008;47:13736–13744. doi: 10.1021/bi8017625. [PMC free article] [PubMed] [Cross Ref]
69. Jiang W, Saleh L, Barr EW, et al. Biochemistry. 2008;47:8477–8484. doi: 10.1021/bi800881m. [PMC free article] [PubMed] [Cross Ref]
70. Vosko SH, Wilk L, Nusair M. Can J Phys. 1980;58:1200–1211.
71. Perdew JP, Chevary JA, Vosko SH, et al. Phys Rev B. 1992;46:6671–6687. doi: 10.1103/PhysRevB.46.6671. [PubMed] [Cross Ref]
72. Perdew JP, Burke K, Ernzerhof M. Phys Rev Lett. 1996;77:3865–3868. doi: 10.1103/PhysRevLett.77.3865. [PubMed] [Cross Ref]
73. Perdew JP, Burke K, Ernzerhof M. Phys Rev Lett. 1997;78:1396. doi: 10.1103/PhysRevLett.78.1396. [Cross Ref]
74. Handy NC, Cohen AJ. Mol Phys. 2001;99:403–412. doi: 10.1080/00268970010018431. [Cross Ref]
75. Swart M, Groenhof AR, Ehlers AW, et al. J Phys Chem A. 2004;108:5479–5483. doi: 10.1021/jp049043i. [Cross Ref]
76. Swart M, Ehlers AW, Lammertsma K. Mol Phys. 2004;102:2467–2474. doi: 10.1080/0026897042000275017. [Cross Ref]
77. Han WG, Liu TQ, Lovell T, et al. J Comput Chem. 2006;27:1292–1306. doi: 10.1002/jcc.20402. [PubMed] [Cross Ref]
78. Martinez-Pinedo G, Schwerdtfeger P, Caurier E, et al. Phys Rev Lett. 2001;87:062701(1–4). [PubMed]
79. Mouesca JM, Chen JL, Noodleman L, et al. J Am Chem Soc. 1994;116:11898–11914. doi: 10.1021/ja00105a033. [Cross Ref]
80. Zhao XG, Richardson WH, Chen JL, et al. Inorg Chem. 1997;36:1198–1217. doi: 10.1021/ic9514307. [PubMed] [Cross Ref]
81. Noodleman L, Chen JL, Case DA, et al. In: Nuclear magnetic resonance of paramagnetic macromolecules. La Mar GN, editor. Kluwer; Netherland: 1995.
82. Sinnecker S, Neese F, Noodleman L, et al. J Am Chem Soc. 2004;126:2613–2622. doi: 10.1021/ja0390202. [PubMed] [Cross Ref]
83. Nordlund P, Eklund H. J Mol Biol. 1993;232:123–164. doi: 10.1006/jmbi.1993.1374. [PubMed] [Cross Ref]
84. Han WG, Tajkhorshid E, Suhai S. J Biomol Struct Dyn. 1999;16:1019–1032. [PubMed]
85. Logan DT, Su XD, Aberg A, et al. Structure. 1996;4:1053–1064. doi: 10.1016/S0969-2126(96)00112-8. [PubMed] [Cross Ref]