Ammonia, as the basic chemical substance for the production of various synthetic chemicals such as chemical fertilizers, dyes and resins, plays an indispensable role in agriculture and industry [1-3]. The traditional preparation of ammonia is mainly based on the Haber-Bosch process, which operates under extremely harsh conditions (400–600 ℃ and 150–300 atm) [4-6], leading to the huge energy consumption and the massive emission of carbon dioxide [7,8]. To overcome such shortage of the Haber-Bosch method, the electrocatalysis of the N₂-to-NH₃ conversion paves a green and sustainable avenue for the ammonia synthesis under the mild reaction [9-12], wherein the performance is significantly depended on the efficiency of the catalyst used. Despite tremendous efforts in recent decades, electrocatalysts that selectively and efficiently reduce nitrogen to ammonia remain elusive. Therefore, the development of the highly-active materials for the ammonia synthesis has attracted great attention [13].

Among the varied catalysts, single-atom catalysts (SACs) with utmost atom-utilization efficiency and distinctive coordination environment have become a research hotspot toward the conversion from N₂ into NH₃ [14-16]. Such as, Zhang et al. demonstrated that FeN₃ embedded graphene significantly promotes the electrocatalysis of nitrogen reduction reaction (NRR) wherein the rate determining step is the first protonation with the free energy change of 0.04 eV [17]. Similarly, Zhao et al. found that WN₃ decorated graphene is a promising NRR electrocatalyst with the limiting potentials of 0.33 V [18]. As the graphene derivative, the VN₃ embedded graphane shows good NRR performance via the enzymatic pathway with the free energy change of 0.28 eV [19]. Additionally, W anchored g-C₃N₄ exhibits the good activity toward NRR with a limiting potential of 0.35 V via associative enzymatic pathway [20]. The mentioned results indicate the fact that the carbon-based SACs with transition metal decoration as the active site delivers good activity toward N₂-to-NH₃ conversion.

Exception with the common carbon-based material mentioned above, the waves of renewed interest have turned to two-dimensional monoelemental structures due to its interesting physical and chemical properties [21,22]. It seems the group-VA elemental monolayer is attractive as an electrocatalyst for ammonia synthesis. As an illustration, Liu et al. found that W embedded black phosphorus have the onset potential of 0.42 V [23]. Furthermore, Xu et al. found the V doped arsenene provides an ultralow overpotential of 0.10 V via the enzymatic pathway [24]. Analogously, Song et al. demonstrated that NbP₃ embedded arsenene has a low overpotential of 0.36 V via the distal mechanism [25]. Besides, Gu et al. discovered that anionic Os and Re metal centers on the defective antimonene can electrochemically catalyze the nitrogen reduction reaction with a limiting potential close to that of stepped Ru(0001) [26]. The similar conclusion has been drawn by Cao et al. that the antimonene increases the efficiency of NH₃ production [27]. For the last group-VA elemental monolayer, bismuthene, it has been successfully synthesized [28,29], which possesses identical atomic structure in comparison with phosphorene, arsenene and antimonene. More encouragingly, Hao et al. discovered that the presence of bismuth nanocrystals and potassium cations in water significantly promotes nitrogen electroreduction to ammonia [30]. Based on the mentioned information, bismuthene might be a promising material for the N₂-to-NH₃ conversion. Therefore, we are naturally interested in whether the pristine bismuthene is active and whether transition metal decoration would boost the efficiency.

To address the mentioned issues, we systematically studied the catalytic activities of bismuthene with and without anchoring transition metal atom in the application of N₂-to-NH₃ electrocatalysis by means of the density functional theory (DFT) calculations. For simplification, it donates as TM@Bis wherein 3d/4d/5d TM are taken into consideration with the exception of Tc and Lu due to toxic or radioactive characteristic. Our results demonstrate that the pristine bismuthene fails to capture the N₂ molecule, indicating its infeasibility to boost the N₂-to-NH₃ conversion. Furthermore, the strategy of transition metal decoration is a valid method to activate the inert bismuthene, wherein that W@Bis exhibits the best activity under the distal pathway with a small limiting potential U_L of 0.26 V. More interestingly, a volcano curve is established between U_L and TM d valence electrons wherein W with 4 electrons in d band has the minimum U_L meanwhile either increase or decrease the number of valence electrons would enlarge the U_L data. Such phenomenon originates from the underlying acceptance-back donation mechanism. Our study provides the theoretical prediction that the functional bismuthene is a feasible to improve the efficiency for the nitrogen electroreduction to ammonia.

All calculations are performed within the DFT framework as implemented in DMol³ code [31,32]. The generalized gradient approximation with the Perdew−Burke−Ernzerhof (PBE) functional is employed to describe exchange and correlation effects [33]. The DFT Semi-core Pseudopots (DSPP) core treat method is implemented for relativistic effects, which replaces core electrons by a single effective potential and introduces some degree of relativistic correction into the core [34]. The double numerical atomic orbital augmented by a polarization function (DNP) is chosen as the basis set [31]. A smearing of 0.005 Ha (1 Ha = 27.21 eV) to the orbital occupation is applied to achieve accurate electronic convergence. In the geometry structural optimization, the convergence tolerances of energy, maximum force and displacement are 1.0 × 10⁻⁵ Ha, 0.002 Ha/Å and 0.005 Å, respectively. The spin-unrestricted method is used for all calculations. A conductor-like screening model (COSMO) was used to simulate water solvent environment for all calculations [35]. COSMO is a continuum model in which the solute molecule forms a cavity within the dielectric continuum. The DMol³/COSMO method has been generalized to periodic boundary cases. The dielectric constant is set as 78.54 for H₂O. As a representative, in order to verify the effect of water solution, the adsorption abilities of W@Bis and Nb@Bis with and without COSMO correction are tabulated in Table S1 (Supporting information) for comparison and the corresponding free energy profiles are presented in Figs. S1 and S2 (Supporting information). In line with the previous reports [36], such correction enhances the adsorption strengths and reduces the thermodynamic barriers, being benefiting for the reactants capture and protonation progress. Furthermore, a first principles molecular dynamics (AIMD) simulation was performed to assess the thermodynamic stability at 300 K. A canonical ensemble (NVT) within Massive GGM thermostat scheme was used to control the temperature [37,38]. During the stimulation processes, a total simulation time was set to 3 ps with the time step of 1 fs. To check the sufficiency of the stimulation time adopted herein, we further perform another AIMD calculation on W@Bis with the stimulation time of 5 ps and the corresponding result is given in Fig. S3 (Supporting information). Consistently, no structural damage is observed for W@Bis, which is in line with the results drawn from the stimulation time of 3 ps. Therefore, we analyze our data on the basis of AIMD with 3 ps in the discussion.

The 15 Å-thick vacuum is added to avoid the artificial interactions between the catalyst and its images. The adsorption energies (E_ads) of NRR intermediates are calculated by

(1)

where E_system, E_catalyst and E_m represent the total energy of the adsorption system, the catalyst and the adsorbates, respectively.

To evaluate the thermodynamic stabilities of the designed TM embedded bismuthene, we calculated the binding energies E_b using the following equation, which has been extensively used in literatures [39]:

(2)

where E_V is the total energy of the catalyst without TM atoms and E_TM is the total energy of the single TM atom.

To evaluate the energetic feasibility of the experimental synthesis, we calculated the formation energies E_form using the following equation [40]:

(3)

where E_TM-bulk is the total energy of the TM bulk.

Gibbs free energy change (ΔG) of the elementary steps was constructed based on the computational hydrogen electrode (CHE) model developed by Nørskov et al., where the chemical potential of (H⁺ + e⁻) pair in solution is equal to the half of the chemical potential of a gas-phase H₂ [41-43]. The corresponding ΔG is determined by the following equation [44]:

(4)

where ∆E is the energy difference, ∆ZPE indicates the change in zero point energy, and ∆S represents the change in the entropy, respectively. Herein, the temperature T is 298.15 K. The ZPE and S are calculated from the vibrational frequencies according to standard methods. Following the suggestion of Lim and Wilcox [45], in order to reduce the cumbersome calculation, the substrates are fully constrained for the frequency calculation. The term of ∆G_U = – eU corresponds to the free energy contribution caused by the variation of electrode potential while ΔG_pH is the pH correction of the free energy, being expressed as ΔG_pH = – kTln10 × pH, where k_B is the Boltzmann constant. In our case, pH is zero under the acid solution. ΔG < 0 corresponds to an exothermic adsorption process vice versa. The methods have been successfully applied for analyzing N₂ fixation process, indicating its robustness and feasibility [46,47].

The limiting potential U_L is determined by the following equation [48]:

(5)

where ΔG_PDS is the free energy change of the elementary step with the maximum endothermic character during the N₂-to-NH₃ conversion.

Fig. 1a describes the atomic configuration of the transition metal atom decorated bismuthene (TM@Bis) with one TM substitution. The binding energies E_b between TM and its surrounding are listed in Table 1. In the same period, the TM element with half-filling or full-filling d band delivers the local minimum E_b due to its stable electron configuration. Therein, the corresponding values are −2.93 and −3.01 eV for Cr and Cu, −3.51 and −2.57 eV for Mo and Ag, −3.96 and −3.25 eV for Re and Au, respectively. Since E_b > −1.5 eV indicates the intensive chemical interaction [17], it reveals that TM atom is firmly fixed in the bismuthene. In order to further reveal the strong interaction, Fig. 1b presents the partial density of states (PDOS) between TM d band and Bi p band of TM@Bis wherein the elements of Ti, V, Nb, Mo, W are taken as the illustration along with the information of pristine bismuthene. Wherein, the intensive pd coupling between TM decoration and its Bi ligand as well as the descended p band of Bi ligand clearly supports the formation of the chemical TM-Bi bonds. Furthermore, the TM dopant enriches the electron occupation at Fermi energy, being beneficial for the electron transfer, another advantage for the application of electrocatalysis [49,50]. In addition, the formation energies E_form of TM@Bis are further evaluated in order to illustrate the synthesis feasibility and the values are tabulated in Table S2 (Supporting information). Exception with the dopants of Sc, Y, Pd, Ag and Au, the rest ones possess the positive E_form values. It indicates that the formation of TM@Bis is an endothermic process and needs the extra energy provided. However, the similar systems, such as Ru/g-C₃N₄, Pt/g-C₃N₄, Fe/N-G and Co/N-G, have been experimentally prepared and their E_form values are 3.12, 2.88, 2.26 and 2.27 eV, respectively [40,51-54]. To a certain extent, the comparable value in the same order of magnitude implies that the synthesis of TM@Bis is feasible from the energetic aspect. Therefore, it is reasonable to believe its realization due to the advanced experimental technology.

Before the activity investigation, the adsorption of N₂ molecule is considered and the corresponding free energies G_ads(N₂) are listed in Table 1. The pristine bismuthene is unable to capture the N₂ molecule due to the endothermic G_ads(N₂), indicating it is inactive toward ammonia production. Furthermore, the similar situations are found for Sc@Bis, Cr@Bis, Mn@Bis, Cu@Bis, Y@Bis, Zr@Bis, Ag@Bis, Pd@Bis, Pt@Bis and Au@Bis, respectively. Therefore, the mentioned systems are ruled out for further consideration. In addition to the good affinity toward N₂ molecule, the relative weaker H adsorption is necessary to address the critical competition of the undesirable hydrogen evolution reaction (HER) in aqueous solutions [55-57]. Herein, the difference of adsorption free energies between N₂ molecule and H is adopted for the first step of quick screening, that is, △△G_ads = G_ads(N₂) – G_ads(H), wherein the negative value indicates the preferred NRR selectivity over HER, visa verse. Fig. 1c reveals that the NRR efficiency of Hf@Bis, Ta@Bis and Ir@Bis would be suppressed by the competition of side HER process. However, the rest configurations would offer good NRR selectivity, involving Ti@Bis, V@Bis, Fe@Bis, Co@Bis, Ni@Bis, Nb@Bis, Mo@Bis, Ru@Bis, Rh@Bis, W@Bis, Re@Bis and Os@Bis. Subsequently, we focus our attention on the free energy changes of the first protonation and the sixth protonation since the potential-determining step (PDS) is generally located at either the first protonation (*N₂ + H⁺ + e^– → *NNH) or the sixth protonation (*NH₂ + H⁺ + e^– → *NH₃) [47,58]. Fig. 1d demonstrates that the elements of Ti, V, Nb, Mo, W and Os would be promising dopants due to the possibility of low thermodynamic barriers if artificially setting the values less than 0.55 eV as our reference according to the previous reports [43,59,60].

To confirm the efficiency of the N₂-to-NH₃ conversion, we perform the analysis of the whole reaction pathways and the free energy changes of elementary steps are given in Table 2. There are three feasible reaction pathways involving distal, alternative and enzymatic pathways [61-63]. For the distal mechanism, three proton-electron pairs (H⁺ + e⁻) first react with the outermost N atom to form the first NH₃ molecule, and with that, the remaining N atom undergoes subsequent hydrogenation until the release of second NH₃ molecule. While for the alternating mechanism, the proton-electron pairs alternatively attack the two N atoms, generating the two NH₃ molecules one after another. The reduction of the side-on adsorbed N₂ proceeds through the enzymatic mechanism, in which two N atoms are hydrogenated by proton-electron pairs alternately, and two NH₃ molecules will be released consecutively. Herein, it is noteworthy that the continuous protonation progress is from the adsorbed N-containing intermediates reacted with the protons in the solution [47,64-67]. It is a universal and valid method for the proton treatment adopted by the calculations. Therefore, we do not consider the co-adsorption of NRR intermediates and *H specie.

As a representative, Fig. 2 shows the free energy profiles of W@Bis via distal mechanism, alternating mechanism and enzymatic mechanism, respectively. For continuous protonation via the distal mechanism, the free energy changes of the elementary steps ΔG are 0.16, –0.32, –1.1, –0.29, 0.25 and 0.26 eV, respectively. The maximum endothermic step is located at the last protonation step of *NH₂ + H⁺ + e^– → *NH₃. Therefore, the PDS is identified as the *NH₃ formation and the corresponding thermodynamic barrier ΔG_PDS is 0.26 eV. Similarly, ΔG of the protonation steps are 0.16, 0.86, –0.59, 0.54, –2.26 and 0.26 eV under the alternating pathway wherein the PDS is the second protonation step of *NNH + H⁺ + e^– → *NHNH and ΔG_PDS is 0.86 eV. Furthermore, the ΔG are 0.24, 0.35, −0.61, 0.46, −2.12 and 0.26 eV along enzymatic mechanism, respectively. Since the step of *NH*NH₂ + H⁺ + e^– → *NH₂*NH₂ encounters the greatest energy barriers of 0.46 eV, the fourth protonation step is regarded as the PDS of enzymatic mechanism. Herein, we notice that the *NH₃ desorption is upshifted energetically and it needs to overcome a barrier of 2.27 eV. In accord with the previous works, the endothermic characteristic of NH₃ desorption is a universal phenomenon. To remove its adsorption and recover the active site, the formed *NH₃ could be further hydrogenated to NH₄⁺ and readily leave the electrode with the aid of strongly acidic aqueous solution [64,68-70]. Therefore, with aid of the solution [71,72], the *NH₃ desorption is not a problematic obstacle in NRR. As a result, distal mechanism with lowest thermodynamic barrier is more favorable for NRR electrocatalysis. In the regard, the N₂-to-NH₃ convention catalyzed by W@Bis material is energetically favored under the limiting potential U_L of 0.26 V. As a supplementary, the NRR activity of W adsorption is further investigated since surface decoration via adsorption is another well-established method to provide the active site for electrocatalysis. According to the free energy profiles presented in Fig. S4 (Supporting information), the limiting potential U_L is 0.51 V, being much higher than 0.26 V offered by the W substitution. Therefore, W substitution is more valid to boost ammonia synthesis.

The free energy profiles of the rest systems are supplemented in Figs. S5–S15 (Supporting information). According to the same analysis, we summarize the ΔG of the protonation steps in Table 2 and present the U_L data of varied TM@Bis in Fig. 3a. It is noteworthy that we do not consider the enzymatic mechanism on Ru@Bis and Rh@Bis due to the unfeasible characteristic of *N₂ side-on adsorption [73]. Herein, the NRR catalyzed by Mo, Nb, W and Re dopants prefer the distal mechanism and the U_L values are 0.41, 0.35, 0.26 and 0.55 V, respectively. Furthermore, the NRR catalyzed by Ni, Ru and Rh sites are through alternating pathway with the large U_L of 1.09, 1.20 and 1.31 V, respectively, indicating the poor activities. Additionally, the NRR catalyzed by Ti, V, Fe, Co and Os dopants prefer enzymatic mechanism with the U_L of 0.43, 0.37, 0.81, 0.74 and 0.55 V, respectively. According to the U_L data, the NRR efficiency is in the order of W@Bis > Nb@Bis > V@Bis > Mo@Bis > Ti@Bis > Re@Bis = Os@Bis > Co@Bis > Fe@Bis > Ni@Bis > Ru@Bis > Rh@Bis. As a supplementary, the U_L comparison (or reduction potential) between NRR and HER are listed in Table S3 (Supporting information). The data further confirms the good efficiency of the W@Bis endowed with a significantly smaller U_L of NRR vs HER. Furthermore, the dopants of Re and Os are also promising candidates for experimental synthesis due to its good selectivity, in despite of the slightly higher U_L of NRR. On the contrary, the NRR activity of TM@Bis with the rest dopant might be suffered from the undesirable side reaction. Besides the activity, Table 2 reveals that the PDS under the energetically preferred pathway is the last protonation for W@Bis mentioned above and it turns to the first protonation for the rest ones with the exception of Os@Bis located at the second protonation. On the basis of Sabatier principle, it implies that the NRR activity of the former slightly suffers from the strong affinity toward NRR intermediates meanwhile the latter provides the insufficient activation ability. Fig. 3b presents the universal linear relation of the adsorption energies E_ads among varied NRR intermediates, being stemmed from the similar adsorption configuration [74,75]. Herein, the adsorption affinity of TM@Bis could be presented by the *N adsorption energy E_ads(*N). Furthermore, Fig. 3c reveals the monotonic scaling between U_L and E_ads(*N) wherein the TM@Bis endowed with the stronger E_ads(*N) offers better performance for N₂ electroreduction into NH₃. Therefore, our results are in accord with our conjecture that the thermodynamic limitation is posed by the correlations among the adsorption of varied NRR intermediates, which is prevalent phenomenon for SAC with TM active site [27,76].

According to the mentioned information, W@Bis performs best toward ammonia generation herein. As previous reported, the W dopant embedded different substrates have been discovered to efficiently catalyze nitrogen fixation. For example, the U_L values are 0.33, 0.35 and 0.42 V for the W supported on the N doped graphene, g-C₃N₄ and black phosphorus, respectively [18,20,23]. It seems a similar good performance of W dopant as the active site. In order to reveal the underlying origin of its universality mentioned above, Fig. 4a schematically shows the integral partial density of states of N₂ adsorption on W@Bis as an illustration. Wherein, the electrons contributed by the d-band of W@Bis disturb the N≡N π bond, which makes the N₂ bond longer and activated. Spontaneously, the electron returns to the d band again. It is well-known as the acceptance-back donation mechanism. It implies that the metal atom in the middle range of Periodic Table would deliver relative better NRR performance, because its partial filling d band has the ability to release electrons and host electrons. On the contrary, the nearly empty d band of early transition metal element limits the electron transfer from the active site to the reactant meanwhile the almost full d band of the late transition metal element suffers from the electron transfer from the reactant to the active site. The mechanism is in line with our data, that is, the W dopant in the middle range offer relative smaller U_L in comparison with the early transition metal elements such as Ti and late transition metal elements such as Fe, Co and Ni. To confirm this view, the fitting between the number of valence electrons and U_L data is performed in Fig. 4b. There is a volcano curve with the W@Bis suited at the summit. That is, W with 4 electrons in d band have the minimum limit potentials of 0.26 V. Either increase or decrease the number of valence electrons would enlarge the U_L data. The final question is why Nb@Bis delivers slightly inferior performance in comparison with W@Bis since W and Nb have the identical valence electrons in d band. As discussed above, the NRR performances of Nb@Bis and W@Bis separately suffer from the thermodynamic barriers of the first and last protonation steps, which originate from the diffident adsorption ability. The reduced reactivity of Nb@Bis is clearly revealed by the PDOS between W/Nb d band and N p band shown in Fig. 4c wherein the pd hybridization of Nb@Bis is relative weaker in comparison with W@Bis, being in line with the trend of the affinity toward NRR intermediates, as reflected by the E_ads listed in Table S1 (Supporting information). It is the reason why the activities of Nb@Bis and W@Bis are different as mentioned above. Furthermore, to gain deep insight into the effect of W@Bis during continuous protonation progress, the electron transfers are analyzed in detail wherein three moieties are artificially divided and displayed in Fig. 4d, involving the adsorbates N_xH_y (moiety 1), the TMBi₃ site (moiety 2) and the bismuthene support (moiety 3). Fig. 4e illustrates the Mulliken charges during the protonation progress, indicating that good conductivity of W@Bis enables the effective electron transfers. More specifically, the electron communications are occurred between the bismuthene support and NRR intermediates meanwhile the electron fluctuation of moiety 2 is negligible [77]. That is, moiety 3 is an electron reservoir and moiety 2 acts as an electron transmitter [19,78].

Before closing our discussion, we pay our attention on the stability of W@Bis, another important factor for catalyst design. Considering the possibility of the structural deformation caused by the adsorbates, Table S4 (Supporting information) presents the variation of W-Bi bond lengths caused by the adsorption. Clearly, the adsorption leads to the elongation of W-Bi bonds. However, the changes of average bond length are within 0.2 Å and the W-Bi bonds keep binding under the adsorption. Furthermore, the PDOS presented in Fig. S16 (Supporting information) provides the evidence to support the strong fixation of W atom wherein the overlapping of pd orbitals is consistently observed below the Fermi energy level, indicating the maintenance of covalent W-Bi bonds, in line with the small variation of bond lengths. Therefore, the W@Bis configuration will not be destroyed by the adsorbates. Moreover, we perform a first principles molecular dynamics simulation with a constant temperature of T = 300 K, in order to ensure the thermodynamic stability of W@Bis at the room temperature. Fig. 4f demonstrates that the atomic structure is well preserved during the simulation, as reflected by the monitoring of the bond length d_W-Bi between W and its surrounding Bi atoms. Despite the observed fluctuation, no breaking of W-Bi bonds is found, indicating the good resistance against the structural collapse. However, a slight deformation of atomic configuration of Nb@Bis is observed in Fig. S17 (Supporting information), indicating its inferior stability referring to W@Bis. Therefore, W@Bis endowed with its high efficiency and good stability is a promising candidate for ammonia synthesis.

In a summary, the catalytic performance of TM@Bis toward N₂ electroreduction to NH₃ is systematically evaluated by means of density functional theory calculations. Based on the free energy profiles, the efficiency is in the order of W@Bis > Nb@Bis > V@Bis > Mo@Bis > Ti@Bis > Re@Bis = Os@Bis > Co@Bis > Fe@Bis > Ni@Bis > Ru@Bis > Rh@Bis, wherein the thermodynamic barriers are 0.26, 0.35, 0.37, 0.41 0.43, 0.55, 0.55, 0.74, 0.81, 1.09, 1.2 and 1.31 eV, respectively. The PDS under the energetically preferred pathway is the last protonation for W@Bis and it turns to the first protonation for the rest ones with the exception of Os@Bis located at the second protonation. It indicates that the NRR performance of the former slightly suffers from the strong affinity toward NRR intermediates meanwhile the latter provides the insufficient activation ability. In addition, the best performance offered by the W dopant stems from the underlying acceptance-back donation mechanism. Therefore, our work offers promising candidates for the experimental synthesis and provides a fundament understanding for the material design for nitrogen reduction electrocatalysis. It represents an efficient and sustainable ammonia synthesis from nitrogen with broad scientific impacts.

Declaration of competing interest

There are no conflicts to declare.

Acknowledgments

The authors greatly acknowledge the financial support from the National Natural Science Foundation of China (Nos. 21503097, 52130101, 51701152, 21806023 and 51702345) and China Scholarship Council (No. 202008320215).

Supplementary materials

The Supporting information is available on line, including details of adsorption energy, formation energy, free energy profile, limiting potential and bond length. Supplementary material associated with this article can be found, in the online version, at doi:10.1016/j.cclet.2022.07.002.