Download PDF
Research Article  |  Open Access  |  14 May 2023

High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction

Views: 730 |  Downloads: 359 |  Cited:   5
Chem Synth 2023;3:23.
10.20517/cs.2023.09 |  © The Author(s) 2023.
Author Information
Article Notes
Cite This Article


Transitional metal single atom (TM1) doped graphene catalysts have been widely applied in electrochemical N2 reduction reaction (NRR). However, it remains a challenge for the rational design of highly active and selective electrocatalysts owing to limited knowledge of structure-activity correlations. Here, we adopted first-principle calculations to high-throughput screen the NRR performance of TM1 coordinated with two boron and two nitrogen atoms in graphene (TM1-B2N2/G). A “five-step” strategy was implemented by progressively considering different metrics such as stability, N2 adsorption, N2 activation, potential-determining step, and selectivity. As a result, a volcano plot of reactivity is established by using the valence electron number of TM1 as the descriptor. Among all catalysts, Cr1-B2N2/G exhibits superior performance with a limiting potential of -0.43 V with high selectivity of NRR interpreted by better spatial symmetry and excellent compatibility in terms of energy when N2 interacts with TM1. Our work reveals the general strategy of computational efforts to predict the next generation of advanced catalytic materials for NRR.


Single-atom catalysts, ammonia synthesis, density functional theory, structure-activity relationship, descriptor


Ammonia (NH3) is an extremely important chemical feedstock in human society and is widely used in the production of industrial chemicals such as fertilizers and pharmaceuticals[1]. The Haber-Bosch method converts hydrogen (H2) and nitrogen (N2) to NH3 under high temperature (> 450 °C) and high pressure (> 100 atm), which consumes tremendous energy and releases plenty of CO2[2-4]. In past decades, the electrocatalytic N2 reduction reaction (NRR) has emerged as a promising technology for the production of NH3 using N2 and H2O as reactants at ambient conditions with much lower energy consumption[5-7]. However, the development of highly active and selective electrocatalysts for NRR remains a challenge due to the complexity of the reaction environment, which makes it difficult to accurately understand the correlation between activity and catalyst structure[8-13].

Single-atom catalysts (SACs) have exhibited outstanding specific activities for many electrochemical reactions[14-19]. Owing to high surface area and intriguing electrical properties, graphene-supported SACs have demonstrated excellent performance for NRR[20-22]. The well-defined structure in graphene-supported SACs offers great potential to modulate the active metal and corresponding coordination environment for establishing the structure-activity relationship and achieving high activity and selectivity[23-27]. The doping of B and N in the graphene structure has been reported to be beneficial for the adsorption and activation of reactants, thus enhancing the electrocatalytic activity of the catalyst. Due to the different electronegativity of B, N, and C, the interaction between active metal and substrate can be improved by electron transfer or charge rearrangement. Different B, N coordination environments can also modulate the chemical natures of active metal sites, giving rise to higher reactivity[28,29].

In this work, we investigated the NRR activity of all 3d transition metals in different B2N2 coordination environments by a high-throughput screening method. The comparative comparisons were performed in terms of stability, N2 adsorption, N2 activation, potential-determining step, and selectivity. Thus a volcano plot of reactivity was established by using the valence electron number of TM1 as the descriptor, which can be understood that the chemical nature of the active center (TM1 coordinated with B and N atoms) determines the electron donation process from TM1 to N2, which highly affect N2 activation and subsequent hydrogenations. Moreover, the scenario demonstrated in this work offers guidance for the synthesis of efficient NRR catalysts by experiment.

Computational details

All spin-polarization calculations were carried out by using the Vienna ab initio simulation package (VASP)[30,31]. The projector augmented wave (PAW) method was adopted to describe the ion-electron interactions[32]. The generalized gradient approximation in the Perdew-Burke-Ernzerhof (PBE) function was used[33,34], and the cut-off energy of 500 eV was set for the plane-wave basis. A 6 × 6 graphene supercell was used as the substrate with a 15 Å vacuum space and the Brillouin zone was sampled with the Monkhorst-Pack 3 × 3 × 1 and 7 × 7 × 1 k-point grids for geometry and electronic structure calculations, respectively. The cut-off energy, size of the graphene supercell, and the selection of k-points were referred to previous literature to ensure the accuracy of the calculation results [35-37]. The van der Waals (vdW) interactions between the adsorbates and catalysts are described by using Grimme’s semiempirical DFT-D3 scheme[38]. The solvation effect was evaluated under an implicit solvent model VASPsol[39]. The energy and force convergence thresholds for the iteration in self-consistent field (SCF) were set to 10-5 eV and 0.02 eV·Å-1, respectively. The calculation of COHP was done with the support of the Lobster program[40].

The adsorption energies (Eads) of the NRR intermediates were obtained by the following equation:

Eads = EtotalEclean catalystsEadsorbate (1)

In the above equation, the Etotal, Eclean catalysts, and Eadsorbate represent the total energies of the species-adsorbed catalysts system, clean catalysts, and adsorbate, respectively. According to this definition, the more negative adsorption energy indicates stronger adsorption of corresponding species.

The Gibbs free energy change (ΔG) for each reaction step was calculated by the computational hydrogen electrode (CHE) model proposed by Nørskov and co-workers[41], which uses one-half of the chemical potential of hydrogen as the chemical potential of the proton-electron pairs. The ΔG was obtained by the following equation:

ΔG = ΔE + ΔEZPETΔS + ΔGpH + eU (2)

In the above equation, the ΔE is reaction energy obtained directly from the DFT calculation. The ΔEZPE and TΔS are the contributions of the zero-point energy and entropy for the ΔG. T represents the temperature, which is set as 298.15 K. The parameters e and U represent the number of electrons and the applied electrode potential, respectively. ΔGpH is the free energy correction of pH, which can be determined as:

ΔGpH = kBT × pH × ln10 (3)

where the pH value was set to zero. The correction of the frequency was done with the support of the VASPKIT program[42]. In this work, the limiting potential (UL) value was used as a descriptor of catalytic activity, which was determined by the potential-determining step (PDS) among six proton-coupled electron transfer (PCET) processes. UL was obtained by:

UL = – ΔGmax/e (4)

Where ΔGmax represents the ΔG of the PDS.


Structure and Stability of TM1-B2N2/G

The 6 × 6 graphene supercell was first used as a substrate, four carbon atoms of which were replaced by two boron and two nitrogen atoms, giving rise to B/N doped graphene (denoted as B2N2/G). Then ten different 3d transition metal single atoms (denoted as TM1) were selected and embedded as active centers (Sc, Ti, V, Cr, Mn, Fe, Co, Ni, Co, Zn) and coordinated with B or N via three different manners including two cis- and one trans- forms, distinguished as TM1-B2N2/G-1 [Figure 1A], TM1-B2N2/G-2 [Figure 1B], and TM1-B2N2/G-3 [Figure 1C], respectively. Hence, a total of 30 different single-atom catalysts were constructed.

High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction

Figure 1. Structural stability of TM1-B2N2/G catalysts. (A-C) The geometric structure of TM1-B2N2/G catalysts; (D) The binding energy (Eb) of TM1-B2N2/G catalysts; (E) The difference between binding energy and cohesion energy (Eb - Ec) of TM1-B2N2/G catalysts.

To assess the stability of the above TM1-B2N2/G structures, we calculated the binding energies (Eb) of TM single atoms with substrates via equation (5) and the cohesion energy (Ec) by following equation (6) which refers to the energy of separating constituent metallic atoms apart from each other and assembling into the bulk:

Eb = EcatalystETMsingleEsupport (5)

Ec = ETMbulk/nETMsingle (6)

where ETMbulk, ETMsingle, and n represent the energy of the bulk crystal cell of the transition metal, the energy of a single transition metal atom in the vacuum, and the number of atoms in the cell, respectively. Thereupon, we compared both the calculated Eb and EbEc values in Figure 1D-E and Supplementary Table 1. The former metric indicates the binding strength of TM single atom with the support. For TM1-B2N2/G catalysts, all Eb < 0 [Figure 1D], which means that all TM single atoms connect strongly with boron and nitrogen in graphene, and more negative values of Eb underline that the cis-structures are more sturdy than trans-type. On the other side, Eb - Ec reflects the aggregation tendency of TM single atom, the more negative the values, the more stable the single atoms in the matrix of B2N2/G. Figure 1E exhibits that Eb - Ec < 0 for most of TM1-B2N2/G (except some Sc, Ti, V SACs), further demonstrating that the TM single atoms are anchored firmly rather than being subject to agglomeration. Additionally, the designing rationality of this work also could be verified by the fact that most of the structures above have been successfully synthesized by experiments and reported in the literature[43,44].

Screen the NRR activities of TM1-B2N2/G

After investigating the geometry and stability of the catalyst, we started to consider the adsorption and activation of N2. According to previous studies, the high-performance NRR process requires: (1) N2 spontaneously adsorbed on catalysts via chemical adsorption; (2) the free energy change of the first proton-coupled electron transfer (PCET) step should be appropriately low; (3) the limiting potential of the whole process should be low enough throughout the NRR process; (4) the generated NH3 on the catalyst surface can be easily desorbed after reaction; (5) the NRR process prevails while hydrogen evolution reaction (HER) is prohibited[45-47]. Hence, we followed the above hints to evaluate the NRR performance of TM1-B2N2/G.

N2 adsorption has a crucial influence on N2 activation and the subsequent reaction process. After a careful summary of the literature, we considered two adsorption modes here: end-on and side-on configurations [Figure 2A][48]. The corresponding energies of N2 adsorption are mapped in Figure 2B (Gad) and Supplementary Table 2 (Ead)-Supplementary Table 3(Gad). For most TM1-B2N2/G catalysts, the negative or near zero values suggest N2 can be spontaneously adsorbed on TM1 [Supplementary Table 2 and Supplementary Table 3]. Meanwhile, the related N≡N lengths of adsorbed *N2 showed significant increases (1.13~1.24 Å) compared to that of N2 in vacuum (1.11 Å) [Supplementary Table 4], confirming that N2 can be easily activated on these catalysts. In contrast, for Cu1-B2N2/G and Zn1-B2N2/G, both the relatively higher adsorption energies and negligible changes of N≡N lengths account for unsatisfactory N2 activation on Cu or Zn single site; thus, Cu1-B2N2/G and Zn1-B2N2/G would not be considered in following sections.

High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction

Figure 2. N2 Adsorption. (A) Adsorption mode of N2; (B) The heatmap of the free energy of adsorption of *N2 (Gad) for all TM1-B2N2/G catalysts.

The previous studies found that there were six PCET steps to accomplish the electrochemical reduction of N2 and produced two NH3 molecules[49-51]. Given that N2 is a thermodynamic much stable molecule (dissociation energy of N≡N as 945 kJ·mol-1), the first PCET plays a vital role among these steps because if the free energy change for this step is more positive, the following reaction would be more difficult to happen. Depending on N2 adsorption behavior, the reaction proceeds through three pathways, as summarized in Supplementary Figure 1, i.e., the so-called distal or alternating mechanism (end-on) and enzymatic mechanism (side-on)[48]. When taking the former way, the first PCET step is *N2 + H+ + e-*NNH, while in the enzymatic mechanism, the first PCET step is *N2 + H+ + e-*N-NH. Taking all the possibilities into account, we calculated and summarized the free energy changes for this step in Figure 3A and Supplementary Table 5. Compared to other catalysts, the hotspots map illustrates that all V1-B2N2/G and Cr1-B2N2/G lay on red areas, suggesting both of the catalysts have lower energies for first PCET and N2 adsorption with side-on pattern more energy-favorable than end-on. Among the three coordination environments of each TM1, V1-B2N2/G-1 and Cr1-B2N2/G-1 had lower energies of the first PCET, 0.06 and 0.41 eV, respectively.

High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction

Figure 3. Catalytic studies of NRR. (A) Heatmap of the free energy change for the first PCET step; (B) Heatmap of the lowest NRR limiting potential (UL) for all TM1-B2N2/G catalysts; (C) Optimal reaction mechanism for all TM1-B2N2/G catalysts; (D) Diagram of the NRR free energy change for the Cr1-B2N2/G-1 catalyst; (E) Selectivity of NRR and HER for all TM1-B2N2/G catalysts; (F) The proposed "five-step" strategy for high-throughput screening of NRR catalysts. NRR: N2 reduction reaction.

Subsequently, we simulated the whole NRR process for all TM1-B2N2/G with three mechanisms, including 72 reaction paths and 576 steps, which were provided in Supplementary Figures 2-4. And 72 limiting potentials [UL, Supplementary Table 6] of NRR can be sorted out and used as metrics to screen the activities or analyze the reaction mechanisms for all catalysts. Thus, the optimal activity and preferential reaction pathway for each catalyst are given in Figure 3B and C, with the UL values listed in Supplementary Table 7. The comparison of UL with the energies of the first PCET in Supplementary Figure 5 clearly demonstrates that the PDS of NRR on near 90 % TM1-B2N2/G is the first PCET, which provides a convenient strategy for the preliminary screening of the catalysts. For all catalysts, the early transitional metal single atoms (e.g., Sc, V, Cr) generally have lower UL than late ones (e.g., Fe, Co, Ni) regardless of the coordination manners. Herein, Cr1-B2N2/G SACs exhibit superior activities to other samples with UL as -0.43, -0.54, and -0.47 V for Cr1-B2N2/G-1, Cr1-B2N2/G-2, and Cr1-B2N2/G-3, respectively. It is worth mentioning that although some TM1-B2N2/G had lower first PCET energies (e.g., Cr1-B2N2/G-1 with 0.41 eV, V1-B2N2/G-1 with 0.59 eV), the UL values were uplifted due to more stable intermediates of *NH-NH2 and *N, which were minority cases but had to be taken into account when using first PCET energy to guide catalyst designing. To get more insights into the catalytic mechanism, Figure 3C screens out the optimal NRR reaction pathways for all catalysts, and it is noted that the enzymatic mechanism triggered by side-on N2 adsorption prevails compared to distal and alternating pathway, caused by modest adsorption energy for N2 activation. However, when N2 binds with the metal single atom strongly (e.g., Ti, V), the hydrogenation of *NH-NH in the enzymatic mechanism will become highly endothermic with an energy consumption of at least 1 eV, then the reaction will prefer to follow the distal or altering pathway. When N2 binds with metal single atom weakly (e.g. Fe, Co), the hydrogenation of the *NH-NH2 intermediate will become highly endothermic, making the reaction will tend to follow the distal or alternate pathway. In particular, NRR processes catalyzed by more active Cr1-B2N2/G catalysts are found to favorably take enzymatic way in Figure 3D and Supplementary Figure 6, which are triggered by first hydrogenations of N2 with lower energies via side-on adsorption mode.

To complete the catalytic cycle, the derived *NH3 species need to be desorbed to restore the catalysts. The free energy changes of second *NH3 desorption (no PCET steps involved) are shown in Supplementary Figure 7. It can be seen that *NH3 desorption requires more energy input on the early transitional metals (e.g., Sc, Ti, V with 0.9~1.1 eV) than the late ones (e.g., Co, Ni with 0.5~0.7 eV). However, Cr1-B2N2/G catalysts have moderate desorption energies to ensure easy desorption of NH3 and lead to the higher reactivity of NRR [Figure 3B]. For instance, the energy on Cr1-B2N2/G-1 is 0.78 eV, lower than V1-B2N2/G-1 (0.95 eV). In addition, the protonation as NH4+ in the acidic electrolyte reaction conditions can further facilitate the NH3 desorption[52,53].

From the above discussions, we can see that Cr1-B2N2/G shows a superior activity towards NRR to other catalysts, and herein Cr1-B2N2/G-1 performs the best, which are contributed by appropriate N2 adsorption energy and adsorption mode that guarantee lower energy inputs requiring for N2 activation, PDS (i.e., first PCET or the hydrogenation of *NH-NH2) and even the desorption of second *NH3.

The hydrogen evolution reaction is the key competing reaction during the NRR process. The presence of *H is facile to block TM single sites, which is deleterious for the NRR selectivity of the catalysts. To study the selectivity of catalysts, we compared the UL of the NRR and HER processes in Figure 3E. When UL (NRR) > UL (HER), i.e., the region above the dashed line, the electrochemical process is dominated by NRR. On the contrary, when UL (HER) > UL (NRR), i.e., the region below the dashed line, the electrochemical process is dominated by HER. According to this criterion, there are six TM1-B2N2/G-1 or 3 (TM = Sc, V, Cr, Fe, Co, Ni) expected to exhibit higher selectivities of NRR over HER, including the most active Cr1-B2N2/G-1.

The descriptor for NRR activity of TM1-B2N2/G

The high-throughput calculations of the TM1-B2N2/G catalysts after a 5-step screening strategy [Figure 3F] allow to establish a volcano plot of the reaction activity [i.e, optimal UL, Figure 4A] as the active centers of 3d TM1 varying from the left to right but with same coordination environments. Cr1-B2N2/G catalysts exhibit better NRR activities regardless of the coordination of Cr1. In order to guide the rational design of NRR catalysts with high performance, we chose the valence electron numbers of TM1 in the catalysts as the descriptor δ that is highly associated with N2 adsorption/activation and could be easily characterized by spectroscopy methods. Here δ is defined as the number of valence electrons of the transition metal. And it is the sum of the d-electron number of the transition metal and the Bader charge that the metal receives from the substrate. As shown in Figure 4B, a similar volcano relationship is constructed between NRR activities and δ for TM1-B2N2/G. The effective catalysts for NRR can be found at δ = 4.5 to 6.0, e.g., Cr1-B2N2/G. Although the estimation of δ descriptor by computation or readily characterization allows predicting the NRR catalysts with outstanding reactivity, the understanding of the above-disputed structure-activity connection remains unambiguous.

High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction

Figure 4. Descriptors of NRR activity. (A) Relationship between d-electron number of metal and NRR activity; (B) Relationship between valence electron number of metal in TM1-B2N2/G and NRR activity. NRR: N2 reduction reaction.

The above comprehensive comparisons in terms of UL and reaction pathways disclose that the NRR activity of TM1-B2N2/G is determined by N2 activation and hydrogenation of key intermediates (e.g., *NH-NH or *NH-NH2), so-called PDS, which are apparently with respect to the interaction between N2 and TM1. Additionally, it is believed that the electronic structure of metallic sites in TM1-B2N2/G is of significant importance in investigating N2 adsorption and activation. Therefore, we chose the most active Cr1-B2N2/G as a probe to interpret the influence of TM1 chemical nature and its corresponding steric coordination on the NRR activity.

Understanding the charge transfer during N2 activation

As shown in the electrostatic potential diagrams [Figure 5A-C] of three Cr1-B2N2/G, isolated Cr sites have positive potentials (blue zone), while B sites show negative potentials (red zone) and N sites nearly remain electrically neutral (green zone), suggesting that when Cr doped in the lattice of B2N2/G, the charge transfer mainly occurred between Cr and B atoms rather than N atoms. Similar phenomena of electron migration can be seen from the charge density difference diagrams of Cr1-B2N2/G [Supplementary Figure 8A-C]. The yellow region represents electron accumulation and the cyan region represents electron depletion. It is obvious that the charges are dominantly located at B sites and only a tiny amount of charges gather at N sites. Correspondingly, the electrons are accumulated at Cr sites. Therefore, the discrete charge distributions point that the interactions between Cr and B/N result in B atoms primarily donor electrons to Cr atoms as evidenced by the charge density distribution [Supplementary Figure 9].

High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction

Figure 5. Electronic structure of Cr1-B2N2/G catalysts. (A-C) The electrostatic potential (EP) of Cr1-B2N2/G catalysts; (D-F) The partial density of states (PDOS) of Cr1-B2N2/G catalysts; (G-I) The electron localization function (ELF) of Cr1-B2N2/G catalysts.

Next, the partial density of states [PDOS, Figure 5D-F] reveals that the interactions between Cr and B/N atoms were achieved via the couplings of d(Cr)-p(B/N) orbitals. The interplay of d(Cr) with p(B) occurs throughout Fermi energy level (EF), whereas the overlaps between d(Cr) and p(N) orbitals are only happened below EF, proving that the Cr-B binding is more intensive than Cr-N, matching with the above analysis of charge transfer process. Continuously, the electron localization function (ELF) was applied to check the bonding characteristics of Cr with B/N; meanwhile, the electron distributions between B and N were also inspected to survey the influence of coordination (i.e., cis and trans) on the chemical nature of Cr. Here the ELF was illustrated by contour plots in real space with values from 0 to 1. The number 1 means the complete localization of electrons, implying covalent bonds may form between Cr and B/N, while 0 means complete delocalization of electrons, implying the formation of ionic bonds. In particular, 0.5, marked as a line, corresponds to the free electron gas, implying metallic bonds may form considering the semimetal properties of boron[54,55]. The calculated ELF values were summarized in Supplementary Table 8 and mapped in Figure 5G-I. It is noted that a similar ELF (Cr-N) of 0.74~0.78 are found for three Cr1-B2N2/G, indicating the electrons are mainly delocalized in between instead of any elemental site. In contrast, the localization of electrons for Cr-B binding is different from each other catalyst, evidenced by three distinct ELF of 0.387, 0.308, 0.645 for Cr1-B2N2/G-1, Cr1-B2N2/G-2, and Cr1-B2N2/G-3 respectively. Such a deviation can be interpreted by the electronic stereo interferences between B···B (or N) site affected by the distances, e.g., 2.90, 2.01, > 3.8 Å for B···B in three catalysts with ELF computed as 0.037, 0.814 and not available for trans sample [Supplementary Table 8], validating that the electron localization associated with the interplay of B···B (or N) follows the order of Cr1-B2N2/G-2 > Cr1-B2N2/G-3 > Cr1-B2N2/G-1, whereas these electrons would not contribute to the charge transfer between Cr and B, resulting in the distinguished Bader charges of Cr as 1.156, 0.862 and 0.943 for three Cr1-B2N2/G catalysts [Supplementary Table 9]. The higher charge densities are expected to deliver better performance of NRR activity for Cr1-B2N2/G-1.

With the electronic structures of Cr1-B2N2/G catalysts in mind, we continued to study the interaction of N2 with the Cr single atom, which was of paramount importance for the following PCET processes and NRR activity. It was uncovered that the N2 molecule preferential adsorbed on Cr with side-on mode [Figure 3]. The prolonged N≡N bond lengths of ~1.21Å (versus 1.11 Å in vacuum) verified the successful activation of N2 rather than physical adsorption [Supplementary Table 2 and Supplementary Table 3]. Figure 6A-C visualizes the charge transfer in N2 activation process on three Cr1-B2N2/G, showing that the charge is migrated from Cr to N.

High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction

Figure 6. Investigation of the interaction between N2 and Cr1-B2N2/G. (A-C) The charge density difference of N2 adsorbed on different Cr1-B2N2/G catalysts. The isosurface value was set to be 0.002 e·Å-3; (D-F) The PDOS of the d-orbitals of Cr and the p-orbitals of N2 in the Cr1-B2N2/G catalysts for N2 adsorption; (G-I) The partial Crystal Orbital Hamiltonian Population (pCOHP) of N≡N bonds in Cr1-B2N2/G catalysts adsorbed with N2 in Side-on mode.

Subsequently, the orbital interactions related to N2 activation were analyzed by PDOS [Figure 6D-F]. Generally, the interaction of N2 with Cr was achieved by the d(Cr)-π*(N2) orbital couplings. As illustrated in Figure 6B, the better spatial symmetry and the excellent compatibility in terms of energy for these two orbitals favors a stronger interaction between N2 and Cr1-B2N2/G-1 at 2.21 eV and 4.35 eV compared to the cases of Cr1-B2N2/G-2 and Cr1-B2N2/G-3, which may be contributed by the more charge densities of Cr in Cr1-B2N2/G-1 [Supplementary Table 9]. The above stronger synergistic effect also effectively eliminates the spin polarization effect of the active site and ensures a higher NRR activity, as reported in the literature[56]. In addition, the intensive interplay of N2 with Cr1-B2N2/G-1 can also be proved by the energy level splitting of d(Cr) orbitals that originally dominated at 2.56 eV [Figure 5D].

To have a quantitative comparison of N2 activation on three Cr1-B2N2/G catalysts, the binding strength of N≡N bonds in the activated N2 molecule was examined by the partial Crystal Orbital Hamiltonian Population (pCOHP)[57]. As depicted in Figure 6G-I, the bonding (-pCOHP > 0) and antibonding (-pCOHP < 0) contributions are shown on the right and left of one panel for each catalyst, respectively. The bonding contributions below the fermi level follow the order of Cr1-B2N2/G-1 (-1.65) < Cr1-B2N2/G-3 (-2.12) < Cr1-B2N2/G-2 (-2.18), agreeing well with the integral COHP values, which are -1.24, -1.35, and -1.32, respectively. Here more negative ICOHP values indicate stronger N≡N bonds. On the contrary, the weaker binding strength of N≡N corresponds to a higher level of N2 activation, as evidenced by the Bader charges of activated N≡N as 0.568, 0.413, 0.458 [Supplementary Table 10], which thus favors the subsequent hydrogenation reaction. This is consistent with the limiting potential order of Cr1-B2N2/G-1 (-0.43 V) < Cr1-B2N2/G-3 (-0.47 V) < Cr1-B2N2/G-2 (-0.58 V) in Supplementary Table 7.

In this regard, the structure-activity relationship revealed for Cr1-B2N2/G can be generalized to other transitional metals. That is to say, the NRR activities of TM1-B2N2/G are mainly governed by the charge transfers between N2 and metal atoms, i.e., more electrons lend from the metal to N2, and higher NRR reactivity is obtained for the catalysts. Therefore, the descriptor of δ is further correlated with the Bader charge of activated N≡N [Supplementary Figures 10-15] for the most active catalyst in each group of TM1-B2N2/G. Interestingly, a volcano trend of N≡N Bader charge depended on the valence electron number (δ) of the transitional metal single atom can be observed in Supplementary Figure 16. Such a connection reveals that the physical nature of TM1 intrinsically dominates the interaction between N2 and the active site, which determines N2 activation and the subsequent consecutive hydrogenation, even the desorption of *NH3. The obtained knowledge paves the way for the development of advanced NRR.


High-throughput DFT calculations were implemented to investigate the NRR performance of TM1-B2N2/G, which was distinguished by different transitional metals and coordination configurations. Among all catalysts, Cr1-B2N2/G-1, with better activity and selectivity, was found by a “five-step” screening strategy, which has an NRR limiting potential of -0.43 V. The volcano shape of the structure-activity relationship highlights the importance of the valence electron number of TM1 as a descriptor, which can accelerate the rational design process of advanced NRR catalysts. In addition, it was revealed that the limiting potentials were highly relevant with electron donation from TM1 to N2, which determined the N2 activation and potential-determining step. The scenario demonstrated in this work offers guidance for the synthesis of efficient NRR catalysts by experiment.


Authors’ contributions

Conceptualization, calculation, data curation, data analysis, and original draft: Cao N

Calculation, data organization, and data analysis: Zhang N

Data analysis: Wang K

Resources: Yan K

Conceptualization, resources, data analysis, review and writing: Xie P

Availability of data and materials

The data that support the findings of this study are available from the corresponding author upon reasonable request.

Financial support and sponsorship

This work was supported by the National Natural Science Foundation of China (22278365), the National Key Research and Development Program of China (2022YFE0128600), the Natural Science Foundation of Zhejiang Province (LR22B060002), and a grant from Shanxi-Zheda Institute of Advanced Materials and Chemical Engineering (2021ST-AT-002).

Conflicts of interest

All authors declared that there are no conflicts of interest.

Ethical approval and consent to participate

Not applicable.

Consent for publication

Not applicable.


© The Author(s) 2023.

Supplementary Materials


1. Comer BM, Fuentes P, Dimkpa CO, et al. Prospects and challenges for solar fertilizers. Joule 2019;3:1578-605.

2. Smil V. Detonator of the population explosion. Nature 1999;400:415-415.

3. Licht S, Cui B, Wang B, Li FF, Lau J, Liu S. Ammonia synthesis. Ammonia synthesis by N2 and steam electrolysis in molten hydroxide suspensions of nanoscale Fe2O3. Science 2014;345:637-40.

4. Qian J, An Q, Fortunelli A, Nielsen RJ, Goddard WA 3rd. Reaction mechanism and kinetics for ammonia synthesis on the Fe(111) surface. J Am Chem Soc 2018;140:6288-97.

5. Guo W, Zhang K, Liang Z, Zou R, Xu Q. Electrochemical nitrogen fixation and utilization: theories, advanced catalyst materials and system design. Chem Soc Rev 2019;48:5658-716.

6. Suryanto BHR, Du H, Wang D, Chen J, Simonov AN, Macfarlane DR. Challenges and prospects in the catalysis of electroreduction of nitrogen to ammonia. Nat Catal 2019;2:290-6.

7. Niu L, Wang D, Xu K, et al. Tuning the performance of nitrogen reduction reaction by balancing the reactivity of N2 and the desorption of NH3. Nano Res 2021;14:4093-9.

8. Pang Y, Su C, Jia G, Xu L, Shao Z. Emerging two-dimensional nanomaterials for electrochemical nitrogen reduction. Chem Soc Rev 2021;50:12744-87.

9. Yang Z, Wang J, Wang J, et al. 2D WO3-x nanosheet with rich oxygen vacancies for efficient visible-light-driven photocatalytic nitrogen fixation. Langmuir 2022;38:1178-87.

10. Giuffredi G, Asset T, Liu Y, Atanassov P, Di Fonzo F. Transition metal chalcogenides as a versatile and tunable platform for catalytic CO2 and N2 electroreduction. ACS Mater Au 2021;1:6-36.

11. Wu T, Melander MM, Honkala K. Coadsorption of NRR and HER intermediates determines the performance of Ru-N4 toward electrocatalytic N2 reduction. ACS Catal 2022;12:2505-12.

12. Majumder M, Saini H, Dědek I, et al. Rational design of graphene derivatives for electrochemical reduction of nitrogen to ammonia. ACS Nano 2021;15:17275-98.

13. Lee CH, Pahari S, Sitapure N, Barteau MA, Kwon JS. DFT–kMC analysis for identifying novel bimetallic electrocatalysts for enhanced NRR Performance by suppressing her at ambient conditions via active-site separation. ACS Catal 2022;12:15609-17.

14. Wang Y, Su H, He Y, et al. Advanced electrocatalysts with single-metal-atom active sites. Chem Rev 2020;120:12217-314.

15. Li X, Liu L, Ren X, Gao J, Huang Y, Liu B. Microenvironment modulation of single-atom catalysts and their roles in electrochemical energy conversion. Sci Adv 2020:6.

16. Li L, Chang X, Lin X, Zhao ZJ, Gong J. Theoretical insights into single-atom catalysts. Chem Soc Rev 2020;49:8156-78.

17. Hu S, Tian N, Li M, et al. Trapezohedral platinum nanocrystals with high-index facets for high-performance hydrazine electrooxidation. Chem Synth 2023;3:4.

18. Kaiser SK, Chen Z, Faust Akl D, Mitchell S, Pérez-Ramírez J. Single-atom catalysts across the periodic table. Chem Rev 2020;120:11703-809.

19. Kim J, Choi S, Cho J, Kim SY, Jang HW. Toward multicomponent single-atom catalysis for efficient electrochemical energy conversion. ACS Mater Au 2022;2:1-20.

20. Novoselov KS, Geim AK, Morozov SV, et al. Electric field effect in atomically thin carbon films. Science 2004;306:666-9.

21. Stoller MD, Park S, Zhu Y, An J, Ruoff RS. Graphene-based ultracapacitors. Nano Lett 2008;8:3498-502.

22. Tan C, Cao X, Wu XJ, et al. Recent advances in ultrathin two-dimensional nanomaterials. Chem Rev 2017;117:6225-331.

23. Zhuo HY, Zhang X, Liang JX, Yu Q, Xiao H, Li J. Theoretical understandings of graphene-based metal single-atom catalysts: stability and catalytic performance. Chem Rev 2020;120:12315-41.

24. Liu P, Nørskov JK. Ligand and ensemble effects in adsorption on alloy surfaces. Phys Chem Chem Phys 2001;3:3814-8.

25. Li X, Shen P, Luo Y, et al. PdFe single-atom alloy metallene for N2 electroreduction. Angew Chem Int Ed Engl 2022;61:e202205923.

26. Ren G, Shi M, Liu S, Li Z, Zhang Z, Meng X. Molecular-level insight into photocatalytic reduction of N2 over Ruthenium single atom modified TiO2 by electronic Metal-support interaction. J Chem Eng 2023;454:140158.

27. Zhang A, Liang YX, Zhang H, Geng ZG, Zeng J. Doping regulation in transition metal compounds for electrocatalysis. Chem Soc Rev 2021; 50:9817-9844.

28. Vineesh TV, Kumar MP, Takahashi C, et al. Bifunctional electrocatalytic activity of boron-doped graphene derived from boron carbide. Adv Energy Mater 2015;5:1500658.

29. Yu X, Han P, Wei Z, et al. Boron-doped graphene for electrocatalytic N2 reduction. Joule 2018;2:1610-22.

30. Kresse G, Furthmüller J. Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set. Comput Mater Sci 1996;6:15-50.

31. Kresse G, Furthmüller J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys Rev B Condens Matter 1996;54:11169-86.

32. Blöchl PE. Projector augmented-wave method. Phys Rev B Condens Matter 1994;50:17953-79.

33. Perdew JP, Chevary JA, Vosko SH, et al. Atoms, molecules, solids, and surfaces: Applications of the generalized gradient approximation for exchange and correlation. Phys Rev B Condens Matter 1992;46:6671-87.

34. Perdew JP, Wang Y. Accurate and simple analytic representation of the electron-gas correlation energy. Phys Rev B Condens Matter 1992;45:13244-9.

35. Yang WJ, Ren JN, Zhang HW, et al. Single-atom iron as a promising low-temperature catalyst for selective catalytic reduction of NOx with NH3: a theoretical prediction. Fuel 2021;302:121041.

36. Yang W, Zhao M, Ding X, et al. The effect of coordination environment on the kinetic and thermodynamic stability of single-atom iron catalysts. Phys Chem Chem Phys 2020;22:3983-9.

37. Yang W, Xu S, Ma K, et al. Geometric structures, electronic characteristics, stabilities, catalytic activities, and descriptors of graphene-based single-atom catalysts. NMS 2020;2:120-31.

38. 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:154104.

39. Gray CM, Saravanan K, Wang G, Keith JA. Quantifying solvation energies at solid/liquid interfaces using continuum solvation methods. Mol Simul 2017;43:420-7.

40. Nelson R, Ertural C, George J, Deringer VL, Hautier G, Dronskowski R. LOBSTER: Local orbital projections, atomic charges, and chemical-bonding analysis from projector-augmented-wave-based density-functional theory. J Comput Chem 2020;41:1931-40.

41. Nørskov JK, Bligaard T, Logadottir A, et al. Trends in the exchange current for hydrogen evolution. J Electrochem Soc 2005;152:J23.

42. Wang V, Xu N, Liu J, Tang G, Geng W. VASPKIT: A user-friendly interface facilitating high-throughput computing and analysis using VASP code. Comput Phys Commun 2021;267:108033.

43. Shi L, Bi S, Qi Y, et al. Anchoring mo single-atom sites on B/N codoped porous carbon nanotubes for electrochemical reduction of N2 to NH3. ACS Catal 2022;12:7655-63.

44. Fan M, Chen X, Zhang M, Cui L, He X, Zou X. Highly dispersed Ru nanoclusters anchored on B,N co-doped carbon nanotubes for water splitting. Inorg Chem Front 2022;9:968-76.

45. Niu H, Wang X, Shao C, Zhang Z, Guo Y. Computational screening single-atom catalysts supported on g-CN for N2 reduction: high activity and selectivity. ACS Sustainable Chem Eng 2020;8:13749-58.

46. Choi C, Back S, Kim N, Lim J, Kim Y, Jung Y. Suppression of hydrogen evolution reaction in electrochemical N2 reduction using single-atom catalysts: a computational guideline. ACS Catal 2018;8:7517-25.

47. Lv X, Wei W, Huang B, Dai Y, Frauenheim T. High-throughput screening of synergistic transition metal dual-atom catalysts for efficient nitrogen fixation. Nano Lett 2021;21:1871-8.

48. Shi L, Yin Y, Wang S, Sun H. Rational catalyst design for N2 reduction under ambient conditions: strategies toward enhanced conversion efficiency. ACS Catal 2020;10:6870-99.

49. Jiao D, Liu Y, Cai Q, Zhao J. Coordination tunes the activity and selectivity of the nitrogen reduction reaction on single-atom iron catalysts: a computational study. J Mater Chem A 2021;9:1240-51.

50. Liu S, Cheng Z, Liu Y, et al. Boosting electrochemical nitrogen reduction reaction performance of two-dimensional Mo porphyrin monolayers via turning the coordination environment. Phys Chem Chem Phys 2021;23:4178-86.

51. Jin H, Guo C, Liu X, et al. Emerging two-dimensional nanomaterials for electrocatalysis. Chem Rev 2018;118:6337-408.

52. Katsounaros I, Figueiredo MC, Chen X, Calle-vallejo F, Koper MTM. Structure- and coverage-sensitive mechanism of NO reduction on platinum electrodes. ACS Catal 2017;7:4660-7.

53. Chun H, Apaja V, Clayborne A, Honkala K, Greeley J. Atomistic insights into nitrogen-cycle electrochemistry: a combined DFT and kinetic monte carlo analysis of NO electrochemical reduction on Pt(100). ACS Catal 2017;7:3869-82.

54. Ma F, Jiao Y, Gao G, et al. Graphene-like two-dimensional ionic boron with double dirac cones at ambient condition. Nano Lett 2016;16:3022-8.

55. Zhang H, Li Y, Hou J, Tu K, Chen Z. FeB6 Monolayers: The graphene-like material with hypercoordinate transition metal. J Am Chem Soc 2016;138:5644-51.

56. Sun X, Zhang H, Wang S, et al. Tuning the coordination microenvironment to boost the electrocatalytic HER activity of M3(C6O3S3)2. J Phys Chem C 2022;126:16606-14.

57. Deringer VL, Tchougréeff AL, Dronskowski R. Crystal orbital Hamilton population (COHP) analysis as projected from plane-wave basis sets. J Phys Chem A 2011;115:5461-6.

Cite This Article

Export citation file: BibTeX | RIS

OAE Style

Cao N, Zhang N, Wang K, Yan K, Xie P. High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction. Chem Synth 2023;3:23.

AMA Style

Cao N, Zhang N, Wang K, Yan K, Xie P. High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction. Chemical Synthesis. 2023; 3(3): 23.

Chicago/Turabian Style

Ning Cao, Nan Zhang, Ke Wang, Keping Yan, Pengfei Xie. 2023. "High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction" Chemical Synthesis. 3, no.3: 23.

ACS Style

Cao, N.; Zhang N.; Wang K.; Yan K.; Xie P. High-throughput screening of B/N-doped graphene supported single-atom catalysts for nitrogen reduction reaction. Chem. Synth. 2023, 3, 23.

About This Article

© The Author(s) 2023. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License (, which permits unrestricted use, sharing, adaptation, distribution and reproduction in any medium or format, for any purpose, even commercially, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

Author Biographies

Ning Cao
Ning Cao received his B.Sc. and M.Sc. degrees from Northeastern Electric Power University in 2019 and 2022, respectively. He is pursuing his Ph.D. degree at the College of Chemical and Biological Engineering, Zhejiang University. His main research interests include computational design and reaction mechanism understanding of heterogeneous catalysis.
Nan Zhang
Nan Zhang received his B.Sc. and M.Sc. degrees from Northeastern Electric Power University in 2019 and 2022, respectively. He is pursuing his Ph.D. degree at the School of Chemistry and Chemical Engineering, Nanjing University. His main research interests include enzyme catalysis and nitrogen conversion.
Ke Wang
Ke Wang received his Ph.D. degree from the College of Chemistry, Fuzhou University. Currently, he is a postdoctor at the College of Chemical and Biological Engineering, Zhejiang University. His research interests focus on developing advanced catalysts for thermal/photo/electro catalysis of CO2 reduction and hydrogen evolution reactions.
Keping Yan
Keping Yan received his B.Sc. and M.Sc. degrees from Beijing Institute of Technology in 1983 and 1986, respectively. He received his Ph.D. degree from Eindhoven University of Technology in 2001. He was a professor at the College of Environmental and Resource Sciences, Zhejiang University from 2006 to 2009. From April 2009 to the present, he is a professor at the College of Chemical and Biological Engineering, Zhejiang University. His research interests include electric dust removal and low temperature plasma for pollutant removal.
Pengfei Xie
Pengfei Xie received his Ph.D. degree in Chemistry from Fudan University in 2016. He worked as a Postdoctoral Researcher in the Chemical and Biomolecular engineering department, at Johns Hopkins University, then promoted to Assistant Research Scientist in 2020. Xie moved to Zhejiang University in 2021 and is currently a professor at the College of Chemical and Biological Engineering. His research interests include atomically dispersed metal catalysts and their applications in the field of energy and environmental science.

Data & Comments




Comments must be written in English. Spam, offensive content, impersonation, and private information will not be permitted. If any comment is reported and identified as inappropriate content by OAE staff, the comment will be removed without notice. If you have any queries or need any help, please contact us at

Download PDF
Cite This Article 6 clicks
Like This Article 7 likes
Share This Article
Scan the QR code for reading!
See Updates
Chemical Synthesis
ISSN 2769-5247 (Online)


All published articles are preserved here permanently:


All published articles are preserved here permanently: