On the Accuracy of Calculated Reduction Potentials of Selected Group 8 (Fe, Ru, and Os) Octahedral Complexes

The theoretical calculations of reduction potentials for the [M­(H₂O)₆]^2+/3+, [M­(NH₃)₆]^2+/3+, [M­(<i>en</i>)₃]^2+/3+, [M­(<i>bipy</i>)₃]^2+/3+, [M­(CN)₆]^4–/3–, and [MCl₆]^4–/3– systems (M = Fe, Os, Ru) were carried out. The DFT­(PBE)/def2-TZVP//DFT­(PBE)/def2-SVP quantum chemical method was employed to obtain presumably accurate ionization energies, whereas the conductor-like screening model for real solvents (COSMO-RS) was selected as the most suitable method for calculations of solvation energies of the oxidized and reduced forms of the studied species. It has been shown that COSMO-RS may overcome problems related to directionality of hydrogen bonds in the second solvation sphere that previously led to errors of ∼1 V for the [Ru­(H₂O)₆]²⁺ complex employing PCM-like models. Thus, most of the values for (2+) → (3+) oxidations are now within 0.1–0.2 V from the experimental data, once the anticipated spin–orbit coupling effects in Os complexes (downshifting the calculated reduction potentials by ∼0.3 V) are taken into account. The robustness of the DFT­(PBE)/COSMO-RS computational protocol is further verified by showing that reduction potentials obtained for selected 2+/3+ redox pairs with and without the inclusion of explicit second-sphere water molecules are almost identical. At the same time, it must be admitted that the calculated values of reduction potentials for systems involving quadruple charged species, exemplified here by [M­(CN)₆]^4–/3– and [MCl₆]^4–/3– redox pairs, might still not be within the grasp of contemporary solvation models, possibly due to the large values of solvation energies of their reduced (4−) forms that are in the range of 700–750 kcal mol^–1 (30–33 eV) and possibly larger errors inherent in their calculations. Finally, a comparison is made with M06-L//SMD computational protocol, which is also shown to correct some of the deficiencies of previous PCM models

Predicting the Stability Constants of Metal-Ion Complexes from First Principles

The most important experimental quantity describing the thermodynamics of metal-ion binding with various (in)­organic ligands, or biomolecules, is the stability constant of the complex (β). In principle, it can be calculated as the free-energy change associated with the metal-ion complexation, i.e., its uptake from the solution under standard conditions. Because this process is associated with the interactions of charged species, large values of interaction and solvation energies are in general involved. Using the standard thermodynamic cycle (in vacuo complexation and solvation/desolvation of the reference state and of the resulting complexes), one usually subtracts values of several hundreds of kilocalories per mole to obtain final results on the order of units or tens of kilocalories per mole. In this work, we use density functional theory and Møller–Plesset second-order perturbation theory calculations together with the conductor-like screening model for realistic solvation to calculate the stability constants of selected complexes[M­(NH₃)₄]²⁺, [M­(NH₃)₄(H₂O)₂]²⁺, [M­(Imi)­(H₂O)₅]²⁺, [M­(H₂O)₃(His)]⁺, [M­(H₂O)₄(Cys)], [M­(H₂O)₃(Cys)], [M­(CH₃COO)­(H₂O)₃]⁺, [M­(CH₃COO)­(H₂O)₅]⁺, [M­(SCH₂COO)₂]^2–with eight divalent metal ions (Mn²⁺, Fe²⁺, Co²⁺, Ni²⁺, Cu²⁺, Zn²⁺, Cd²⁺, and Hg²⁺). Using the currently available computational protocols, we show that it is possible to achieve a <i>relative</i> accuracy of 2–4 kcal·mol^–1 (1–3 orders of magnitude in β). However, because most of the computed values are affected by metal- and ligand-dependent systematic shifts, the accuracy of the “absolute” (uncorrected) values is generally lower. For metal-dependent systematic shifts, we propose the specific values to be used for the given metal ion and current protocol. At the same time, we argue that ligand-dependent shifts (which cannot be easily removed) do not influence the metal-ion selectivity of the particular site, and therefore it can be computed to within 2 kcal·mol^–1 average accuracy. Finally, a critical discussion is presented that aims at potential caveats that one may encounter in theoretical predictions of the stability constants and highlights the perspective that theoretical calculations may become both competitive and complementary tools to experimental measurements

Predicting the Stability Constants of Metal-Ion Complexes from First Principles

The most important experimental quantity describing the thermodynamics of metal-ion binding with various (in)­organic ligands, or biomolecules, is the stability constant of the complex (β). In principle, it can be calculated as the free-energy change associated with the metal-ion complexation, i.e., its uptake from the solution under standard conditions. Because this process is associated with the interactions of charged species, large values of interaction and solvation energies are in general involved. Using the standard thermodynamic cycle (in vacuo complexation and solvation/desolvation of the reference state and of the resulting complexes), one usually subtracts values of several hundreds of kilocalories per mole to obtain final results on the order of units or tens of kilocalories per mole. In this work, we use density functional theory and Møller–Plesset second-order perturbation theory calculations together with the conductor-like screening model for realistic solvation to calculate the stability constants of selected complexes[M­(NH₃)₄]²⁺, [M­(NH₃)₄(H₂O)₂]²⁺, [M­(Imi)­(H₂O)₅]²⁺, [M­(H₂O)₃(His)]⁺, [M­(H₂O)₄(Cys)], [M­(H₂O)₃(Cys)], [M­(CH₃COO)­(H₂O)₃]⁺, [M­(CH₃COO)­(H₂O)₅]⁺, [M­(SCH₂COO)₂]^2–with eight divalent metal ions (Mn²⁺, Fe²⁺, Co²⁺, Ni²⁺, Cu²⁺, Zn²⁺, Cd²⁺, and Hg²⁺). Using the currently available computational protocols, we show that it is possible to achieve a <i>relative</i> accuracy of 2–4 kcal·mol^–1 (1–3 orders of magnitude in β). However, because most of the computed values are affected by metal- and ligand-dependent systematic shifts, the accuracy of the “absolute” (uncorrected) values is generally lower. For metal-dependent systematic shifts, we propose the specific values to be used for the given metal ion and current protocol. At the same time, we argue that ligand-dependent shifts (which cannot be easily removed) do not influence the metal-ion selectivity of the particular site, and therefore it can be computed to within 2 kcal·mol^–1 average accuracy. Finally, a critical discussion is presented that aims at potential caveats that one may encounter in theoretical predictions of the stability constants and highlights the perspective that theoretical calculations may become both competitive and complementary tools to experimental measurements

Accurate Prediction of One-Electron Reduction Potentials in Aqueous Solution by Variable-Temperature H‑Atom Addition/Abstraction Methodology

A robust and efficient theoretical approach for calculation of the reduction potentials of charged species in aqueous solution is presented. Within this approach, the reduction potential of a charged complex (with a charge |<i>n|</i> ≥ 2) is probed by means of the reduction potential of its neutralized (protonated/deprotonated) cognate, employing one or several H-atom addition/abstraction thermodynamic cycles. This includes a separation of one-electron reduction from protonation/deprotonation through the temperature dependence. The accuracy of the method has been assessed for the set of 15 transition-metal complexes that are considered as highly challenging systems for computational electrochemistry. Unlike the standard computational protocol(s), the presented approach yields results that are in excellent agreement with experimental electrochemical data. Last but not least, the applicability and limitations of the approach are thoroughly discussed

Reduction Pathways of 2,4,6-Trinitrotoluene: An Electrochemical and Theoretical Study

The reduction pathways of trinitrotoluene are studied using electrochemical and computational methods. The electrochemical reduction of three nitro groups in 2,4,6-trinitrotoluene (TNT) is characterized by three major reduction peaks in cyclic voltammograms at the peak potentials of −0.310, −0.463, and −0.629 V vs a normal hydrogen electrode (NHE). The second and third peaks coincide with the two peaks observed for the 2-amino-4,6-dinitrotoluene (at the potentials of −0.475 and −0.627 V vs NHE), whereas the two peaks in the 4-amino-2,6-dinitrotoluene voltammograms appear at −0.537 and −0.623 V and deviate more significantly from the corresponding two peaks in 2,4,6-trinitrotoluene. It suggests that the first NO₂ group reduced in the overall process is the one in <i>ortho</i> position with respect to the CH₃ group. Analogously, the 2,6-diamino-4-nitrotoluene exhibits a reduction peak at −0.629 V, almost identical to the third and second reduction peaks of 2,4,6-trinitrotoluene and 2-amino-4,6-dinitrotoluene, respectively. Since the other isomer, 2,4-diamino-6-nitrotoluene, exhibits a reduction peak at −0.712 V, we conclude that the second reduction occurs also in the <i>ortho</i> position with respect to the methyl group. Most of these observations are corroborated by quantum chemical calculations, which yielded reduction potentials in a good agreement with the experimental values (in relative scale). Thus, studying in detail all of the possible protonation and redox states in the reduction of the first nitro group and the key steps in the reduction of the second and third nitro groups, we have obtained a comprehensive and detailed picture of the mechanism of the full 18<i>e</i>^–/18H⁺ reduction of TNT. Last but not least, the calculations have shown that the thermodynamic stabilities of (isomeric) neutral radical species (<b>X </b>+<b> </b><i><b>e</b></i>^<b>–</b><b> </b>+<b> H</b>^<b>+</b>)presumably the regioselectivity-determining steps in the 6<i>e</i>^–/6H⁺ reductions of the individual NO₂ groupsare within 2 kJ·mol^–1 (i.e., comparable to RT). Therefore, the course of the reduction can be governed by the effect of the surroundings, such as the enzymatic environment, and a different regioselectivity can be observed under biological conditions

Computational Electrochemistry as a Reliable Probe of Experimentally Elusive Mononuclear Nonheme Iron Species

Despite the growing number of reported Fe^IVO complexes, an unambiguous experimental characterization of their redox properties, such as one-electron reduction potentials, remains a challenging task. To this aim, we describe an efficient and straightforward theoretical protocol for accurate calculations of redox potentials and calibrate the protocol on a set of diverse 37 mononuclear nonheme iron (NHFe) redox couples. It is shown that the methodology, further applied to a set of 10 Fe^IVO species, not only serves for near-quantitative predictions of reduction potentials, but also is an elegant tool for interpretation of the experimental electrochemical data. The general need for such a computational methodology is illustrated on the prototypical example of the (N4Py)­Fe^IVO complex, whose electrochemistry has been studied for many years and still raises many questions

Accurate Prediction of One-Electron Reduction Potentials in Aqueous Solution by Variable-Temperature H‑Atom Addition/Abstraction Methodology

A robust and efficient theoretical approach for calculation of the reduction potentials of charged species in aqueous solution is presented. Within this approach, the reduction potential of a charged complex (with a charge |<i>n|</i> ≥ 2) is probed by means of the reduction potential of its neutralized (protonated/deprotonated) cognate, employing one or several H-atom addition/abstraction thermodynamic cycles. This includes a separation of one-electron reduction from protonation/deprotonation through the temperature dependence. The accuracy of the method has been assessed for the set of 15 transition-metal complexes that are considered as highly challenging systems for computational electrochemistry. Unlike the standard computational protocol(s), the presented approach yields results that are in excellent agreement with experimental electrochemical data. Last but not least, the applicability and limitations of the approach are thoroughly discussed

Computational Electrochemistry as a Reliable Probe of Experimentally Elusive Mononuclear Nonheme Iron Species

Despite the growing number of reported Fe^IVO complexes, an unambiguous experimental characterization of their redox properties, such as one-electron reduction potentials, remains a challenging task. To this aim, we describe an efficient and straightforward theoretical protocol for accurate calculations of redox potentials and calibrate the protocol on a set of diverse 37 mononuclear nonheme iron (NHFe) redox couples. It is shown that the methodology, further applied to a set of 10 Fe^IVO species, not only serves for near-quantitative predictions of reduction potentials, but also is an elegant tool for interpretation of the experimental electrochemical data. The general need for such a computational methodology is illustrated on the prototypical example of the (N4Py)­Fe^IVO complex, whose electrochemistry has been studied for many years and still raises many questions

Reduction Pathways of 2,4,6-Trinitrotoluene: An Electrochemical and Theoretical Study

The reduction pathways of trinitrotoluene are studied using electrochemical and computational methods. The electrochemical reduction of three nitro groups in 2,4,6-trinitrotoluene (TNT) is characterized by three major reduction peaks in cyclic voltammograms at the peak potentials of −0.310, −0.463, and −0.629 V vs a normal hydrogen electrode (NHE). The second and third peaks coincide with the two peaks observed for the 2-amino-4,6-dinitrotoluene (at the potentials of −0.475 and −0.627 V vs NHE), whereas the two peaks in the 4-amino-2,6-dinitrotoluene voltammograms appear at −0.537 and −0.623 V and deviate more significantly from the corresponding two peaks in 2,4,6-trinitrotoluene. It suggests that the first NO₂ group reduced in the overall process is the one in <i>ortho</i> position with respect to the CH₃ group. Analogously, the 2,6-diamino-4-nitrotoluene exhibits a reduction peak at −0.629 V, almost identical to the third and second reduction peaks of 2,4,6-trinitrotoluene and 2-amino-4,6-dinitrotoluene, respectively. Since the other isomer, 2,4-diamino-6-nitrotoluene, exhibits a reduction peak at −0.712 V, we conclude that the second reduction occurs also in the <i>ortho</i> position with respect to the methyl group. Most of these observations are corroborated by quantum chemical calculations, which yielded reduction potentials in a good agreement with the experimental values (in relative scale). Thus, studying in detail all of the possible protonation and redox states in the reduction of the first nitro group and the key steps in the reduction of the second and third nitro groups, we have obtained a comprehensive and detailed picture of the mechanism of the full 18<i>e</i>^–/18H⁺ reduction of TNT. Last but not least, the calculations have shown that the thermodynamic stabilities of (isomeric) neutral radical species (<b>X </b>+<b> </b><i><b>e</b></i>^<b>–</b><b> </b>+<b> H</b>^<b>+</b>)presumably the regioselectivity-determining steps in the 6<i>e</i>^–/6H⁺ reductions of the individual NO₂ groupsare within 2 kJ·mol^–1 (i.e., comparable to RT). Therefore, the course of the reduction can be governed by the effect of the surroundings, such as the enzymatic environment, and a different regioselectivity can be observed under biological conditions

Macrocycle Conformational Sampling by DFT-D3/COSMO-RS Methodology

To find and calibrate a robust and reliable computational protocol for mapping conformational space of medium-sized molecules, exhaustive conformational sampling has been carried out for a series of seven <i>macrocyclic</i> compounds of varying ring size and one acyclic analogue. While five of them were taken from the MD/LLMOD/force field study by Shelley and co-workers (Watts, K. S.; Dalal, P.; Tebben, A. J.; Cheney, D. L.; Shelley, J. C. Macrocycle Conformational Sampling with MacroModel. J. Chem. Inf. Model. 2014, 54, 2680−2696), three represent potential macrocyclic inhibitors of human cyclophilin A. The free energy values (<i>G</i>_{DFT/COSMO‑RS}) for all of the conformers of each compound were obtained by a composite protocol based on <i>in vacuo</i> quantum mechanics (DFT-D3 method in a large basis set), standard gas-phase thermodynamics, and the COSMO-RS solvation model. The <i>G</i>_{DFT/COSMO‑RS} values were used as the reference for evaluating the performance of conformational sampling algorithms: standard and extended MD/LLMOD search (simulated-annealing molecular dynamics with low-lying eigenvector following algorithms, employing the OPLS2005 force field plus GBSA solvation) available in MacroModel and plain molecular dynamics (MD) sampling at high temperature (1000 K) using the semiempirical quantum mechanics (SQM) potential SQM­(PM6-D3H4/COSMO) followed by energy minimization of the snapshots. It has been shown that the former protocol (MD/LLMOD) may provide a more complete set of initial structures that ultimately leads to the identification of a greater number of low-energy conformers (as assessed by <i>G</i>_{DFT/COSMO‑RS}) than the latter (i.e., plain SQM MD). The CPU time needed to fully evaluate one medium-sized <i>compound</i> (∼100 atoms, typically resulting in several hundred or a few thousand conformers generated and treated quantum-mechanically) is approximately 1,000–100,000 CPU hours on today’s computers, which transforms to 1–7 days on a small-sized computer cluster with a few hundred CPUs. Finally, our data sets based on the rigorous quantum-chemical approach allow us to formulate a composite conformational sampling protocol with multiple checkpoints and truncation of redundant structural data that offers superior insights at affordable computational cost