Understanding potential-dependent competition between electrocatalytic dinitrogen and proton reduction reactions

A key challenge to realizing practical electrochemical N2 reduction reaction (NRR) is the decrease in the NRR activity before reaching the mass-transfer limit as overpotential increases. While the hydrogen evolution reaction (HER) has been suggested to be responsible for this phenomenon, the mechanistic origin has not been clearly explained. Herein, we investigate the potential-dependent competition between NRR and HER using the constant electrode potential model and microkinetic modeling. We find that the H coverage and N2 coverage crossover leads to the premature decrease of NRR activity. The coverage crossover originates from the larger charge transfer in H+ adsorption than N2 adsorption. The larger charge transfer in H+ adsorption, which potentially leads to the coverage crossover, is a general phenomenon seen in various heterogeneous catalysts, posing a fundamental challenge to realize practical electrochemical NRR. We suggest several strategies to overcome the challenge based on the present understandings.

In the manuscript titled "Understanding Potential-dependent Competition Between Electrocatalytic Dinitrogen and Proton Reduction Reactions", Choi et al. reported a computational study of nitrogen reduction reaction (NRR). By combining density functional theory (DFT) calculations and kinetic model, the authors concluded that the reason for the early loss of NRR efficiency attributes to the loss of N2 adsorption when H adsorption becomes stronger at more negative applied potential. Thus, hydrogen evolution reaction (HER) becomes dominant. Overall, it looks to me that this manuscript is worth publishing. But it seems that the evidence currently provided by the authors cannot fully support their conclusions as they expected. Based on my expertise and understanding of NRR and HER. I raise several comments for the authors to consider as follows. 1. I am more skeptical of the calculation results shown in Figure 5. According to the hypothesis proposed by the authors, competitive adsorption between N2 and H* is the major cause of the competition between NRR and HER. More specifically, N2 adsorption is less dependent on voltage. Therefore, the relationship between the adsorption energy and the voltage change is not obvious. Instead, the adsorption of H highly depends on voltage. Therefore, its adsorption energy has a significant correlation with voltage. The above analysis can explain the calculation results in Figure  4, but not Figure 5. According to the above assumptions, Figure 5 should also appear similar trend as shown in Figure 4. Although the author emphasizes that CHE has limitations in describing charge transfer for transition states, the CHE model should capture the above-mentioned dependence of adsorption energy and voltage. I can't understand the significant differences between Figure 4 and Figure 5. I suggest the authors add more discussions and explanations.
2. What is the physical meaning of dG/U in Figure 6? Is it a capacitance?
3. I cannot agree with the authors in the discussions of Figure 1. Both of the above reactions are multiple-step reactions (NRR is 6 e-and HER is 2e-). The overall reduction potential, in most cases, is known to be not related to the performance. Therefore, I can not be convinced that "NRR should have more negative potential than CO2RR" as the authors concluded.
Reviewer #2 (Remarks to the Author): Choi et al used electronic structure calculations in order to understand the competition between electrochemical ammonia production and hydrogen evolution on a single-site iron catalyst. This is an interesting and important research topic and I believe the use of grand-canonical DFT methods can provide significant new insights to processes like these. Unfortunately, I found this work to suffer from many severe flaws (discussed below) and to not provide meaningful insights to the literature.
Major flaws: -The free energy diagrams found in Figure SI-3 form the basis for most of this manuscript. Here, they compare the diagrams produced from the well-known computational hydrogen electrode (CHE) model to what they are referring to as the constant-electrode potential (CEP) model. It is shocking how different the results are between the two models, especially since the CHE model is thermodynamically exact in the limit of integer electron transfer (or alternatively, of zero change in surface-normal dipole). If their CEP model is thermodynamically consistent, I would expect deviations from the CHE model to be small, perhaps on the order of 0.1 eV (for example, if an electron transfer was 0.9 instead of 1). Yet we can observe dramatic differences: for example, the step *NHNH -> *NHNH2 appears uphill by ~0.5 eV in the CHE model and downhill by about the same amount in the CEP model. This implies that CHE made a 1-electron error in this step! This is hard to believe; if this is really the case, the authors need to fully explain the difference and the physical meaning of such.
-The mechanism they have employed is not justified. There are many possible pathways that this reaction could follow, even on a single-site catalyst, and the authors have only listed one, without justification.
-Related, the free energy diagrams have a step that appears to liberate NH2. (*NH2NH2 -> *NH2 + NH2) This is a free radical! This seems very unlikely. Why was only this single mechanism explored? There are many ways that N-N bond could be broken, and I do not expect the liberation of a radical from the surface would be the most favorable way. (Additionally, they calculate the liberation of this radical to be exergonic, which seems suspect.) -In Figure 2, the transition states corresponding to part (c), (e) and (f) look quite similar to the final states. I.e., I do not see any significantly stretched or distorted bonds. How do the authors confirm that TS obtained are the true transition states? For clarity and better understanding, the authors should include the charge/electron transfer associated between the IS, TS and FS, as well as potential-energy diagrams of the minimum energy pathway.
-If I understand it correctly, Figure 3(a) is pure thermodynamics and should align with CHE in the limit of integer electron transfer; that is, deviations from integer slope should be caused by noninteger electron transfer. To what do the authors attribute a ~0.2 e transfer for *N2 adsorption? I would have suspected that this elementary step doesn't involve significant charge transfer. Similarly, why are the charge transfers associated with the other two steps not closer to 1? [As a more minor point, how are the slopes calculated? The data points are not connected by a single straight line.] -Similarly, Figure 6 does not make sense to me. (The authors are plotting nearly the same quantity in the x and y directions.) It is bizarre to see a huge window in the charge transfer associated with *N2 adsorption (ranging from 0 to almost 0.6 e). How do the authors justify a charge transfer of 0.6 e for a step that doesn't involve explicit charge transfer? Similarly, for *H adsorption, the charge-transfer window varies over a huge window from roughly 1 to 2 e; how could the authors justify the 2-electron case (almost equivalent to the total electrons transferred in the entire HER)?
-The microkinetic results seem are suspect, since they didn't calculate barriers for many of their steps. (It seems they only have barriers for the couple of steps that they declared to be "potential determining" in the manuscript; I do not see others reported.) I personally would say that all barriers should be included in such a model since they are now rather straightforward with constant-potential methodologies; although I agree that it's likely that setting some to zero won't affect the kinetics. However, I can see a couple of apparently significant omissions just from looking at the free energy diagrams: * At 0 V, the highest-energy intermediate is *NHNH2 or *NHNH, depending on if you examine the CHE or CEP model. (As noted earlier, the difference in elementary thermodynamics between these two approaches is troubling.) Therefore, it's a good guess that the largest *overall* barrier (measured from the most-abundant reactive intermediate, MARI) would be from either the formation or destruction of these high-energy species. The largest *overall* barrier, not the largest individual barrier, is what typically best approximates the kinetics, so the omission of these barriers is significant. * There are non-electrochemical steps involved, and notably *NH2NH2 -> *NH2 + NH2. Since this step will not be affected by potential (to a first order), its barrier will not diminish with potential as the researchers are implicitly assuming other steps will. Therefore, this barrier and ones like it must be calculated. It is a bit strange to examine N2 reduction without examining breaking of the N-N bond.
-They seem to give two different explanations for the H2 production overwhelming N2 reduction (which is advertised to be the significant conclusion of this manuscript). On page 8, they attribute it to the higher coverage of H as the potential gets more negative (apparently blocking sites for N2 adsorption). Later on the same page they seem to be making the argument that it is the barrier to form *NNH, as compared to the barrier for *H, that limits. Which is it? I don't expect two different steps to be equally responsible.

Significant issues:
There are many other significant issues with this manuscript, which are too numerous to mention. I will list some of them here, but I do not mean to imply this is a complete list.
-The free energy diagrams should be in the manuscript itself; it's too confusing to have to look to ill-explained SI figures to follow the logic of the discussion in the manuscript.
-In Supplementary figure 6, why is the red plot (-0.5 V) not half-way between the black (0 V) and blue (-1 V)? -It appears that for the competitive H2 evolution reaction, the authors only considered the Volmer-Heyrovsky mechanism and ignored Volmer-Tafel. Volmer-Tafel is thought to dominate on Pt, for example. Were both mechanisms calculated? If there's reason to believe that Volmer-Tafel is inactive, please provide it.
-I found it confusing that the authors were including their results (for the limiting potential, U_L) already in their introduction, without telling us how they calculated these quantities or what they assumed in doing so.
-The authors spend a decent amount of time in the introduction saying why the CHE model is not good enough, yet then on page 8 they make claims about what the potential-determining steps are in N2 and proton reduction, which as far as I can tell are based solely on CHE logic. (That is, the "last step to turn downhill".) Since they can potentially make a much better analysis with the CEP model, why are they relying on CHE here? -In line 106, the authors cite 6 articles for the CEP model and explain it very briefly. There is no detailed description of the CEP model in the "methods" section as well. It is confusing to readers from which of these six articles the CEP model is coming from and how exactly it works. It is important to cite the original paper that developed the CEP model and explain in detail if authors used any modification of that method in this paper. Also, it is a good idea to discuss the method in the "method" section of the paper as it is the core of all calculations done; i.e., only the expert reader will recognize that the paragraph on "Implicit solvation...of the electrolyte" is referring to the CEP model.
-In Line 118 and elsewhere, ("Further calculations on other catalysts reveal that larger charge transfer..."). This is a general comment to the authors. There is no data available that is used to make plots in the manuscript. While charge transfer is discussed as an important parameter to decide the slope of delta_G vs U plot, there is no availability of the charge transfer values for any reaction in the manuscript or supporting information. It is important that all the data that is used to make any plot is available either in the form of mathematical equation (e.g. coverage) or tables (e.g. Gibbs free energy, charge transfer etc.).
-In microkinetic modeling, it is crucial to have a table specifying all rate constants and Gibbs free energies at all three pH and different potentials; or some way for the interested reader to reproduce the results.
-I found the comparison of CO2RR with NRR to be confusing. It is unclear what the authors reasoning towards the premature decrease in NRR activity which is not seen for CO2RR (when in competition with HER). In line 73, the authors suggest that the premature decrease in NRR activity should be attributed to the intrinsic properties of catalysts. However, in the rest of the manuscript, they do not mention what properties of the catalysts that could result in such a behavior. In fact, the rest of the manuscript only compares the elementary steps of NRR and HER to support their reasoning. I do not really understand what the CO2RR comparison was doing for this manuscript.
-The authors use the term coverage (theta) for N2 and H in the micro-kinetic modeling section. Also, in Figure 4, the y-axis represents a continuously varying coverage. However, several systems tested in this manuscript are single-site catalyst (for example Fe@N4). How is coverage defined for such a system? -How is pH incorporated when making the free energy surfaces? It is very clear that the bulk pH and local pH are very different. How do the authors treat the local pH? Also, based on my understanding the competition of NRR with HER should be more rigorous under acidic conditions provided the abundance of H+ to block the active sites. However, this should be less sever under alkaline conditions.
-The CHE definition in Supplementary Note 2 seems problematic. It is not the chemical potential of H+, but rather of the (H+ and e-) pair that can be related to the mu(H2) at certain conditions. Also, credits should be given to the correct reference for CHE (i.e., Norskov 2004).
Reviewer #3 (Remarks to the Author): The authors have written a very nice manuscript and should be commended for producing a very clean, concise, and well-reasoned work of science. As I outline below, the results of this paper are somewhat expected. The main contributions of this paper are 1. The insufficiency of the CHE model in describing the potential dependence of the species relevant to NRR and 2. Per microkinetic modeling too much importance is put on the limiting potentials from DFT. This work clearly outlines and corrects a flaw in the current thinking and literature around NRR electrocatalysis by testing the CHE model and introducing microkinetic modeling. I thus recommend this manuscript for publication with minor revisions.

Specific Comments:
The authors discuss the general trend of charge transfer being larger for H* adsorption than N2* adsorption. While the demonstration of this trend is valuable, the trend itself is not surprising. N2 is known to be very inert, and thus typically physisorbs on surfaces and thus has very limited charge transfer. Additionally, the greater dependence of H* on applied potential relative to N2* is also expected. The adsorption energies of N2* and H* are typically modeled using the computational hydrogen electrode model. In this model the adsorption of N2 is assumed to be independent of applied potential, necessitating a crossover between H* and N2* adsorption energies at some applied potential. The authors find this crossover as expected. While authors show in their manuscript that the CHE model is not sufficient to describe the real physics I present it here just to explain the expectations a reader would have. However, the differing potential dependence of N2H* relative to H* as well as N2* having significant potential dependence at all is interesting and not expected. Much of this lengthy discussion is outlined in [1] below (reference 56 in the manuscript.) I am somewhat surprised this work is not referenced more extensively, as this piece cleanly picks up the thread of this previous work.
The authors examine the three surface species that are important for NRR: N2*, N2H*, and H*. However, there is a 4th important species: NH2*. NH2* is often seen control the other half of the volcano plot when scaling relations are generated [2]. Thus, the authors cannot claim to have fully examined the species relevant for NRR. This limitation should be noted.
The section titled "Origin and descriptor for different slopes" is an interesting discussion of why there are differing slopes in the adsorption vs potential of various species with respect to their extent of charge transfer. However, as the section lays out this is more a demonstration of the nerst equation than a new discovery. This is not really a flaw, however it should not be posed as a major discovery.

Reviewer #1 (Remarks to the Author):
In the manuscript titled "Understanding Potential-dependent Competition Between Electrocatalytic Dinitrogen and Proton Reduction Reactions", Choi et al. reported a computational study of nitrogen reduction reaction (NRR). By combining density functional theory (DFT) calculations and kinetic model, the authors concluded that the reason for the early loss of NRR efficiency attributes to the loss of N 2 adsorption when H adsorption becomes stronger at more negative applied potential. Thus, hydrogen evolution reaction (HER) becomes dominant. Overall, it looks to me that this manuscript is worth publishing.
But it seems that the evidence currently provided by the authors cannot fully support their conclusions as they expected. Based on my expertise and understanding of NRR and HER. I raise several comments for the authors to consider as follows.
Reply: We are greatly grateful to the reviewer for positive comments on our work as well as very helpful suggestions. We respond to each question the reviewer brings one by one as below. The replies in this letter are marked in blue, while the revised texts in the manuscript and supporting information are highlighted in yellow.
1. I am more skeptical of the calculation results shown in Figure 5. According to the hypothesis proposed by the authors, competitive adsorption between N 2 and H* is the major cause of the competition between NRR and HER. More specifically, N 2 adsorption is less dependent on voltage. Therefore, the relationship between the adsorption energy and the voltage change is not obvious. Instead, the adsorption of H highly depends on voltage.
Therefore, its adsorption energy has a significant correlation with voltage. The above analysis can explain the calculation results in Figure 4, but not Figure 5. According to the above assumptions, Figure 5 should also appear similar trend as shown in Figure 4. Although the 2 author emphasizes that CHE has limitations in describing charge transfer for transition states, the CHE model should capture the above-mentioned dependence of adsorption energy and voltage. I can't understand the significant differences between Figure 4 and Figure 5. I suggest the authors add more discussions and explanations.
Reply: We thank the reviewer for this insightful comment. As the reviewer points out, the H adsorption energy becomes more negative with more negative potential in both the CEP and CHE models. However, the key difference lies in the different slopes with which the energy (and activation energy) lowering occurs for Volmer reaction (* + H 2 O → *H + OH -) vs.
Specifically, in the CEP model in which the charge transfer is rigorously calculated, the charge transfer in the TS of Volmer reaction and that of the Heyrovsky reaction (e.g., activation energy slope) is 0.96 and 0.56, respectively. Also, ∆G(*H) decreases with the slope of 1.34, higher than that of the CHE model (1.0). Thus, with increasingly more negative U, the rate for the Volmer reaction (*H formation) becomes faster than that for the Heyrovsky reaction (*H elimination), which leads to the potential-dependent accumulation of the *H coverage.
However, in the CHE model, since the same potential-dependence of the reaction energetics is assumed as long as the reactions involve the same number of electrons, both Volmer (*H formation) and Heyrovsky (*H elimination) reactions formally involving 1-e transfer, have the same potential-dependent behaviors. As a result, the accumulation of the H coverage is not observed with potential. To verify it, we plotted the MKM results using the CHE model for a hypothetical case (inspired by the CEP results) in which the activation energy of Volmer and that of Heyrovsky reactions have different slopes of 0.9 and 0.5, respectively, as well as the identical slopes of 0.5 (a CHE case shown in the original manuscript). As the comparison is shown in Supplementary Fig 16 and below, due to the identical potential-dependence 3 assumed in CHE model, the accumulation of H coverage with U is not observed in usual treatment of the same slopes, and the H-coverage crossover with potential is only observed in a hypothetical case in which the slopes for Volmer and Heyrovsky reactions are different, a situation that only the CEP model can treat rigorously.
We now added such discussions in the main text (in Page 19) with Supplementary Fig. 16.

Supplementary Fig. 16 | The change in
, and by U obtained by the microkinetic modeling (MKM) using the CHE model with different charge transfer in TS of *H formation. The relative and coverage are shown in upper and lower panels, respectively. Electrode potential (U) is in RHE scale at pH = 13.
2. What is the physical meaning of dG/U in Figure 6? Is it a capacitance?
Reply: The slope (∆G/U) is change of ∆G with respect to the change of potential and its physical meaning (∆G (in eV) / ∆U (in V)) is the amount of transferred electrons (e -) during chemical reaction.
3. I cannot agree with the authors in the discussions of Figure 1. Both of the above reactions 4 are multiple-step reactions (NRR is 6 e-and HER is 2e-). The overall reduction potential, in most cases, is known to be not related to the performance. Therefore, I cannot be convinced that "NRR should have more negative potential than CO 2 RR" as the authors concluded.
Reply: We thank and agree with the reviewer's comment. We now revised the discussion of "The experimental CO 2 RR activity increases with more negative potential and the maximum CO 2 RR activity is observed at around -0.7 V (Fig. 1b) 24 . At U = -0.7 V which is more negative than the theoretical limiting potential U L (-0.32 V), CO2RR can be sufficiently facilitated and its activity begins to decrease due to approaching the mass-transfer limit. Thus, the potential-dependent CO 2 RR activity can be qualitatively explained by conventional DFT calculations. NRR activity also increases with more negative potential at first, however, it begins to decrease quickly at -0.4 V (pH = 7.2) or -0.05 V (pH = 13) ( Fig. 1b), much earlier before reaching its U L (-1.29 V). This result indicates that NRR activity prematurely decreases with increasing reduction potential while its kinetics has not reached its expected theoretically maximum. NRR shows an unusual potential-dependent behavior which is unexplained by the conventional DFT calculations. Thus the premature decrease in the NRR activity should be attributed to the intrinsic properties of catalysts."

Reviewer #2 (Remarks to the Author):
Choi et al used electronic structure calculations in order to understand the competition between electrochemical ammonia production and hydrogen evolution on a single-site iron catalyst. This is an interesting and important research topic and I believe the use of grand-canonical DFT methods can provide significant new insights to processes like these.
Unfortunately, I found this work to suffer from many severe flaws (discussed below) and to not provide meaningful insights to the literature.
Reply: We are greatly grateful to the reviewer for many very helpful suggestions to enrich the perspectives of our paper. As per the reviewer's suggestions, we included the activation energies of all the other reaction steps in the MKM, and extensively clarified the physical origin for the large difference between the CEP and CHE models. We respond to each question the reviewer brings one by one as below. The replies in this letter are marked in blue, while the revised texts in the manuscript and supporting information are highlighted in yellow.
Major flaws: 1. The free energy diagrams found in Figure  Here, they compare the diagrams produced from the well-known computational hydrogen electrode (CHE) model to what they are referring to as the constant-electrode potential (CEP) model. It is shocking how different the results are between the two models, especially since the CHE model is thermodynamically exact in the limit of integer electron transfer (or alternatively, of zero change in surface-normal dipole). If their CEP model is thermodynamically consistent, I would expect deviations from the CHE model to be small, perhaps on the order of 0.1 eV (for example, if an electron transfer was 0.9 instead of 1). Yet 6 we can observe dramatic differences: for example, the step *NHNH -> *NHNH 2 appears uphill by ~0.5 eV in the CHE model and downhill by about the same amount in the CEP model. This implies that CHE made a 1-electron error in this step! This is hard to believe; if this is really the case, the authors need to fully explain the difference and the physical meaning of such.
Reply: We thank the reviewer for this comment. The reaction energy of *NHNH + (H + + e -) → *NHNH 2 the reviewer is referring to, we are guessing, is -0.30 eV in the CEP model and 0.33 eV in the CHE model (instead of 0.5 eV). The physical origin for such a large discrepancy (0.63 eV) is the change in the potential of zero charge (U PZC ) during reaction.
The U PZC represents the electrode potential of the slab model in neutral state (net charge = 0).
For example, the U PZC of Fe@N 4 (denoted as U PZC (*)) is -0.83 V and that of *H (U PZC (*H)) is -0.54 V. To set the U of * and *H to 0 V, electrons amounting to 0.83 V is extracted in * while electrons amounting to 0.54 V is extracted in *H. Consequently, extra electrons corresponding to 0.29 V (0.83 -0.54) are engaged to compensate the change of U PZC . In the CHE model, all intermediates are assumed to be at same U, and thus, the charge due to the change of U PZC is not included. We compared the change in U PZC (∆U PZC ) with ∆∆G CEP-CHE .
A linear correlation of ∆U PZC with ∆∆G CEP-CHE was obtained, indicating that the discrepancy mainly originates from the change in U PZC during reaction. We found that ∆U PZC (*NHNH) is -0.03 V, while ∆U PZC (*NHNH 2 ) is -0.26 V. Thus, the ∆∆G CEP-CHE (*NHNH) (-0.09 eV) is much closer to 0 than ∆∆G CEP-CHE (*NHNH 2 ) (-0.71 eV), leading to the large discrepancy in *NHNH → *NHNH 2 step between the CEP model and CHE model. 2. The mechanism they have employed is not justified. There are many possible pathways that this reaction could follow, even on a single-site catalyst, and the authors have only listed 8 one, without justification.
Reply: We thank the reviewer for this highly relevant and helpful comment to enrich our results. We now considered all possible reaction intermediates listed in Fig. 2. Combined with the Comment 7, we added activation free energies. The lowest activation energy pathway is highlighted with red line in Fig. 3. Consideration of these additional energetics of intermediates do not change the conclusions. The number of transferred protons is listed in the first row. For *NH + NH 3 , *NH 2 + NH 3 , *NH 3 + NH 3 , NH 3 is omitted for the clarity. 3. Related, the free energy diagrams have a step that appears to liberate NH 2 (*NH 2 NH 2 -> *NH 2 + NH 2 ). This is a free radical! This seems very unlikely. Why was only this single mechanism explored? There are many ways that N-N bond could be broken, and I do not expect the liberation of a radical from the surface would be the most favorable way.
(Additionally, they calculate the liberation of this radical to be exergonic, which seems suspect.) Reply: We thank the reviewer for pointing this out. The "*NH 2 + NH 2 " was a typo. The correct mechanism which we considered is *NH 2 NH 2 + (H + + e -) → *NH 2 + NH 3 , liberating NH 3 molecule. We now corrected the typos in all free energy diagrams. Figure 2, the transition states corresponding to part (c), (e) and (f) look quite similar to the final states. I.e., I do not see any significantly stretched or distorted bonds. How do the authors confirm that TS obtained are the true transition states? For clarity and better 10 understanding, the authors should include the charge/electron transfer associated between the IS, TS and FS, as well as potential-energy diagrams of the minimum energy pathway.

In
Reply: We verified the TS to be a first order saddle point by performing vibrational analysis.
The TS structure of *NNH formation is especially difficult to visually differentiate from its FS since the *NNH formation is highly endothermic, and hence, the geometry as well as charge transfer in the TS are very similar to those of FS. This agrees with the Hammond postulate [3], that is, TS of endothermic reaction resembles FS (and vice-versa). We now added the discussion for the geometry of TS (in Page 12), images for reaction pathway (Fig.   4), charge transfer and energy diagram ( Supplementary Fig. 4, 6 and Supplementary Table 4) for the clarity.
[3] A Correlation of Reaction Rates. J. Am. Chem. Soc. 1955, 77, 334 Fig. 6 in below). This result indicates that more negative potential promotes charge transfer to *N 2 and *N 2 activation.
Due to the geometry change by U (which will be discussed during the response to Comment 10 below), the data points of *N 2 and *NNH are not completely in a single line in a wide potential range (0 V ~ -1 V). The slopes are obtained by a simple first order linear regression.
In Fig. 3a (Fig. 5a in the revised manuscript), we focused on the trend in ∆G vs. U rather than comparing the exact slopes. We note that more narrow range of U (e.g. 0 V ~ -0.5 V) is used to calculate accurate slope to construct the potential-dependent energetics in the MKM.
These discussions are now added in Page 10-13 with Supplementary Fig. 5   ∆G(*N 2 ) and ∆G(*N 2 →*NNH) on various catalysts. We found that the deviation of the charge transfer well agrees with ∆U PZC (Supplementary Fig. 19 in below). This result indicates that catalysts whose U PZC easily changes during chemical reaction, requires extra (or deficient) electrons during electrochemical reaction.
We also note that such a large deviation in charge transfer has also been reported on SAC and formation or destruction of these high-energy species. The largest *overall* barrier, not the largest individual barrier, is what typically best approximates the kinetics, so the omission of these barriers is significant. * There are non-electrochemical steps involved, and notably *NH 2 NH 2 -> *NH 2 + NH 2 . Since this step will not be affected by potential (to a first order), its barrier will not diminish with potential as the researchers are implicitly assuming other steps will. Therefore, this barrier and ones like it must be calculated. It is a bit strange to examine N 2 reduction without examining breaking of the N-N bond.
Reply: We thank the reviewer for helpful suggestions to strengthen our results. Except for the adsorption/desorption of N 2 and NH 3 , activation energies are now calculated for all reaction steps when using H 2 O as a proton donor. We note that highly exothermic reactions such as *NH + H 2 O → *NH 2 + OHand *NH 2 + H 2 O → *NH 3 + OHproceed barrierlessly. In the lowest activation energy pathway at 0 V (red line in Fig. 3a below), the *NHNH formation has a higher apparent activation energy than the *NNH formation as in the reviewer's expectation. Thus, we now included the activation energies of all reaction steps in the MKM simulations under neutral and alkaline conditions. The revised MKM results at pH 7.2 and 13 are shown in Fig. 6 as below.
Since the apparent activation of *NNH formation decreases slower than that of *NHNH formation with more negative U, the *NNH formation begins to show the highest apparent activation energy (i.e. RDS) at -0.12 V ( Fig. 3b and 3c). In this study, we mainly discuss NRR under pH = 13 and pH = 7.2, the same conditions with the reported experiments. Thus, we simplified the MKM simulations under acidic conditions by considering only *NNH formation in activation energy, assuming that *NNH formation determines overall kinetics of NRR under negative potential. At all three conditions considered in this study, the premature decrease of ammonia yield rate and N 2 coverage still hold in the revised MKM results (Fig.   19 6a) with the precise U at maximum NH 3 yield rate or N 2 coverage changing only slightly by ~ 0.05 V. 8. They seem to give two different explanations for the H 2 production overwhelming N 2 reduction (which is advertised to be the significant conclusion of this manuscript). On page 8, they attribute it to the higher coverage of H as the potential gets more negative (apparently blocking sites for N 2 adsorption). Later on the same page they seem to be making the argument that it is the barrier to form *NNH, as compared to the barrier for *H, that limits.

21
Which is it? I don't expect two different steps to be equally responsible.
Reply: The coverage crossover between H and N 2 is mainly responsible for the premature decrease of NRR activity (i.e. decreasing NRR activity with more negative U). However, higher *NNH formation barrier than *H formation explains why the NRR kinetics is slower than HER in absolute magnitudes. We now clarified the role of coverage crossover and *NNH formation vs. *H formation.

Significant issues:
There are many other significant issues with this manuscript, which are too numerous to mention. I will list some of them here, but I do not mean to imply this is a complete list.
9. The free energy diagrams should be in the manuscript itself; it's too confusing to have to look to ill-explained SI figures to follow the logic of the discussion in the manuscript.
Reply: We now moved the potential-dependent free energy diagram of NRR at U = 0, -0.23 and -0.5 V (vs. RHE in pH 13) to the main text as Fig. 3. 10. In Supplementary figure 6, why is the red plot (-0.5 V) not half-way between the black (0 V) and blue (-1 V)?
Reply: The reason for that the reaction energy at -0.5 V is not half-way between that at 0 V and -1V, is the change of *NNH geometry with U. We optimized all structures with the proper number of electrons in the slab model and found that *NNH optimizes to a more bent structure with negative potential. As a comparison, we also obtained the potential-dependent *NNH formation energy when the geometry of * and *NNH is fixed at their optimized geometries at neutral state. Without the needed potential-dependent geometry relaxation, the 22 slope (∆G vs. U) is constant ( Supplementary Fig. 7 in below). This result indicates that the geometry relaxation by U leads to the potential-dependent charge transfer. We now added such discussion in the Page 12-13 with Supplementary Fig. 7 as below.   Supplementary Fig. 7 | The potential-dependent change of geometry and formation energy of *NNH. (a) The optimized geometry of *NNH at neutral state (-0.7 V), 0 V, -0.25 V, -0.5 V and -1 V (vs. SHE). (b) Potential-dependent *NNH formation energy. Red line represents formation energy obtained by single point calculation using optimized structure at neutral state, while black line represents *NNH formation energy at fully optimized structure. At -0.5 V to -1 V, noticeable changes of geometry and *NNH formation energy are not observed.
11. It appears that for the competitive H 2 evolution reaction, the authors only considered the Volmer-Heyrovsky mechanism and ignored Volmer-Tafel. Volmer-Tafel is thought to dominate on Pt, for example. Were both mechanisms calculated? If there's reason to believe that Volmer-Tafel is inactive, please provide it.
Reply: The Tafel reaction can proceed with two different pathways as illustrated in Supplementary Fig. 13. One is the reaction of H at Fe (*H) with H at adjacent N atom (*-NH) and the other is the reaction of two H's at Fe atom (*2H). We first compared the formation free energy of *-NH and *H by using the CEP model and found that *H is more stable than 23 *-NH by 0.92 eV (Supplementary Fig. 14a). Thus, we exclude adjacent N atom as an active site of HER. Next, we compared the activation energy of *2H formation (*H + H 2 O → *2H + OH -) and that of Heyrovsky reaction (*H + H 2 O → H 2 (g) + OH -) ( Supplementary Fig. 14b).
The activation energy for *2H formation is higher than that of Heyrovsky reaction by 0.34 eV at 0 V. Thus, we conclude that the Volmer-Tafel mechanism is less active, and thus considered the Volmer-Heyrovsky reaction. Above discussion are now added in the Supplementary Note 5 with Supplementary Fig. 14-15. 12. I found it confusing that the authors were including their results (for the limiting potential, U L ) already in their introduction, without telling us how they calculated these quantities or what they assumed in doing so.
Reply: The limiting potential (U L ) represents electrode potential where all individual electrochemical reactions become exothermic. In the CHE model, U L is obtained from the free energy change at potential determining step (denoted as ∆G PDS ) as U L = -∆G PDS / e (V).
We now added the definition of the limiting potential in the introduction as below.
"Here, we calculated U L by the computational hydrogen electrode (CHE) model, which has been the most widely used method in estimating the energetics of electrochemical reactions.
The U L is equal to the ∆G PDS / e, where ∆G PDS is the free energy change at the most uphill individual step (i.e. potential-determining step)." 13. The authors spend a decent amount of time in the introduction saying why the CHE model is not good enough, yet then on page 8 they make claims about what the potentialdetermining steps are in N 2 and proton reduction, which as far as I can tell are based solely on CHE logic. (That is, the "last step to turn downhill".) Since they can potentially make a much better analysis with the CEP model, why are they relying on CHE here?
Reply: The PDS is the highest overpotential requiring individual step among electrochemical reactions, which determines the overall thermodynamic overpotential. Thus, the PDS is the same in the CHE model and CEP model for both NRR and HER. To clarify, we now modified the sentence to "Next, we analyze the trend of NRR and HER activity by comparing G a (*N 2 →*NNH) and G a (*H), a RDS under negative U as well as PDS in NRR and HER on Fe@N 4 (Fig. 3 and Supplementary Fig. 4), respectively." 14. In line 106, the authors cite 6 articles for the CEP model and explain it very briefly. There is no detailed description of the CEP model in the "methods" section as well. It is confusing to readers from which of these six articles the CEP model is coming from and how exactly it works. It is important to cite the original paper that developed the CEP model and explain in detail if authors used any modification of that method in this paper. Also, it is a good idea to discuss the method in the "method" section of the paper as it is the core of all calculations done; i.e., only the expert reader will recognize that the paragraph on "Implicit solvation...of the electrolyte" is referring to the CEP model.
Reply: As per the reviewer's suggestion, we now clarified the explanation of the CEP model used here in the Introduction as well as method sections.

26
In the Introduction, "In this work, we use the constant electrode potential (CEP) model which treats the electrode-electrolyte interface as a polarizable continuum with implicit solvation model 34,35 . In this model, the number of electrons is adjusted to guarantee different states to have the same workfunction in the grand canonical states. This method has recently been used to understand many electrochemical reactions [36][37][38][39][40][41][42][43][44][45][46] ." In the method section, "We treat the electrode-electrolyte interface as a polarizable continuum via the linearized Poisson-Boltzmann equation which is implemented in the VASPsol 34,35 , extension of VASP. In this method, an ionic counter-charges are implicitly placed at the interface. The net dipole originates from adding (or extracting) electrons is screened by the ionic counter-charges." 15. In Line 118 and elsewhere, ("Further calculations on other catalysts reveal that larger charge transfer..."). This is a general comment to the authors. There is no data available that is used to make plots in the manuscript. While charge transfer is discussed as an important parameter to decide the slope of delta_G vs U plot, there is no availability of the charge transfer values for any reaction in the manuscript or supporting information. It is important that all the data that is used to make any plot is available either in the form of mathematical equation (e.g. coverage) or tables (e.g. Gibbs free energy, charge transfer etc.).
Reply: The charge transfer at each potential and the slope are now added as Supplementary Gibbs free energies at all three pH and different potentials; or some way for the interested reader to reproduce the results.
Reply: The potential-dependent energetics at all three pH used for the MKM, are now listed in the Supplementary Table 2. The formation energy (or activation energy) at certain U is expressed as ∆G(U) (or G a (U)) = αU + β, where the α and β represents the slope (or the amount of charge transfer in the reaction) and ∆G (or G a ) at 0 V, respectively. Thus, α and β of all reactions included in the MKM are provided. We used the transition state theory, thus the rate constants can be derived with just Gibbs free energy of activation. Reply: The main motivation of the comparison of NRR with CO 2 RR is to highlight that the unusual potential-dependent NRR activity is unexplained by the conventional DFT calculations for the theoretical limiting potential. The discussion of CO 2 RR vs. NRR is now revised as below in the Introduction.

Supplementary
"The CO 2 RR activity increases with more negative potential and the maximum CO 2 RR activity is observed at around -0.7 V (Fig. 1b). At U = -0.7 V which is more negative than U L (-0.32 V), CO 2 RR can be sufficiently facilitated and begins to decrease due to the masstransfer limit. Thus, the potential-dependent CO 2 RR activity can be qualitatively explained by conventional DFT calculations. NRR activity also increases with more negative potential at first, however, it begins to decrease quickly at -0.4 V (pH = 7.2) or -0.05 V (pH = 13) ( Fig.   1b), before reaching its U L (-1.29 V). This result indicates that NRR activity prematurely decreases with increasing reduction potential while its kinetics has not reached its expected theoretically maximum. NRR shows an unusual potential-dependent behavior which is unexplained by the conventional DFT calculations. Thus the premature decrease in the NRR activity should be attributed to the intrinsic properties of catalysts." 18. The authors use the term coverage (theta) for N 2 and H in the micro-kinetic modeling section. Also, in Figure 4, the y-axis represents a continuously varying coverage. Reply: Chemical potential of proton can be related to the that of H 2 as (H ) = (H ) + Φ − 0.059pH). We omitted the derivation of (H ) from the CHE model in original submission. The (H ) can be derived from the CHE model as below.
The 0 V (vs. SHE) is defined when the 1 bar of gas H 2 is equilibrium with the pair of proton and electron at pH = 0 (i.e. 1/2 H 2 ↔ H + + e -). At 0 V (vs. SHE), the chemical potential of a pair of proton and electron is equal to half of the chemical potential of H 2 (1), which is the well-known CHE model.
(H ) = (H ) + Φ − 0.059pH (5) We now added the above discussion in Supplementary Note 2 and cited the original paper for the CHE model.

Reviewer #3 (Remarks to the Author):
The authors have written a very nice manuscript and should be commended for producing a very clean, concise, and well-reasoned work of science. As I outline below, the results of this paper are somewhat expected. The main contributions of this paper are 1. The insufficiency of the CHE model in describing the potential dependence of the species relevant to NRR and 2. Per microkinetic modeling too much importance is put on the limiting potentials from DFT.
This work clearly outlines and corrects a flaw in the current thinking and literature around NRR electrocatalysis by testing the CHE model and introducing microkinetic modeling. I thus recommend this manuscript for publication with minor revisions.
Reply: We are very grateful to the reviewer for positive view on our work. We respond to each question the reviewer brings one by one as below. The replies in this letter are marked in blue, while the revised texts in the manuscript and supporting information are highlighted in yellow.
Specific Comments: 1. The authors discuss the general trend of charge transfer being larger for H* adsorption than N 2 * adsorption. While the demonstration of this trend is valuable, the trend itself is not surprising. N 2 is known to be very inert, and thus typically physisorbs on surfaces and thus has very limited charge transfer. Additionally, the greater dependence of H* on applied potential relative to N 2 * is also expected. The adsorption energies of N 2 * and H* are typically modeled using the computational hydrogen electrode model. In this model the adsorption of N 2 is assumed to be independent of applied potential, necessitating a crossover between H* and N 2 * adsorption energies at some applied potential. Reply: We thank the reviewer for this comment. Catalysts which lie on the left-leg of the NRR volcano plot (i.e. strongly N binding), the reaction *NH 2 + (H + + e -) → * + NH 3 has been known as the potential-determining step. For these catalysts, as the reviewer correctly points out, *NH 2 rather than *NNH should be considered. We now added this discussion in the "Origin and descriptor for different slopes" section (in Page 21) as below.
"Here, we estimate the overall kinetics of NRR by using *NNH formation which is the PDS of Fe@N 4 as well as various catalysts. However, catalysts with strongly N binding strength, which lies on the left-leg of the NRR volcano plot, are limited by *NH 2 + (H + + e -) → NH 3 6,62 .
Thus, we note the potential-dependent energetics and charge transfer associated with *NH 2 should be considered for strongly N binding catalysts." 3. The section titled "Origin and descriptor for different slopes" is an interesting discussion of why there are differing slopes in the adsorption vs potential of various species with respect to their extent of charge transfer. However, as the section lays out this is more a demonstration of the Nernst equation than a new discovery. This is not really a flaw, however it should not be posed as a major discovery.
[  Black, red, and blue colors represent *H formation, *N 2 formation and *NNH formation from *N 2 , respectively.

REVIEWERS' COMMENTS
Reviewer #1 (Remarks to the Author): In the revised manuscript, the authors have fully addressed the concerns as the reviewer raised. I have no more comments.
Reviewer #4 (Remarks to the Author): In this work, Choi et al use DFT to analyze the competition between HER and NRR on a singleatom FeN4 catalyst. They demonstrate that it is essential to consider the explicit charge state of the catalyst site using a grand canonical formalism to establish a constant electrode potential (CEP), which they compare to the computational hydrogen electrode (CHE) model. The CEP model captures non-integer charge transfer effects impacting the coverage dependence of H* and N2* on potential (i.e., capturing deviations from predictions made with the simpler CHE model). The balance between H* and N2* coverage varies with potential, which explains (in part) the experimentally observed competition between HER and N2R. I see that the manuscript has already been reviewed and extensively revised. I believe that the authors have thoroughly answered the questions and comments raised by the other reviewers. For this reason, this manuscript should be accepted once the following minor points are considered.
1) The authors state on Pg. 3 that "Theoretical studies have suggested that the theoretical limiting potential (UL) where the electrochemical reaction becomes exothermic…" This phrasing is vague and should be restated to be more precise. UL is the potential at which all elementary reaction steps become energetically downhill (as opposed to the potential where the overall reaction becoming exothermic).
2) The authors state that "The lowest energy pathway based on the apparent activation energy is represented by red line in Fig. 3." How is the "apparent activation energy" defined in this case? E.g., is it defined as it is in the energetic-span model of Kozuch and Shaik (Acc. Chem. Res. 2011, 44, 2, 101-110)? The authors use the term "apparent activation energy" several times without defining what they mean by it. Many different definitions for "apparent activation energy" have been used in the literature, so this is a point that needs to be clarified.
3) Fig. 3 is difficult to read and interpret because many lines are overlapping, and it is not clear which label corresponds to which line. This figure should be revised.

1
Reviewer #4 (Remarks to the Author): In this work, Choi et al use DFT to analyze the competition between HER and NRR on a single-atom FeN4 catalyst. They demonstrate that it is essential to consider the explicit charge state of the catalyst site using a grand canonical formalism to establish a constant electrode potential (CEP), which they compare to the computational hydrogen electrode (CHE) model. The CEP model captures non-integer charge transfer effects impacting the coverage dependence of H* and N2* on potential (i.e., capturing deviations from predictions made with the simpler CHE model). The balance between H* and N2* coverage varies with potential, which explains (in part) the experimentally observed competition between HER and N2R. I see that the manuscript has already been reviewed and extensively revised. I believe that the authors have thoroughly answered the questions and comments raised by the other reviewers. For this reason, this manuscript should be accepted once the following minor points are considered.
Reply: We are very grateful to the reviewer for the positive review on our work. We respond to each question of the reviewer one-by-one below. The replies in this letter are marked in blue, while the revised texts in the manuscript and supporting information are highlighted in yellow.
1) The authors state on Pg. 3 that "Theoretical studies have suggested that the theoretical limiting potential (UL) where the electrochemical reaction becomes exothermic…" This phrasing is vague and should be restated to be more precise. UL is the potential at which all elementary reaction steps become energetically downhill (as opposed to the potential where the overall reaction becoming exothermic).
Reply: We thank the reviewer for this correction. We now edited the sentence to "Theoretical 2 studies have suggested that the theoretical limiting potential (U L ) where the all electrochemical elementary steps become exothermic…" 2) The authors state that "The lowest energy pathway based on the apparent activation energy is represented by red line in Fig. 3." How is the "apparent activation energy" defined in this case? E.g., is it defined as it is in the energetic-span model of Kozuch and Shaik (Acc. Chem. Res. 2011, 44, 2, 101-110)? The authors use the term "apparent activation energy" several times without defining what they mean by it. Many different definitions for "apparent activation energy" have been used in the literature, so this is a point that needs to be clarified.
Reply: We thank the reviewer for this comment and agree that there are many different definitions for the apparent activation energy. Here, the apparent activation energy is defined as the energy difference between the highest transition state (TS) energy and the lowest energy intermediate in the catalytic cycle from the energetic span model. We now added such a sentence with the reference suggested by the reviewer in Page 9 and clarified the definition of apparent activation energy in this study.
3) Fig. 3 is difficult to read and interpret because many lines are overlapping, and it is not clear which label corresponds to which line. This figure should be revised.
Reply: We appreciate the reviewer's suggestion. We changed the color of a black line to gold when three lines overlap and clarified which label corresponds to which line.
4) The authors state in the rebuttal letter in response to comment 7 of reviewer 2 that "we included the activation energies of all the other reaction steps in the MKM." In Fig. 3 and Table S2 the authors report that some steps "proceed without a barrier." I am surprised to see