Exploring Hydrogen Evolution Accompanying Nitrogen Reduction on

May 24, 2019 - The Supporting Information is available free of charge on the ACS ... and activation Gibbs free energies, and optimized geometries (PDF...
1 downloads 0 Views 2MB Size
Article Cite This: Inorg. Chem. 2019, 58, 7969−7977

pubs.acs.org/IC

Exploring Hydrogen Evolution Accompanying Nitrogen Reduction on Biomimetic Nitrogenase Analogs: Can Fe−NxHyIntermediates Be Active Under Turnover Conditions? Zsolt Benedek,† Marcell Papp,† Julianna Oláh,*,† and Tibor Szilvási*,‡

Downloaded via GUILFORD COLG on July 25, 2019 at 04:07:01 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.



Department of Inorganic and Analytical Chemistry, Budapest University of Technology and Economics, Szent Gellért tér 4, 1111 Budapest, Hungary ‡ Department of Chemical and Biological Engineering, University of WisconsinMadison, 1415 Engineering Drive, Madison, Wisconsin 53706, United States S Supporting Information *

ABSTRACT: Nitrogen reduction reaction (N2RR) carried out on biomimetic catalytic systems is considered to be a promising alternative for the traditional Haber−Bosch ammonia synthesis. Unfortunately, the selectivity of the currently known biomimetic catalysts is poor, as they also catalyze the unproductive hydrogen evolution reaction (HER). In the present computational study, we examine the HER activity of early N2RR intermediates in EP3 (E = B, Si) ligated single-site biomimetic iron complexes by calculating and comparing the activation Gibbs free energies of HER and N2RR elementary steps. We find that, in contrast to previous suggestions, early N2RR intermediates are not likely sources of HER under turnover conditions, as the barriers of the competing N2RR steps are significantly lower. Consequently, future research should focus on preventing other potential HER mechanisms, e.g., hydride formation, rather than accelerating the consumption of early N2RR intermediates as proposed earlier to design more efficient biomimetic catalysts.



INTRODUCTION Even though more than a century has elapsed since the invention of the Haber−Bosch ammonia synthesis, this process is still vital for modern agriculture and industry. In the past few years, the annual global production of Haber−Bosch plants exceeded 140 million tons of NH3 gas,1 which accounted for roughly 1% of the world’s energy consumption.2 In spite of the extensive efforts to develop more effective and environmental friendly ammonia synthesis methods, the industrial application of alternative synthesis procedures is at present insignificant. One potential approach for nitrogen fixation is the design of biomimetic transition metal (TM) complexes.3 These compounds mimic the active site of nitrogenase enzymes,4 which catalyze the conversion of N2 to biologically usable NH3 in microorganisms under ambient conditions. In recent years, numerous “artificial nitrogenase” families were developed,5 among which the single-site Fe complexes with EP3 ligand (E = B, C, Si)6−8 (Scheme 1, left) resemble the natural enzymatic structure most closely.9 These complexes are apt to bind dinitrogen and catalyze the N2 reduction reaction (N2RR) in the presence of proton and electron sources (Scheme 1, top right).10 The catalytic cycle of EP3-ligated Fe nitrogenases proceeds under atmospheric pressure, which is a significant advantage compared to the traditional high-pressure synthesis of NH3 from N2 and H2. Nonetheless, biomimetic molecular catalysts are far from practical application in their current state, partly due to the fact that these systems are not selective to ammonia formation: a © 2019 American Chemical Society

Scheme 1. Structure (Left) and Catalytic Reactions (Right) of Biomimetic EP3 Ligated (E = B, C, Si) Single-Site Fe Complexes with Net Chemical Equation of N2RR (Top Right) and HER (Bottom Right)

significant part of protons and electrons provided by the reagents is converted to H2 in the catalytic mixture, as the Fe complexes also function as catalysts for the hydrogen evolution reaction (HER; Scheme 1, bottom right).10 Thus, one of the main challenges in developing more efficient catalysts is finding a way to increase the selectivity toward the N2RR. The rational design of future artificial nitrogenases requires the knowledge of the complete mechanism of both N2RR and HER. In contrast to the N2RR process, whose elementary steps have been thoroughly studied by experimental10−14 and Received: March 12, 2019 Published: May 24, 2019 7969

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977

Article

Inorganic Chemistry computational works,15,16 the reaction mechanism that can contribute to HER are still under debate. There are two main HER pathways proposed in the literature as summarized in Scheme 2. One possible HER source is the early intermediates of

mechanism (route I) that competes with the productive catalytic process.17 On the other hand, the gradual accumulation of iron hydride resting states in the catalytic mixture (Scheme 2, bottom left)7,10 raises the possibility of the hydride-mediated pathway (route II), as metal hydride intermediates are typically invoked in TM-catalyzed HER.23 Furthermore, taking into account that HER is also observed with natural nitrogenases and that the enzymatic H2 formation can be traced back to labile iron hydride compounds,24 this hypothesis also seems reasonable. In parallel with experimental works, computational studies have also touched upon the question of selectivity. Recently, Matson and Peters published a DFT study, which posits the early N2RR intermediates, EP3Fe−NNH, EP3Fe−NNH2+, and EP3FeN-NH2 (E = B, C, Si), as the major sources of HER.17 Their suggestion is based on the following arguments: (i) the calculated bond dissociation free energy (BDFE) of the N−H bonds present in these structures is remarkably small (17−50 kcal/mol); (ii) bimolecular HER is thermodynamically favorable for the mentioned early N2RR intermediates; (iii) computed BDFEs follow the same trend as the experimentally determined catalytic selectivity for two out of three intermediates (see Table S9 in the Supporting Information). We note however that we very recently analyzed the Gibbs free energy profile of Fe-mediated N2RR15 and our calculations revealed that the early N2RR intermediates, including EP3Fe− NNH, EP3Fe−NNH2+, and EP3Fe−N−NH2, are favorable downhill processes: the Gibbs free energy gradually decreases by 10−30 kcal/mol in each reaction step. This is in sharp contrast to the rate-determining late N2RR intermediates of the catalyst regeneration phase, which involve thermoneutral or even endergonic elementary reaction steps. Our results suggest a very short lifetime for early N2RR intermediates that may mitigate bimolecular HER under turnover conditions while long lifetime late N2RR intermediates can participate in HER by hydride formation. We highlight here that both above presented studies15,17 rely solely on thermodynamic considerations, and only the analysis of the kinetics could reveal reliably whether early N2RR intermediates are involved in HER or not. The importance of HER vs N2RR competition motivated us to examine the HER-potential of EP3Fe−NxHy intermediates by focusing on the reaction kinetics. Therefore, in the current study, we present a computational model for quantifying the kinetic preference between HER and N2RR for the N−H bonded intermediates of the most selective catalyst (BP3Fe−N2−) and the least selective analogue (SiP3Fe−N2−). We calculate the activation Gibbs free energies of ET and PT reaction steps of N2RR and bimolecular HER elementary reaction steps starting from the most probable HER sources (i.e., the most N−H labile structures). The obtained data allow us to evaluate whether bimolecular HER can be kinetically competent with N2RR and selectivity may be improved for biomimetic catalysts in the future.

Scheme 2. Reaction Mechanisms Proposed for Biomimetic Nitrogen Fixation with EP3 (E = B, C, Si) Ligated Single Site Fe Complexesa

a

[Fe] denotes for the complete EP3 ligated single site Fe complexes shown in Scheme 1. Catalytic cycle of N2RR (green frame) highlighting the elementary steps involving [Fe]−NxHy intermediates. Catalyst deactivation pathways (orange frame) resulting in the formation of catalytically inactive and thermally stable resting states, presumably via multiple [Fe]-H bonded species. Y refers to all groups coordinating to the Fe atom apart from the EP3 ligand and the hydride moiety (e.g. N2 or a hydronitrogen group, additional hydride(s), etc.). Proposed HER mechanisms (red frame) involving [Fe]−NxHy intermediates (route I), hydride side products (route II), or only proton/electron sources (route III).

the N2RR cycle (Scheme 2, top left frame) where intermediates noted as [Fe]−NxHy (x = 1−2, y = 1−5; [Fe] = EP3Fe, E = B, C, Si) of the N2RR cycle are suggested to be able to react in bimolecular HER due to their relatively weak N−H bonds (Scheme 2, right panel, route I).17 Alternatively, H2 might be released by iron hydrides (Scheme 2, route II),18,19 which are formed as intermediates of catalyst deactivation (Scheme 2, bottom left frame),7,10 probably following the thermodynamically highly favorable protonation of a sterically accessible Fe center.20 On the basis of the HER mechanisms observed with TM-based proton reduction catalysts,21,22 homolytic (route II, pathway a) and heterolytic (route II, pathway b) HER pathways can be considered. A special case of hydride-mediated HER is the experimentally detected catalyst regeneration from catalytically inactive resting states through a series of proton transfer (PT) and electron transfer (ET) processes (route II, pathway c).10 We also note that the noncatalytic HER (Scheme 2, route III) is a significant external factor that influences the overall selectivity−nevertheless, the rate of this process is independent of the properties of the molecular catalyst. Rational catalyst and ligand design should address how to suppress the prevalent HER mechanism(s); thus, it is a fundamental question whether route I or route II is more dominant (Scheme 2). Previous experimental observations can, however, support both scenarios. On the one hand, the rapid loss of H2 from EP3Fe−NNH intermediates, formed upon the stoichiometric protonation of the catalyst anions (EP3Fe− N2−),6,10 is considered as evidence for a bimolecular HER



COMPUTATIONAL DETAILS All computations were carried out using the ORCA 4.0.1.2. program.25 See the Supporting Information for sample input files for the calculations described below. As a general consideration, the properties of all studied structures were computed in the most stable spin state, as determined by experiments10,24,26 andin the absence of experimental data our previously published theoretical analysis.15 Geometry Optimizations. Similarly to our recent DFT study where we extensively tested our computational method7970

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977

Article

Inorganic Chemistry ology,15 we use BP86 density functional with resolution of identity (RI) approximation27,28 with ZORA-def2-SVP basis set29,30 and def2/J auxiliary basis set31 for the optimizations. Aiming to obtain reasonable geometries for the studied transition metal complexes, we introduced zeroth order relativistic approximation (ZORA)32 and Grimme-type dispersion correction with the Becke−Johnson damping function (D3(BJ))33,34 to the calculations. Local minima found on the potential energy hypersurface (PES) were confirmed by harmonic vibrational frequency calculations at the same level of theory (no imaginary mode was observed in the vibrational analysis of initial states). Single-Point Energy Calculations. To calculate reliable energies, we benchmarked the chosen functional and basis set using ab initio reference calculations. As canonical coupled cluster (CC) methods are computationally unaffordable for the bulky complexes studied in this work, a compromise needed to be found between accuracy and computational cost. After deliberation, we chose the recently developed DLPNO− CCSD(T)35,36 method for computing reference energies (see Supporting Information for the detailed description of these calculations), as it has been successfully applied to study the reactivity of transition metal complexes akin to single-site Fe nitrogenases.37,38 We carried out single-point calculations on the optimized geometries of early intermediates of nitrogen reduction. The selected DFT protocol (B3LYP density functional39−42 including relativistic (ZORA) and dispersion corrections (D3(BJ)), used with ZORA-def2-TZVP basis set29,30 on all atoms) was carefully validated against available experimental data in our previous study.15 Additionally, we found that our DFT protocol provides accurate energy values compared to DLPNO−CCSD(T). See the Supporting Information for comparison of ab initio and DFT results. In order to mimic the experimental conditions, DFT single-point calculations were performed in diethyl ether (DEE) solution, using the solvation model based on density (SMD).43 The RIJCOSX density fitting approximation44 with def2/J auxiliary basis set was applied to reduce computational time. Calculation of Activation Barriers. Activation Gibbs free energies of ET and PT depend on the identity of acid and reductant reagents present in the catalytic mixture. In the present study, we calculated the barriers with [(Et2O)2H]+ as proton source and KC8 as electron source. In this way, the catalytic conditions of ref 10 are modeled. Transition State Theory. In the case of PT and bimolecular HER, we located the first order saddle point on the PES that corresponds to the transition state (TS) of the reaction. The obtained TS structures were confirmed by the presence of an imaginary frequency in the calculated vibrational spectrum, belonging to a proton migration (in the case of PT) or the formation/dissociation of an H2 molecule (in the case of HER), respectively. Handling Spin State Changes during PT and HER. In the case of PT and bimolecular HER steps, the ground spin state of initial states might differ from that of the product. In the reactions studied herein both spin-allowed and spin-forbidden routes are conceivable.45 The former case can be properly addressed by traditional transition state theory, but the precise treatment of spin-forbidden transitions requires the location of the minimum energy crossing points (MECPs) between the potential energy surfaces (PESs) of the involved spin states. These MECPs serve as TSs for spin-forbidden reactions and

nonadiabatic transition state theory can be used to obtain effective activation energies for the reaction. The predicted barriers will strongly depend on the spin−orbit coupling term (SOC) of the system’s Hamiltonian.46 As the calculation of the MECPs and the SOC matrix elements requires an enormous computational effort for the large systems studied in this work, we decided to approximate barriers by determining the transition states in all conceivable spin states and the Gibbs free energy of the activation barrier was taken from the lowest energy TS (see Tables S12 and S18 in the Supporting Information). Marcus Theory. The activation barrier of ET reactions was calculated using Marcus theory,47 according to which the activation Gibbs free energy of an ET process (ΔG‡) can be written as ΔG‡ =

(λ + ΔG)2 4λ

where ΔG is the total Gibbs free energy change and λ denotes the total reorganization energy. The latter can be expressed as the sum of inner-sphere and outer-sphere reorganization energies: λ = λis + λos

Given that KC8 shows electrode-like behavior toward the electron acceptor iron complexes, the inner-sphere contribution is48,49 λis = 0.5[Eox (R red) − Eox (R ox ) + Ered(R ox ) − Ered(R red)]

where Eox and Ered denote the gas-phase energy of the iron complex to be reduced by KC8, computed in the oxidized and reduced state, respectively. (As DEE solvation is handled here as an outer-sphere effect, the SMD energy contributions were removed from these calculations, see details in the Supporting Information) Rox and Rred stand for reduced and oxidized equilibrium geometries, respectively, for which the single-point energies are to be calculated. The outer-sphere reorganization energy can be approximated as49,50 2 2 ij 1 ij 1 1 yz Δq 1 yz Δq − zzz − jjj − zzz λos = jjj j ε∞ j ε∞ ε0 z{ 2r ε0 z{ 4d k k

where Δq is the charge transferred during the reduction, ε∞ and ε0 refer to the optical and static dielectric constants of the solvent, respectively, r is the radius of the cavity formed by the electron acceptor transition metal complex in the solvent (which is handled as spherical in this model), and d is the distance between the center of the complex and the surface of the reductant. Substituting the characteristic constants of DEE (ε∞ = 1.99, ε0 = 5.56 relative to the permittivity of vacuum, estimated at the temperature of the catalytic process (195 K)51) and Δq = 1e to the equation above and assuming that the iron complex moieties touch the outer graphene layer of KC8 (d = r), the term for λos gives ÄÅ É ÅÅ kcal ÑÑÑ Å ÑÑ = 336 λos ÅÅ ÅÅÇ mol ÑÑÑÖ r [Å] The r parameters were calculated from the solvent excluded surface (SES) based molecular cavity volume, as provided by the GEPOL algorithm52 implemented in ORCA. 7971

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977

Article

Inorganic Chemistry

Figure 1. Computationally determined N−H bond dissociation free energy (BDFEN−H in kcal/mol) of BP3Fe−NxHy intermediates (left) and SiP3Fe−NxHy intermediates (right), formed in an N2RR cycle. HER sources proposed earlier17 are marked with blue lines and fonts. The borders between the ranges of HER-inactive structures, potential HER sources, and possible HER sources were drawn based on the experimentally determined BDFEN−H of the proven HER-inactive intermediate SiP3Fe−CNH+ (62 kcal/mol) and that of the proven HER-active SiP3FeCNH (41 kcal/ mol) under noncatalytic conditions.

Thermochemistry. As a general protocol, Gibbs free energies were computed at 195 K and 1 atm, in accordance with the experimental conditions reported for the studied reactions.10 In the case of BDFE calculations, however, the temperature was set to 298 K as the reference free energies available in the literature26 were determined at room temperature. In order to simplify the comparison of the thermodynamics and kinetics of the different reaction types that can occur in the modeled catalytic mixture, the energy changes (ΔG, Gibbs free energy of reaction; ΔG‡, activation Gibbs free energy) calculated according to (1) and (2) from the computed individual Gibbs free energies (G) of the starting materials, products, and transition states were postprocessed. The standard-state correction formulas (3) and (4) were applied to obtain the Gibbs free energies relevant for a 1 M solution;53 in these equations, the molar concentrations of the starting materials and products are given in parentheses, c = 0.0625 M stands for the gas-phase molar concentration (c = p/RT) at experimental temperature (T) and pressure (p), and Q denotes the reaction quotient. On the basis of (3) and (4), we added 1.08 kcal/mol to the raw Gibbs free energy if the number of moles increases (bimolecular HER, protonation from [Et2O)2H]+ during the reaction) and subtracted the same value if the number of moles decreases (formation of the transition states of HER and PT). ΔG =

∑ Gproducts − ∑ Gstarting materials

ΔG‡ =

∑ Gtransition state − ∑ Gstarting materials

correction was applied here as the electron transfer does not change the number of molecules, as a result of which ΔG equals to ΔG(c = 1 M) in eq 3. Then, the resulting ΔG value was adjusted according to the lower reduction potential (E0) of KC8 (E0 ≤ −3.0 V compared to the ferrocene electrode) relative to Cp*2Co (E0 = −1.96 V compared to the ferrocene electrode):13 a correction of F(E0KC8,max. − E0Cp*2Co) = −24.00 kcal/mol was added to ΔG, where F stands for the Faraday constant. Calculation of Bond Dissociation Free Energies (BDFEs). BDFEs of the N−H bonds present in EP3Fe−NxHy intermediates were calculated based on the computed Gibbs free energy (ΔG) of the following hypothetic hydrogen atom transfer (HAT) reaction: SiP3Fe −CN+ + EP3Fe −NxHy → SiP3Fe −CNH+ + EP3Fe −NxHy − 1

The complex SiP3Fe−CNH+ was used here as a reference because its BDFEN−H (62 kcal/mol) has been determined experimentally.26 On the basis of this value, the bond dissociation free energy of a given EP3Fe−NxHy complex can be calculated as BDFE N − H(EP3Fe− NxHy) = ΔG(HAT) + 62 kcal/mol

The experimental BDFE of SiP3Fe−CNH+ was measured in THF solution at 298 K and 1 atm. Thus, we also performed the single point calculations taking these conditions (instead of 195 K and DEE solvent as in other cases) into account. The THF environment was simulated using the CPCM solvation model.54

(1) (2)



ΔG(c = 1 M) = ΔG + RT ln(Q (c = 1 M)/Q (c = 0.0625 M))

RESULTS AND DISCUSSION We elaborated a two-step methodology to explore the HER activity of all N2RR intermediates with reasonable computational cost. In the first step, we prescreen the intermediates based on N−H bond energy, which should be cleaved during HER, and identify the most probable sources of bimolecular HER with small N−H bond energy. Then, the N2RR and potential bimolecular HER reactivity of the selected intermediates are examined in detail by comparing the activation Gibbs free energies of the competing elementary steps. In case the barrier of at least one of the conceivable bimolecular HER

(3)

ΔG‡(c = 1 M) = ΔG‡ + RT ln(Q (c = 1 M)/Q (c = 0.0625 M))

(4)

In the case of ET processes, the Gibbs free energy was determined indirectly using decamethylcobaltocene (Cp*2Co) as standard because it is difficult to model potassium graphite by DFT directly. In the first step, the Gibbs free energy of reduction by Cp*2Co was calculated at 195 K. No concentration 7972

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977

Article

Inorganic Chemistry

be generalized: in the case of EP3FeNH+, EP3FeNH2+, and EP3FeNH−NH2+, E = Si gives higher BDFEN−Hs, which confirms that no systematic effect of E can be determined. It has not escaped our attention that our data point to different conclusions than the results presented in ref 17, where the general stabilizing effect of the boron atom and surprisingly large differences between BDFEN−Hs (13 kcal/mol in average; consult with Table S9 of the Supporting Information for details) were reported.17 We think that these differences stem from the different density functionals used in the two studies. Our benchmark calculations compared to coupled cluster results (see details in the computational methods section and the Supporting Information) revealed consistently small errors for our DFT method while inconsistent and large errors for TPSS functional used in ref 17 (see details in the Supporting Information, section 1.2). Following the prescreening process, we continued our examinations handling EP3Fe−NNH (E = B, Si) as the only potential source of H2, andfor completenessthe intermediates EP3Fe−NNH2+ (E = B, Si), EP3FeN-NH2 (E = B, Si), SiP3Fe-NH3, and SiP3Fe-NH2−NH2 as conceivable additional HER sources. We note that BP3Fe−NNH2+ was characterized here as a HER-inactive compound, but we studied its reactivity further because an earlier study suggested it to be a bimolecular HER source,17 and unequivocal classification of such borderline cases is not possible as DFT errors are inevitably present. Comparison of Activation Energies of HER and N2RR Steps. The N−H bond dissociation free energy alone is not satisfactory to characterize bimolecular HER under turnover conditions, as it is possible that the concurrent N2RR step occurs significantly faster because of lower activation barrier. Thus, the previously identified possible and ambiguous HER sources can only be considered as the actual sources of H2 side product observed experimentally if the kinetic competence of the hydrogen evolution can be demonstrated. In order to gain a quantitative picture of the relation of the reaction rates of N2RR and HER, our next goal was to calculate the activation Gibbs free energies. We began our investigations with the earliest N2RR intermediates (EP3Fe−NNH, EP3Fe−NNH2+, EP3FeN−NH2) which are of special interest as they are currently viewed as the most probable source intermediates of HER.17 We aimed to compare the activation barriers of the conceivable bimolecular HER steps to those obtained for the proton additions from [(Et2O)2H]+ to neutral structures (EP3Fe−NNH + H+ → EP3Fe−NNH2+, EP3FeN−NH2 + H+ → EP3FeN−NH3+ or EP3FeN− NH2 + H+ → EP3FeNH−NH2+) and for the reduction of the positively charged EP3Fe−NNH2+ by KC8. These choices are motivated by recent computational15 and experimental10 results showing that the formation of negatively charged EP3Fe−NxHy intermediates during N2RR, as well as the proximal protonation of EP3Fe−NNH, is incompetent under catalytic conditions. The search for the transition state of HER steps is especially challenging due to the size and complexity of the bimolecular activated complex. Furthermore, in contrast to PT steps from [(Et2O)2H]+ to early N2RR intermediates, in the case of which certain transition structures have already been identified with E = B,16 TS calculations for the bimolecular HER process are unprecedented. At this stage of our research, we realized that several simplifications need to be introduced to our theoretical model (described in the Computational Details) in order to complete all calculations within reasonable computational time.

reactions is calculated to be lower than that for the competing N2RR step, it can be regarded as a strong evidence for the presence of HER route I (see Scheme 2) under turnover conditions. In the opposite case, however, the hypothesis of facile H-transfer between N2RR intermediates becomes highly questionable. Identifying the Potential HER Sources Based on BDFEN−H. As the initial step of our investigations, we calculated the N−H bond dissociation free energies of all EP3Fe−NxHy intermediates suspected to be formed during biomimetic nitrogen fixation. On the basis of the reasonable assumption that low BDFEN−H indicates high affinity to bimolecular HER (and vice versa),55 this procedure enabled us to classify the intermediates according to their potential HER activity and screen out the implausible sources of HER. Nevertheless, it is not evident how to choose the cutoff value of BDFEN−H, above which an EP3Fe−NxHy species can be considered to be stable against HER. In order to interpret the obtained BDFEs in an experimentally relevant way, we compared the N−H bond energy of the studied EP3Fe−NxHy intermediates to that of two reference structures: SiP 3 Fe−CNH + (experimental BDFEN−H: 62 kcal/mol) and SiP3FeCNH (experimental BDFEN−H: 41 kcal/mol). The former iron complex is known to be permanently stable at 195 K, while the latter rapidly decomposes at the same temperature to SiP3Fe−CN and H2 due to bimolecular HER.26 In accordance with these reference data, we divide the studied intermediates into three categories, as shown in Figure 1: (i) HER-inactive structures (BDFEN−H ≥ 62 kcal/mol; green range), which cannot participate in bimolecular HER; (ii) possible HER sources (BDFEN−H ≤ 41 kcal/mol; red range), which may participate in bimolecular HER under catalytic conditions; (iii) ambiguous HER sources (62 kcal/mol > BDFEN−H > 41 kcal/mol; yellow range), the HER activity of which cannot be excluded based on the calculated BDFEN−Hs. According to Figure 1 (see numerical data in the Supporting Information), we find that the N−H bonds in the intermediates of N2RR are in general far from labile. Our results indicate that the vast majority of N−H bonded structures is clearly HER inactive and only one N2RR intermediate, EP3Fe−NNH (the bimolecular HER activity of which is experimentally confirmed) falls into the category of possible HER sources (BDFEN−H = 29 kcal/mol with both E = B and E = Si). Additionally, we identify a few ambiguous HER sources (BP3FeN−NH2, SiP3FeN− NH2, SiP3Fe−NNH2+, SiP3Fe−NH3 and SiP3Fe−NH2− NH2), however, the BDFEN−Hs calculated for these structures are all higher than 52 kcal/mol thus located in the upper half of the corresponding range of ambiguous HER sources (41−62 kcal/mol). Strikingly, our results do not confirm previous computational suggestions17 that EP3Fe−NNH2+/0 intermediates are more likely to release H2 than EP3Fe−NNH; what is more, BP3Fe−NNH2+ is characterized as a HERinactive intermediate in the light of our data. Our results also show that the effect of the axial ligand atom (E) on BDFEN−H is negligible. In Figure 1, it can be observed that the analogous structures, which differ only in the atom E, are located close to each other: the average energy difference between BDFEN−Hs is as low as 4 kcal/mol, and even the highest deviation remains below 9 kcal/mol. Thus, it can be concluded that the distance of two or three chemical bonds between E and the N−H bond to be cleaved is long enough to minimize axial ligand effects. Additionally, while the N−H bond is slightly more stable with E = B in the majority of the cases, this relation cannot 7973

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977

Article

Inorganic Chemistry

While the activation Gibbs free energy of symmetrical HER steps was calculated directly from the energy level of the TS structures above, the barrier for HER steps involving two different intermediates were approximated. We applied the Brønsted−Evans−Polányi (BEP) principle,15,58,59 according to which the barrier (ΔG‡) can be expressed as a linear function of the Gibbs free energy change (ΔG): ΔG‡ = αΔG + β. The BEP parameters (α, β) were determined from ΔG−ΔG‡ data pairs obtained from symmetrical reactions (see Figure 3).

First, the activation barriers for bimolecular HER reactions were determined with methyl groups in lieu of the isopropyl substituents of the full ligand scaffold. We note that this simplification causes the underestimation of bimolecular HER activation barriers because of the partial neglect of steric repulsion effects. Second, transition states of HER were calculated only for the symmetrical cases (i.e., reaction between two identical iron complexes), while the barriers for asymmetrical HER steps (i.e., reaction between two different intermediates) were estimated by assuming the Brønsted− Evans−Polányi (BEP) relation, a linear correlation between Gibbs free energies of reaction (ΔG), and activation Gibbs free energies (ΔG‡).56,57 BEP relations are widely used to reliably estimate activation energies for structurally similar activation barriers.58,59 Additionally, we show that barriers obtained from BEP relations are generally very high; thus, they cannot play any role in bimolecular HER even considering very large 10 kcal/ mol error bars. The obtained geometries for the transition states, shown together with the initial and final complexes for comparison, are visualized in Figure 2 in the case of E = Si. An analogous figure for E = B, showing highly similar distances and angles, can be found in the Supporting Information as Figure S1.

Figure 3. Brønsted−Evans−Polányi (BEP) parameters determined for bimolecular HER, based on ΔG−ΔG‡ data pairs obtained for symmetrical HER reactions (see Table 1 for numerical data).

The resulting activation barriers for the HER, PT and ET processes specified above are summarized in Table 1. In accordance with the high stability of N−H bonds in EP3Fe− NNH2+/0 intermediates (BDFEN−H= 52−63 kcal/mol), huge differences can be observed between HER and N2RR activation barriers. The reduction of EP3Fe−NNH2+ by KC8 proceeds with a very small activation barriers (ΔG‡ = 4−6 kcal/mol), while the activation Gibbs free energy of HER between two such intermediates exceeds 50 kcal/mol and even the smallest barrier is 32 kcal/mol. The situation is similar with the neutral EP3Fe N−NH2: the obtained smallest barrier for HER is 29 kcal/mol for the reaction of SiP3FeN−NH2 with SiP3Fe−NNH, but the PT barrier is only 19 kcal/mol which is kinetically significantly more favorable. We emphasize that the ΔG‡ of HER is strictly an underestimation as the iPr groups were reduced to methyl groups, which decreased the steric repulsion. Thus, we think that the 10 kcal/mol difference between the most favorable HER and PT barrier cannot be attributed to the error of the computational model or DFT method alone. We also note that our conclusion that EP3Fe−NNH2+/0 intermediates are stable against HER under turnover conditions is also supported by some previous experiments as it has been shown that SiP3Fe−NNH2+ is “persistent for hours in solution” at the temperature of the catalytic reaction.12 As expected from BDFEN−H calculations, we find the activation Gibbs free energies of bimolecular HER to be the lowest if EP3Fe−NNH species are involved. The symmetrical HER reaction between two EP3Fe−NNH complexes proceeds with a barrier of 12 and 18 kcal/mol for SiP3Fe− NNH and BP3Fe−NNH, respectively. These small numbers are in agreement with the experimentally observed spontaneous HER from in situ formed EP3Fe−NNH.6,10,17 Nevertheless, the concurrent PT reaction can still occur under turnover conditions. In the case of BP3Fe−NNH, the kinetic hindrance of HER is clearly indicated by the 11 kcal/mol difference in the activation barriers (ΔG‡(PT) = 7 kcal/mol, ΔG‡(HER) = 18 kcal/mol). As for SiP3Fe−NNH, the activation Gibbs free energy of concurrent PT and HER steps seems to be very close to each other at first sight (ΔG‡(PT) = 13

Figure 2. Structures of initial states (a), transition states (b), and final states (c) of relevant HER reactions starting from SiP3Fe−NNH (left column) SiP3Fe−NNH2+ (middle column), and SiP3FeN− NH2 (right column). As the transition states of HER were localized with methyl (instead of isopropyl) substituted phosphine ligands, we also show methyl substituted complexes in rows a and c for consistency. Hydrogen atoms (except for those adjacent to N) are omitted for clarity. Color code: yellow, Fe; blue, N; white, H; cyan, Si; orange, P. The central atoms of the complexes are shown in ball and stick, for which selected geometrical parameters (bond distances in ångstroms, angles in degrees) are provided.

Considering the bond distances and angles provided for the presented structures, in can be observed that except for the simultaneous elongation of the N−H bonds to be cleaved (from 1.03 to 1.04 to 1.23−1.42 Å), no significant distortion of the central part of initial complexes occur while they approach each other (Fe−N bond lengths, N−N bond lengths, and N−N−H angles remain almost intact). In the TS, where the reacting H atoms are only separated by a distance of 1.02−1.11 Å, formation of the H−H bond takes place concurrently with the eventual breakage of both N−H bonds. 7974

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977

Article

Inorganic Chemistry

Table 1. Gibbs Free Energy (ΔG) and Activation Gibbs Free Energy (ΔG‡) of Concurrent HER and N2RR Elementary Steps Starting from Early EP3Fe−NxHy Intermediates (E = B, Si)c reaction types

BP3Fe−NNH

BP3Fe−NNH2+

BP3FeN−NH2

BDFEN−H = 29 kcal/mol

BDFEN−H = 63 kcal/mol

BDFEN−H = 56 kcal/mol

ΔG (kcal/mol) bimolecular HER with

N2RR

BP3Fe−NNH BP3Fe−NNH2+ BP3FeN−NH2 PT to distal N PT to proximal N ET

reaction types

−65 18 −36 34a −40 31a −16 7 − − − − SiP3Fe−NNH

N2RR

SiP3Fe−NNH SiP3Fe−NNH2+ SiP3FeN−NH2 PT to distal N PT to proximal N ET

−67 −39 −43 −21 − −

ΔG (kcal/mol)

ΔG‡ (kcal/mol)

ΔG (kcal/mol)

−36 34 −7 59 −11 50a − − − − −37 4b SiP3Fe−NNH2+ BDFEN−H = 60 kcal/mol

ΔG‡ (kcal/mol)

ΔG (kcal/mol)

12 32a 29a 13 − −

−39 −11 −15 − − −30

ΔG‡ (kcal/mol)

−40 31a −11 50a −15 42 −7 13 −27 14 − − SiP3FeN−NH2

a

BDEFN−H = 29 kcal/mol ΔG (kcal/mol)

Bimolecular HER with

ΔG‡ (kcal/mol)

ΔG‡ (kcal/mol)

BDFEN−H = 52 kcal/mol ΔG (kcal/mol)

ΔG‡ (kcal/mol)

−43 −15 −20 −9 −36 −

29a 47a 40 19 20 −

a

32 53 47a − − 6b

a

Calculated from the Brønsted−Evans−Polányi relation. bCalculated from Marcus theory. cPT: proton transfer. ET: electron transfer.

kcal/mol, ΔG‡(HER) = 12 kcal/mol). However, as we already mentioned earlier, ΔG‡(HER) is strictly an underestimation because of the smaller steric groups used for the calculations. Because this TS is the only one where we could not determine the preference for N2RR or for bimolecular HER, we calculated the bimolecular HER step using iPr groups on SiP3Fe−NNH, which was also an important check whether iPr groups indeed induce steric exclusion or not. Our calculated barrier with iPr groups is 30 kcal/mol, which is significantly larger than that we obtained with methyl groups (12 kcal/mol, Table 1). We note that, for such a large system, finding the most favorable conformation for the lowest energy bimolecular TS is a daunting problem and computationally it is not feasible to fully tackle it nowadays. We built on our previous experience finding TSs with methyl groups and chemical intuition to obtain the lowest energy TS that we could. Nevertheless, because of the large 18 kcal/mol difference between the iPr and methyl-substituted bimolecular HER TSs, we think that our original assumption holds, that is, iPr groups indeed lead to increased steric hindrance compared to methyl groups even if there might be somewhat lower energy TSs with more favorable conformations for the iPr-substituted system. Additionally, because of the 17 kcal/mol difference between the iPr-substituted bimolecular HER TS and the concurrent PT step, we think that the PT step with a barrier of only 12 kcal/mol is clearly favored even if we consider DFT error or lower energy TS of bimolecular HER. For the rest of the ambiguous HER source intermediates (SiP3Fe−NH3 and SiP3Fe−NH2−NH2), we also found that the barriers of the corresponding elementary steps of N2RR are sufficiently small to avoid the undesired hydrogen evolution (Table 2). As has been pointed out earlier, the smallest activation Gibbs free energy of bimolecular HER is obtained when the reaction partner is SiP3Fe−NNH with the most labile N−H bond. In this case, the change in Gibbs free energy during HER is −37 kcal/mol with both SiP3Fe−NH2−NH2 and SiP3Fe−NH3), which corresponds to barriers of 33 kcal/mol according to the BEP relationship (Figure 3). In contrast, the activation Gibbs free energy of ammonia release from SiP3Fe−

Table 2. Activation Gibbs Free Energy (ΔG‡) of Concurrent HER and N2RR Elementary Steps Starting from SiP3Fe− NH2−NH2 and SiP3Fe−NH3 Intermediatesb reaction types Bimolecular HER with SiP3Fe−NNH N2RR PT NH3↑

SiP3Fe−NH2−NH2

SiP3Fe−NH3

33a

33a

19 −

− 14

a

Calculated from the Brønsted−Evans−Polányi relation. bPT: proton transfer. NH3↑: ammonia release.

NH3 and that of PT to SiP3Fe−NH2−NH2 are 14 and 19 kcal/ mol, respectively.



CONCLUSIONS The present theoretical study aimed to determine whether the bimolecular HER mechanism (2 EP3Fe−NxHy → 2 EP3Fe− NxHy−1 + H2) can be competitive under catalytic turnover conditions, which has been previously suggested to be an important selectivity-limiting factor of biomimetic single-site Fe nitrogenases. According to the results of our calculations, all possible HER steps starting from EP3Fe−NxHy intermediates are kinetically hindered in a typical catalytic mixtureeither due to the high stability (high BDFE) of the N−H bonds to be cleaved for H2 evolution or due to the significantly lower activation barrier of the competing productive N2RR steps. Consequently, we propose the H2 side product observed experimentally originates from other sources, presumably from iron hydrides (as observed with proton reduction catalysts21−23). Therefore, the recently suggested strategy of accelerating the initial N2RR steps is not expected to improve catalytic selectivity. We recommend further experimental and theoretical studies on the mechanism of hydride formation, as this knowledge can be crucial for increasing the selectivity of future catalysts. Finally, we note that the exact mechanism of hydride -mediated H2 release is still unknown and should be the subject of future research. 7975

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977

Article

Inorganic Chemistry



(11) Anderson, J. S.; Cutsail, G. E.; Rittle, J.; Connor, B. A.; Gunderson, W. A.; Zhang, L.; Hoffman, B. M.; Peters, J. C. Characterization of an Fe≡N-NH2Intermediate Relevant to Catalytic N2Reduction to NH3. J. Am. Chem. Soc. 2015, 137 (24), 7803−7809. (12) Rittle, J.; Peters, J. C. An Fe-N2Complex That Generates Hydrazine and Ammonia via FeNNH2: Demonstrating a Hybrid Distal-to-Alternating Pathway for N2Reduction. J. Am. Chem. Soc. 2016, 138 (12), 4243−4248. (13) Chalkley, M. J.; Del Castillo, T. J.; Matson, B. D.; Roddy, J. P.; Peters, J. C. Catalytic N2-to-NH3Conversion by Fe at Lower Driving Force: A Proposed Role for Metallocene-Mediated PCET. ACS Cent. Sci. 2017, 3 (3), 217−223. (14) Thompson, N. B.; Green, M. T.; Peters, J. C. Nitrogen Fixation via a Terminal Fe(IV) Nitride. J. Am. Chem. Soc. 2017, 139 (43), 15312−15315. (15) Benedek, Z.; Papp, M.; Oláh, J.; Szilvási, T. Identifying the RateLimiting Elementary Steps of Nitrogen Fixation with Single-Site Fe Model Complexes. Inorg. Chem. 2018, 57 (14), 8499−8508. (16) Kaczmarek, M. A.; Malhotra, A.; Balan, G. A.; Timmins, A.; de Visser, S. P. Nitrogen Reduction to Ammonia on a Biomimetic Mononuclear Iron Centre: Insights into the Nitrogenase Enzyme. Chem. - Eur. J. 2018, 24 (20), 5293−5302. (17) Matson, B. D.; Peters, J. C. Fe-Mediated HER vs N2RR: Exploring Factors That Contribute to Selectivity in P3EFe(N2) (E = B, Si, C) Catalyst Model Systems. ACS Catal. 2018, 8 (2), 1448−1455. (18) Deegan, M. M.; Peters, J. C. Electrophile-Promoted Fe-to-N2 hydride Migration in Highly Reduced Fe(N2)(H) Complexes. Chem. Sci. 2018, 9 (29), 6264−6270. (19) Shima, T.; Hou, Z. Dinitrogen Fixation by Transition Metal Hydride Complexes. Top. Organomet. Chem. 2017, 60, 23−43. (20) Creutz, S. E.; Peters, J. C. Exploring Secondary-Sphere Interactions in Fe-NxHy complexes Relevant to N2 fixation. Chem. Sci. 2017, 8 (3), 2321−2328. (21) Dempsey, J. L.; Brunschwig, B. S.; Winkler, J. R.; Gray, H. B. Hydrogen Evolution Catalyzed by Cobaloximes. Acc. Chem. Res. 2009, 42 (12), 1995−2004. (22) Besora, M.; Lledós, A.; Maseras, F. Protonation of TransitionMetal Hydrides: A Not so Simple Process. Chem. Soc. Rev. 2009, 38 (4), 957−966. (23) Bullock, R. M.; Appel, A. M.; Helm, M. L. Production of Hydrogen by Electrocatalysis: Making the H-H Bond by Combining Protons and Hydrides. Chem. Commun. 2014, 50 (24), 3125−3143. (24) Khadka, N.; Milton, R. D.; Shaw, S.; Lukoyanov, D.; Dean, D. R.; Minteer, S. D.; Raugei, S.; Hoffman, B. M.; Seefeldt, L. C. Mechanism of Nitrogenase H2Formation by Metal-Hydride Protonation Probed by Mediated Electrocatalysis and H/D Isotope Effects. J. Am. Chem. Soc. 2017, 139 (38), 13518−13524. (25) Neese, F. Software Update: The ORCA Program System, Version 4.0. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2018, 8 (1), e1327. (26) Rittle, J.; Peters, J. C. N-H Bond Dissociation Enthalpies and Facile H Atom Transfers for Early Intermediates of Fe-N2 and Fe-CN Reductions. J. Am. Chem. Soc. 2017, 139 (8), 3161−3170. (27) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38 (6), 3098−3100. (28) Perdew, J. P. Density-Functional Approximation for the Correlation Energy of the Inhomogeneous Electron Gas. Phys. Rev. B: Condens. Matter Mater. Phys. 1986, 33 (12), 8822−8824. (29) Weigend, F.; Ahlrichs, R. Balanced Basis Sets of Split Valence, Triple Zeta Valence and Quadruple Zeta Valence Quality for H to Rn: Design and Assessment of Accuracy. Phys. Chem. Chem. Phys. 2005, 7 (18), 3297−3305. (30) Pantazis, D. A.; Chen, X.-Y.; Landis, C. R.; Neese, F. All-Electron Scalar Relativistic Basis Sets for Third-Row Transition Metal Atoms. J. Chem. Theory Comput. 2008, 4 (6), 908−919. (31) Weigend, F. Accurate Coulomb-Fitting Basis Sets for H to Rn. Phys. Chem. Chem. Phys. 2006, 8 (9), 1057−1065. (32) van Wüllen, C. Molecular Density Functional Calculations in the Regular Relativistic Approximation: Method, Application to Coinage

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.inorgchem.9b00719.



Sample ORCA input files, ab initio benchmark calculations, details of calculation of Gibbs free energies and activation Gibbs free energies, and optimized geometries (PDF)

AUTHOR INFORMATION

Corresponding Authors

*(J.O.) E-mail: [email protected]. *(T.S.) E-mail: [email protected]. ORCID

Zsolt Benedek: 0000-0002-0174-9676 Marcell Papp: 0000-0002-2464-4818 Julianna Oláh: 0000-0003-4608-3115 Tibor Szilvási: 0000-0002-4218-1570 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We acknowledge NIIF for awarding us access to supercomputing resources based in Hungary at Debrecen and Szeged. Z.B. and M.P. express thanks for the support of the New National Excellence Program of the Ministry of Human Capacities of Hungary. J.O. is also thankful for the financial support of NKFIH Grant No. 115503 and for a KGYNK Fellowship.



REFERENCES

(1) Mineral Commodity Summaries; U. S. Geological Survey: Reston, VA, 2018. (2) Renner, J. N.; Greenlee, L. F.; Ayres, K. E.; Herring, A. M. Electrochemical Synthesis of Ammonia: A Low Pressure, Low Temperature Approach. Electrochem. Soc. Interface 2015, 24 (2), 51− 57. (3) Shaver, M. P.; Fryzuk, M. D. Activation of Molecular Nitrogen: Coordination, Cleavage and Functionalization of N2 Mediated By Metal Complexes. Adv. Synth. Catal. 2003, 345 (910), 1061−1076. (4) MacLeod, K. C.; Holland, P. L. Recent Developments in the Homogeneous Reduction of Dinitrogen by Molybdenum and Iron. Nat. Chem. 2013, 5 (7), 559−565. (5) Stucke, N.; Flöser, B. M.; Weyrich, T.; Tuczek, F. Nitrogen Fixation Catalyzed by Transition Metal Complexes: Recent Developments. Eur. J. Inorg. Chem. 2018, 2018 (12), 1337−1355. (6) Anderson, J. S.; Rittle, J.; Peters, J. C. Catalytic Conversion of Nitrogen to Ammonia by an Iron Model Complex. Nature 2013, 501 (7465), 84−87. (7) Creutz, S. E.; Peters, J. C. Catalytic Reduction of N2 to NH3 by an Fe-N2 complex Featuring a C-Atom Anchor. J. Am. Chem. Soc. 2014, 136 (3), 1105−1115. (8) Lee, Y.; Mankad, N. P.; Peters, J. C. Triggering N2 uptake via Redox-Induced Expulsion of Coordinated NH3 and N2 silylation at Trigonal Bipyramidal Iron. Nat. Chem. 2010, 2 (7), 558−565. (9) Broda, H.; Tuczek, F. Catalytic Ammonia Synthesis in Homogeneous Solution-Biomimetic at Last? Angew. Chem., Int. Ed. 2014, 53 (3), 632−634. (10) Del Castillo, T. J.; Thompson, N. B.; Peters, J. C. A Synthetic Single-Site Fe Nitrogenase: High Turnover, Freeze-Quench 57Fe Mössbauer Data, and a Hydride Resting State. J. Am. Chem. Soc. 2016, 138 (16), 5341−5350. 7976

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977

Article

Inorganic Chemistry

computation of a solvent-excluding surface. J. Comput. Chem. 1994, 15 (10), 1127−1138. (53) Ashcraft, R. W.; Raman, S.; Green, W. H. Ab Initio Aqueous Thermochemistry: Application to the Oxidation of Hydroxylamine in Nitric Acid Solution. J. Phys. Chem. B 2007, 111 (41), 11968−11983. (54) Barone, V.; Cossi, M. Quantum Calculation of Molecular Energies and Energy Gradients in Solution by a Conductor Solvent Model. J. Phys. Chem. A 1998, 102 (11), 1995−2001. (55) Bezdek, M. J.; Guo, S.; Chirik, P. J. Coordination-induced weakening of ammonia, water, and hydrazine X−H bonds in a molybdenum complex. Science 2016, 354 (6313), 730−733. (56) Brønsted, J. N. Acid and Basic Catalysis. Chem. Rev. 1928, 5 (3), 231−338. (57) Evans, M. G.; Polányi, M. Inertia and driving force of chemical reactions. Trans. Faraday Soc. 1938, 34, 11−24. (58) Wang, S.; Temel, B.; Shen, J.; Jones, G.; Grabow, L. C.; Studt, F.; Bligaard, T.; Abild-Pedersen, F.; Christensen, C. H.; Nørskov, J. K. Universal Brønsted-Evans-Polanyi Relations for C−C, C−O, C−N, N−O, N−N, and O−O Dissociation Reactions. Catal. Lett. 2011, 141 (3), 370−373. (59) Sutton, J. E.; Vlachos, D. G. A Theoretical and Computational Analysis of Linear Free Energy Relations for the Estimation of Activation Energies. ACS Catal. 2012, 2 (8), 1624−1634.

Metal Diatomics, Hydrides, Fluorides and Chlorides, and Comparison with First-Order Relativistic Calculations. J. Chem. Phys. 1998, 109 (2), 392−399. (33) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu. J. Chem. Phys. 2010, 132 (15), 154104. (34) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the Damping Function in Dispersion Corrected Density Functional Theory. J. Comput. Chem. 2011, 32 (7), 1456−1465. (35) Riplinger, C.; Neese, F. An Efficient and near Linear Scaling Pair Natural Orbital Based Local Coupled Cluster Method. J. Chem. Phys. 2013, 138 (3), 034106. (36) Riplinger, C.; Pinski, P.; Becker, U.; Valeev, E. F.; Neese, F. Sparse MapsA Systematic Infrastructure for Reduced-Scaling Electronic Structure Methods. II. Linear Scaling Domain Based Pair Natural Orbital Coupled Cluster Theory. J. Chem. Phys. 2016, 144 (2), 024109. (37) Mondal, B.; Neese, F.; Ye, S. Control in the Rate-Determining Step Provides a Promising Strategy To Develop New Catalysts for CO2 Hydrogenation: A Local Pair Natural Orbital Coupled Cluster Theory Study. Inorg. Chem. 2015, 54 (15), 7192−7198. (38) Mondal, B.; Neese, F.; Ye, S. Toward Rational Design of 3d Transition Metal Catalysts for CO2 Hydrogenation Based on Insights into Hydricity-Controlled Rate-Determining Steps. Inorg. Chem. 2016, 55 (11), 5438−5444. (39) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. Ab Initio Calculation of Vibrational Absorption and Circular Dichroism Spectra Using Density Functional Force Fields. J. Phys. Chem. 1994, 98 (45), 11623−11627. (40) Becke, A. D. Density-functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98 (7), 5648−5652. (41) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37 (2), 785−789. (42) Vosko, S. H.; Wilk, L.; Nusair, M. Accurate Spin-Dependent Electron Liquid Correlation Energies for Local Spin Density Calculations: A Critical Analysis. Can. J. Phys. 1980, 58 (8), 1200− 1211. (43) Marenich, A. V.; Cramer, C. J.; Truhlar, D. G. Universal Solvation Model Based on Solute Electron Density and on a Continuum Model of the Solvent Defined by the Bulk Dielectric Constant and Atomic Surface Tensions. J. Phys. Chem. B 2009, 113 (18), 6378−6396. (44) Neese, F.; Wennmohs, F.; Hansen, A.; Becker, U. Efficient, Approximate and Parallel Hartree-Fock and Hybrid DFT Calculations. A “chain-of-Spheres” Algorithm for the Hartree-Fock Exchange. Chem. Phys. 2009, 356 (1−3), 98−109. (45) Poli, R.; Harvey, J. N. Spin forbidden chemical reactions of transition metal compounds. New ideas and new computational challenges. Chem. Soc. Rev. 2003, 32 (1), 1−8. (46) Harvey, J. N. Understanding the kinetics of spin-forbidden chemical reactions. Phys. Chem. Chem. Phys. 2007, 9 (3), 331−343. (47) Marcus, R. A. On the Theory of Oxidation-Reduction Reactions Involving Electron Transfer. I. J. Chem. Phys. 1956, 24 (5), 966−978. (48) Klimkans, A.; Larsson, S. Reorganization Energies in Benzene, Naphthalene, and Anthracene. Chem. Phys. 1994, 189 (1), 25−31. (49) Fernandez, L. E.; Horvath, S.; Hammes-Schiffer, S. Theoretical Analysis of the Sequential Proton-Coupled Electron Transfer Mechanisms for H2 oxidation and Production Pathways Catalyzed by Nickel Molecular Electrocatalysts. J. Phys. Chem. C 2012, 116 (4), 3171−3180. (50) Liu, Y. P.; Newton, M. D. Reorganization Energy for Electron Transfer at Film-Modified Electrode Surfaces: A Dielectric Continuum Model. J. Phys. Chem. 1994, 98 (29), 7162−7169. (51) CRC Handbook of Chemistry and Physics, 90th ed. (Internet Version 2010); Lide, D. R., Haynes, W. M., Eds.; CRC Press/Taylor and Francis, Boca Raton, FL, 2010. (52) Pascual-Ahuir, J. L.; Silla, E.; Tuñon, I. GEPOL: An improved description of molecular surfaces. III. A new algorithm for the 7977

DOI: 10.1021/acs.inorgchem.9b00719 Inorg. Chem. 2019, 58, 7969−7977