Air-stable phosphorus-doped molybdenum nitride for enhanced electrocatalytic hydrogen evolution

Molybdenum-based electrocatalysts for hydrogen evolution have been investigated extensively in recent years. However, unlike other non-oxides, molybdenum nitride generally shows a weak preference for hydrogen evolution and low performance owing to surface oxidation and the strong Mo–H bond. Here, we prepare an air-stable molybdenum nitride through a multi-step solid-state reaction. We find that a uniformly dispersed mixture of the precursors is optimal for preparation of the electrocatalyst. To further enhance hydrogen evolution performance towards practical device applications, phosphorus doping is carried out, using a few layered black phosphorus source. The phosphorus-doped molybdenum nitride (P–Mo–N) sample catalyzes hydrogen evolution with potentials of 105, 145, and 157 mV at the current densities of 10, 50, and 100 mA/cm2, respectively, in 0.5 M H2SO4 solution with a small Tafel slope of 43 mV/dec. Thus it outperforms many of the state-of-art molybdenum-based hydrogen evolution catalysts reported to date. Molybdenum-based electrocatalysts are widely investigated as water splitting catalysts but they suffer from low efficiency due to surface oxidation. Here the authors report the synthesis of air-stable molybdenum nitride which displays enhanced hydrogen evolution activity.

H ydrogen evolution reaction (HER) from electrocatalytic water splitting has been regarded as one of the important strategies to solve the world's increasing energy and environment issues [1][2][3][4][5][6][7][8][9][10] . Water electrolysis requires a stable electrocatalyst which has the ability to overcome high overpotential and provide efficient activity. To date, Pt-group materials are still the most commonly used catalysts in the commercial HER and there is a demand for a necessary replacement due to their limited supply and high cost [11][12][13][14][15][16][17][18][19][20] . For the sustainable and large-scale application of HER, earth-abundant materials have been under intense investigation in an effort to reach a suitable Pt-like HER performance [21][22][23][24][25][26][27][28][29][30] . For the typical HER process, two accepted reaction mechanisms are described, namely the Volmer-Heyrovsky route and Volmer-Tafel route [31][32][33][34] . In the Volmer reaction, the formation of M-H takes place by proton adsorption on the surface sites of the catalyst. In the subsequent Heyrovsky reaction, one M-H undergoes proton-coupled electron transfer reaction to release one H 2 . Whereas, for the Tafel step, two nearby adsorbed *H (H-M) combine into one H 2 . Based on the abundant research in recent years, a good Pt-like HER catalyst in acidic solution should exhibit suitable proton adsorption ability and a fast carrier migration [21][22][23][24][25][26][27][28][29][30] .
Mo-based materials are widely explored for improved HER performance [21][22][23][24][25][35][36][37][38][39][40] . In particular, undoped Mo 2 N has a metallic electronic structure and an excellent Pt/Pd-like hydrogenation ability. Undoped Mo 2 N are frequently used in HER, but exhibits moderate performance in comparison with Mo 2 C or MoP [41][42][43] . This is caused by the high electronegativity of N species, which can induce a strong Mo-H bond and then weaken the H 2 generation rate 40 . To improve the electrocatalytic performance, doping strategies have been widely used [44][45][46][47][48][49][50][51] . The most significant advantage of doping has been to obtain a suitable M-H bond with moderate strength 44,46,50,51 . Moreover, N doping is reported to improve the electrocatalytic hydrogen performance of MoP by reducing the electrical resistance, increasing the number and quality of active sites 52 . The non-metal constituents, nitrogen (N) and phosphorous (P) belong to the same group family and thus have a similar electronic structure. However, phosphorus is more electronegative (2.19) than nitrogen (3.04), and thus, the replacement of N by P may lead to a weakened M-H bond 50,51 . However, thus far, P-doped molybdenum nitride for improving HER performance has not been studied in detail.
One key challenge for obtaining a high HER activity is to locate the electroactive sites of catalysts. The edge sites exposed doped elements or vacancies are reported to be responsible for the enhanced reaction performance [50][51][52][53][54][55] . Increasing the number of reaction sites has been one major explored strategy. For the synthesis of Mo 2 N, ammonia reduction of molybdenum oxide is widely used through high-temperature calcination, which could induce the expected aggregation, thus reducing the exposed reaction sites and weakening the corresponding electrocatalytic activity. However, the releasing of the in situ generated gas can promote the formation of holes or disordered lattice structures, thereby avoiding the aggregation and increasing active sites 53,55 .
Herein, we report a two-step synthesis strategy to prepare Pdoped Mo 2 N. The mixture of MoO 3 and ethanediamine are first treated by high-strength solid-phase grinding and then calcinated in an N 2 atmosphere. Mo 2 N is obtained with abundant exposed grain-boundary sites (Mo-N). Few layered black phosphorus (BP) is used for P doping via vacuum stirring of the BP and Mo 2 N mixture, followed by calcination to obtain P-Mo 2 N (P-Mo-N). The obtained samples of Mo-N and P-Mo-N display notable electrocatalytic HER performance and stability, even at a high current density of 160 mA/cm 2 .

Results
Synthesis of Mo-N and P-Mo-N. First, Mo 2 N (Mo 2 N-r) was synthesized from the reduction reaction of MoO 3 by the release of in situ NH 3 similar to previous work conducted using melamine as the NH 3 source 39 . Solid-state reaction has been used previously to generate samples with numerous holes and vacancies due to the in situ release of gases 53,55 . Owing to its active character, Mo 2 N can be oxidized easily at ambient conditions. Supplementary Figure 1a shows the photograph of the autoignition of Mo 2 N-r. Such an active material owns so many surface oxidation states that its catalytic activity is relatively low [41][42][43] . More exposed Mo sites may be responsible for the above active phenomenon of autoignition. To obtain a relatively stable sample, ethanediamine was used to generate in situ NH 3 . This was also expected to change the crystal faces via the release of gases during its decomposition. The layered MoO 3 and ethanediamine (mass ratio: 1:2) were mixed thoroughly by a ball grinder (Fig. 1a). The mixture was then calcined in an N 2 atmosphere. After that, a layered graphite-like sample (Mo-N-2) without autoignition characteristics was obtained as shown in Supplementary Figure 1b. The crystal structure was confirmed to be Mo 2 N from the X-ray diffraction (XRD) measurement ( Supplementary Figure 1c). The direct appearance of the two samples was tested by scanning electron microscopy (SEM). As shown in Supplementary Figure 2, Mo 2 N-r and Mo-N-2 give the aggregated and layered structures, respectively. The samples from the two different processes were measured by transmission electron microscope (TEM) to further check their morphologies. Mo 2 N-r and Mo-N-2 both show the layered structure; however, the ethanediamine-induced sample presents some obvious particles (Supplementary Figure 3). The high resolution TEM in Supplementary Figure 4 shows the detailed difference of the two samples. Mo 2 N-r shows the typical continuous lattice structure while Mo-N-2 gives obvious crossed lattice fringes, that is, for the Mo-N-2 sample, the formed nanoparticles combine together to maintain the apparent nanosheet morphology (Supplementary Figure 2). This suggests that the ethanediamine-based solid-state reaction can induce the dismemberment of layered MoO 3 to form nanoparticles, and then generate more grain-boundaries under the releasing of gases. To get the suitable grain-boundary and prevent aggregation, the amount of ethanediamine was tested. The mass ratios of MoO 3 to ethanediamine were increased from 1:3 to 1:4 to 1:5. Supplementary Figure 5 shows the corresponding XRD results of obtained samples (Mo-N-n). Only the Mo 2 N phase can be found in three samples indicating no noticeable change to the crystal lattice structure even at the highest ratios (Supplementary Figure 5a). The TEM images show obvious grain-boundaries, as well as some detectable lattice distortion (Supplementary Figure 5b-d). Please note that all the synthesized samples of Mo-N-x present the appearance of nanosheet, which is composed of nanoparticles.
The electrocatalytic hydrogen generation from water splitting in 0.5 M H 2 SO 4 solution was carried out for detecting the activity of the obtained samples. The Mo-N-n samples give an enhanced HER performance obviously compared with Mo 2 N-r (Supplementary Figure 6), suggesting that the ethanediamine-driven Mo 2 N samples have a dominant Mo-H adsorption-desorption performance. More grain-boundary sites, i.e. more nanoparticle formation from the sufficient reduction reaction between ethanediamine and MoO 3 is the internal origin 53,55 . With increasing ethanediamine content, the HER performance increases accordingly. For the samples obtained from MoO 3 and ethanediamine mass ratio: 1:4 and 1:5 (Mo-N-4 and Mo-N-5) show a similar activity. Therefore, Mo-N-4 (Mo-N) was selected for further research.
At the second step, P doping was carried out to further modulate their HER performance as show in Fig. 1b-d. The Mo-N sample and layered black phosphorus (BP) 56,57 were used as the precursors. Few layered BP is not stable in air and can be oxidized to PO x species, which can be reduced and doped into Mo 2 N lattices. Moreover, few-layer BPs still keep the structural integrity of Mo 2 N during the P doping process. Supplementary Figure 7 shows the typical TEM images of BP, the layered structure can be detected. The mixture of the Mo-N sample and BP was collected after about 3 h vacuum-stirring with the BP being loaded onto the surface of Mo-N (Fig. 1c). After storage in air for 1 day, the sample was treated by calcination with the adding of ethanediamine in an N 2 atmosphere to finalize the P doping (Fig. 1d). Figure 1 also shows the corresponding TEM images of the Mo-N samples. Clearly, the image of the P doped sample (Supplementary Figure 8) shows a disordered lattice structure, because of the large ionic radius of P. To confirm the negligible impact of secondary calcination on the HER performance, the reference sample without BP addition was synthesized. As given in Supplementary Figure 9, the XRD pattern does not show an obvious change compared with the initial one (Supplementary Figure 9a). The HER activity also remains unchanged (Supplementary Figure 9b). The obtained sample was named as P-Mo-N for the following discussion, in which the P content was 10 wt.% (Please note that the P content of P-Mo-N samples was detected by X-ray photoelectron spectroscopy (XPS) measurements and is semi-quantitative). The specific surface areas of the samples were further measured by N 2 adsorption-desorption method, calculated by multiple spot Brunauer-Emmett-Teller (BET) theory. As shown in Supplementary Figure 10, the BET surface areas of Mo 2 N-r, Mo-N, and P-Mo-N are 4.48, 5.54, and 5.38 m 2 /g, respectively, further confirming that the ethanediamine-driven method can avoid the aggregation.
Characterization of P-Mo-N. Figure 2a shows the XRD pattern of P-Mo-N sample, four typical peaks are located at 2 theta of 37°, 41°, 62°, and 74.5°, which can be attributed to (112), (200), (220), and (312) of tetragonal Mo 2 N (PDF: 25-1368). Compared with the initial Mo-N sample, no obvious change of peak position and intensity suggests that the secondary calcination of P doping hardly changes the crystalline structure. From the SEM images in  Fig. 2c. Clearly, the nanoparticles can be found along the nanosheet edge. The inset shows the corresponding electron diffraction spectrum (EDS), revealing the typical four crystal faces, which is consistent with the XRD pattern (Fig. 2a). We also presented the TEM images of the reference P-Mo-N samples with 9.5 wt.% and 10.5 wt.% P doping. As shown in Supplementary Figure 12, TEM images similar to that of the 10 wt.% sample P-Mo-N can be found, with both of them giving the rough sheet surface. The corresponding EDS results suggest that the Mo 2 N phase was obtained and the P content in this range hardly affects the final product formation. The HRTEM image of the P-Mo-N sample in Fig. 2d further suggests the existence of nanoparticles. The obvious grain boundaries can be found and two interplanar crystal spacing values are measured to be 0.24 and 0.21 nm, which can be assigned to the (112) and (200) crystal faces, respectively. The disordered lattice fringes can also be detected. Supplementary Figure 13 gives the corresponding HRTEM images of the reference P doped samples with similar results. Figure  However, the PO x /P-Mo 2 N sample, which is not a completely reduced sample, shows a different result compared with the initial sample. At the binding energy of ca. 230.0 and 232.6 eV, two weak peaks, which can be assigned to the formation of MoP x , are found, suggesting that a proportion of the used P has been doped into the Mo 2 N lattice. For the N 1 s XPS spectra, the similar results compared with the patterns in Supplementary Figure 16b can be found. However, the N 1 s binding energy of PO x /P-Mo 2 N sample is ca. 298.28 eV, slightly smaller than that of PO x /Mo 2 N (398.3 eV) and higher than that of P-Mo-N (398.26 eV), further confirming the successful P doping into the Mo 2 N lattice. For the P 2p XPS spectrum in Supplementary Figure 16c, the sample of PO x /Mo 2 N exhibits a strong peak located at 133.7 eV, which comes from the surface of PO x species. Meanwhile, for PO x / P-Mo 2 N, the peak centered at 133.7 eV is weakened and two new peaks (130.0 and 131.1 eV) appear, suggesting the doped P element exists at the low valance state. Based on the above analysis, Supplementary Figure 16d gives the models of the P doping process: at the first stage, the PO x species is only adsorbed onto the surface of Mo 2 N; during calcination, the P element of PO x will be reduced into the low valance P and doped into the lattice of Mo 2 N; after the calcination (2 h), all the P was reduced and in situ incorporated into the Mo 2 N phase. As discussed earlier, to further confirm the present samples are not formed by simple physical aggregation, the Mo-N and P-Mo-N were treated by 12 h ultrasound. After the centrifugation, the supernate was chosen for the TEM measurement. As shown in Supplementary Figure 17, the two samples both show unchanged morphology. The crystal boundaries can be still found, suggesting that the different nanoparticles are fastened together with each other strongly, not a physical mixture. During the synthesis process, the ethanediamine and MoO 3 were grinded thoroughly and then a uniform mixture was obtained. Under the calcination, ethanediamine released NH 3 and other gases, which reacted with MoO 3 to form Mo 2 N. During the reaction, the layered structure of MoO 3 was destroyed with the formation of Mo 2 N nanoparticles. The calcination further combined the particles with many crystal boundaries with the disordered crystal lattices obtained. To the best of our knowledge, it is the first time to report one stable Mo 2 N sample in air, and the exposed different crystal surfaces may be responsible for the above novel phenomena i.e. air-stable and high HER performance.
HER performance of Mo-N and P-Mo-N. Most of electrocatalytic HER among Mo-based samples reported to date investigate their carbides or phosphates in acidic media. Owing to the relative strong electronegativity of N, the molybdenum nitride samples exhibited a low activity 37,39 . Herein, the capability of HER electrocatalysis of P-Mo-N sample was explored in 0.5 M H 2 SO 4 solution at ambient temperature. The samples were prepared into a homogeneous ink with Nafion solution and isopropanol and then dropped onto the glassy carbon electrodes for the HER measurements. The potentials obtained from the measurements were converted to values relative to a reversible hydrogen electrode (RHE). Figure 3a shows the relevant polarization curves of different electrodes, including P-Mo-N and the reference samples of Mo 2 N-r and Mo-N. The purchased 20% Pt/ C electrode was also used as a standard sample. Obviously, the activity trend is as follows: Mo 2 N-r < Mo-N < P-Mo-N < Pt/C, with the onset potentials (defined as the overpotential at current density of −1 mA cm −2 ) of 74, 73, 28, and 17 mV, respectively. Moreover, compared to the reference Mo 2 N-r and Mo-N, the P doped sample presents a sharply enhanced current response at low applied potentials, suggesting that the P doping can boost a better proton transfer and adsorption-desorption processes. Then the potentials at four current density values of 10, 50, 100, and 150 mA/cm 2 were compared as shown in Fig. 3b. Mo 2 N-r and Mo-N both present relatively high potentials, larger than 200 mV. For the P doped Mo 2 N, a sharp decrease can be detected in comparison to the two undoped samples. In details, the P-Mo-N sample gives the potentials of 105, 145, 157, and 164 mV at 10, 50, 100, and 150 mA/cm 2 , respectively. As a contrast, the Pt/C sample exhibits a potential of 47 mV at the 10 mA/cm 2 condition, similar to the reported values in the literature 8,42,44 . Figure 3c gives the corresponding Tafel slope results of the samples under study. The referenced commercial Pt/C sample gives the smallest slope of 37 mV/dec, followed by P-Mo-N of 43 mV/dec, Mo-N of 76 mV/dec, and Mo 2 N-r of 71 mV/dec. A slope larger than 30 mV/dec for the P-Mo-N sample suggests a Volmer-Heyrovsky route for HER, 10,11,14 . The exchange current density (j 0 ), the most inherent measure of HER activity, was also calculated based on the Tafel equations (Supplementary Table 2). The j 0 values of P-Mo-N, Mo-N and Mo 2 N-r are 0.02, 0.0052, and 6.67 × 10 −6 mA/cm 2 , respectively. The P doping can facilitate a large exchange current density, increasing by a factor of 4 in comparison with the undoped samples. Supplementary Table 3 gives the direct activity comparison of P-Mo-N with the Mobased samples recently reported in the literature. Obviously, our P-Mo-N sample gave the best HER performance in 0.5 M H 2 SO 4 solution in terms of current density and Tafel slope. To further confirm the reaction sites to be the P sites on the surface of Mo 2 N, the MoP and N doped MoP were also synthesized 52 .
Supplementary Figure 18a gave the XRD patterns of MoP (PDF: 24-0771) and Supplementary Figure 18b showed the corresponding HER performance. Clearly, although Mo 2 N shows the lowest HER activity, the P doped Mo 2 N gives the best.
The electrochemical impedance spectroscopy (EIS) was carried out to probe charge transfer processes on the catalysts in 0.5 M H 2 SO 4 solution. Figure 3d shows the Nyquist plots of the samples of Mo 2 N-r, Mo-N and P-Mo-N under a 250 mV potential condition. Clearly, the P doped sample gives the smallest circular arc, followed by Mo-N, Mo 2 N-r, suggesting that the P-Mo-N has a relative small carrier migration resistance. Supplementary Figure 19 shows the Nyquist plot of reference Pt/C sample and the corresponding fitted circuit diagram. The R Ω and R ct mean the electrode resistance and carrier transfer resistance across the catalyst/solution interface, respectively. The Warburg and CPE are diffusion resistance and the constant phase element, both are out of consideration owing to the limited applying frequency in our present work. Supplementary Table 4 shows the corresponding fitted results of R Ω and R ct of the three samples. Typically, all the samples give a similar R Ω value (0.15 ohm), suggesting the resistance of catalyst/electrode is almost the same. Regarding R ct , the three samples of Mo 2 N-r, Mo-N and P-Mo-N give decreased values from 111.2, to 50.02 to 14.69 ohm, respectively, confirming that the P doping can induce a suitable proton adsorption-desorption, and then decrease the resistance of the solid/solution interface. To get more information about the catalyst/solution interface, double-layer capacitance (C dl ) calculations were carried out using cyclic voltammetry (CV) to roughly evaluate their effective electrochemical active surface area on Mo 2 N-r, Mo-N and P-Mo-N. The CV was tested with a potential range of −0.2 to 0 V vs RHE at different scan rates, and the Δj at −0.1 V vs scan rate was plotted. The slope was supposed to be twice that of the C dl 44,60 . Our results (Supplementary Figure 20) reveal a considerably larger C dl of the P-Mo-N sample (275 mF cm −2 ) compared with bare Mo-N (115 mF cm −2 ) and the referenced Mo 2 N-r (60 mF cm −2 ), indicating more accessible electroactive sites created by P doping in the Mo 2 N network. The P doping was achieved by BP loading first and then calcination, thus more P element was confined into the surface lattice, which may contribute to the small carrier migration resistance across the catalyst/solution interface and a large C dl .
Moreover, the operating stability of the catalyst is of great concern in the real industrial application. The polarization curve after 10,000 cycles was performed on P-Mo-N catalyst in the 0.5 M H 2 SO 4 solution. As given in Fig. 3e, negligible change can be found before and after the continuous cycle test. The chronoamperometry (j-t) response was further carried out under three potentials of 145, 157, and 164 mV. The j responses with the corresponding values of 50, 100, and 150 mA/cm 2 give no obvious decrease or increase during the continuous 12 h measurement (Fig. 3f). The EIS spectrum and XRD pattern of the P-Mo-N catalyst were further tested as given in Supplementary Figure 21. No change to the resistance and phase can be found, further suggesting a good stability of P-Mo-N catalyst. Supplementary Figure 22 gives the contrastive HRTEM images with no change. Moreover, the electrocatalytic HER activity after 3 months (storing in a sample tube without vacuum treatment) was also conducted. Almost no shift of polarization curves after one and three months (Supplementary Figure 23). The effect of P content for HER performance was verified. As given in Supplementary Figure 24a, the 9.5 wt%, 10 wt%, and 10.5 wt% sample show the similar XRD pattern, confirming that the P content in these three values does not change the phase of Mo 2 N. The HER activity in Supplementary Figure 24b confirms that the 10 wt.% P doping can induce the best electrocatalytic performance.
COMMUNICATIONS CHEMISTRY | https://doi.org/10.1038/s42004-018-0097-9 ARTICLE COMMUNICATIONS CHEMISTRY | (2018) 1:95 | https://doi.org/10.1038/s42004-018-0097-9 | www.nature.com/commschem HER mechanism of P-Mo-N. It is found that the as-prepared Mo 2 N, P-Mo-N with 10 wt% P, catalyst owns numerous exposed crystal faces, which could reduce its activity in air and lower the Mo-H bonds. This is proven to effectively catalyze water splitting for H 2 production whose performance is superior over most of state-of-art catalysts reported in the literature 25,37,41 . Furthermore, the P element doping boosts the HER performance greatly and the enhancement should be clarified. In realistic water electrocatalysis, two vital parameters, i.e. interfacial property of catalyst/solution and the intrinsic nature of the catalyst, determine the final performance 18,21,37,61 . Thus, the interface properties were firstly studied to investigate the improvement shown by P doping. Four solutions with different pH values were chosen for this purpose.
As shown in Fig. 4a, the polarization curves of the P-Mo-N sample in the four varying pH solutions present different patterns. Typically, in the 0.5 M H 2 SO 4 media (pH~0), the catalyst exhibits the best electrocatalytic performance, followed by that in 0.05 M H 2 SO 4 (pH~1), 0.5 M Na 2 SO 4 (pH~6.8) and 0.1 mM NaOH (pH~10), suggesting that more free H + can help obtain a high H 2 evolution performance. Moreover, at the current density of 10 mA/cm 2 , we measured the corresponding potential values: 135, 128, and 100 mV for 0.1 mM NaOH, 0.5 M Na 2 SO 4 , and 0.05 M H 2 SO 4 , respectively. There is a double electrode layer (DEL) at the solid|liquid interface, this can be attributed to different H + concentrations will inducing a varied H + adsorption or arrangement 62 . The 0.1 mM NaOH condition for P-Mo-N HER presented in Fig. 4a demonstrated a better performance than bare Mo 2 N in 0.5 M H 2 SO 4 (Fig. 3a), suggesting the P doping can induce a favorable H + arrangement for HER in the DEL even in a low H + concentration solution. The corresponding EIS spectrum  Fig. 4b in different pH solutions shows that similar trends can be found. For studying the adsorption behavior of reaction species in the DEL, a bias potential of −100 mV (after this potential, the H 2 evolution will start 62 ) was used according to the polarization curves. Supplementary Table 5 shows the corresponding fitted resistance values based on the fitting circuit (Supplementary Figure 19b). Obviously, the R ct in the 0.5 M H 2 SO 4 solution gives the smallest value, and it increases with the elevated pH, suggesting that more free H + in the electrolyte will be suitable for DEL formation with a low resistance. The EIS spectrum of bare Mo 2 N catalyst was also conducted under the same condition with an R ct value of 62.3 ohm, further confirming that the P doping can induce a suitable H + adsorption and then a low resistance DEL. To get more information about the DEL, the electrocatalytic reaction area was further carried out according to the CV curves under the potential range of −0.2 to 0 V as shown in Supplementary Figure 25. The C dl values were calculated to be 275, 115, 80, and 29 mF/cm 2 for 0.5 M H 2 SO 4 , 0.05 M H 2 SO 4 , 0.5 M Na 2 SO 4 , and 0.1 mM NaOH, respectively, further validating the above results, that is, more free H + can be adsorbed by the P doped surface of Mo 2 N. The intrinsic nature of the P doped Mo 2 N was next investigated for understanding the HER enhancement. Density functional theory (DFT) calculations were carried out. Owing to the more exposed facets of the synthesized Mo 2 N sample, four main crystal face models, (112), (200), (220), and (312) with the N replaced by P were established as shown in Supplementary  Figure 26 and Fig. 5a-b. For studying the reaction centers, the elements for proton adsorption were also labeled (Supplementary Figure 26). The computational H adsorption free energy values on different reaction centers in different surfaces of the bare Mo 2 N and P doped Mo 2 N samples are summarized as given in Supplementary Table 6-7 and Fig. 5c. Typically, for the undoped sample, the (200) facet gives a good H adsorption result either on exposed N center or Mo center. However, the Mo atoms in (200) surface can be oxidized easily and then the bare Mo 2 N always gives a low HER performance (Fig. 3a). Exposing more facets with lots of grain boundaries can stabilize (200) facet and induce more reaction sites, and then the Mo 2 N sample with a grain boundary displays a better HER performance than the one with a single exposed surface. For the P doped sample, the P sites in two facets of (112) and (312) give suitable HER reaction sites with the relative small H adsorption free energy values of 0.04 and 0.08 eV, respectively. Figure 5c clearly shows that the hydrogen adsorption free energy (G(*H)) gives a very small G(*H) value for the P-Mo-N sample at P sites, whereas the Mo sites or N sites in doped or bare samples show the strong water adsorption ability and consequent a weak desorption behavior. Supplementary  Figure 27 and Fig. 5d, e give the contrastive N-p, P-p and Mo-d orbital density of states (DOS) in four facets of bare and doped Mo 2 N samples. Obviously, the p orbitals of P atom in all conditions show a right shift compared to N atoms, suggesting that the more anti-bonding orbital of the P atom will hybridize with the s orbital of the H atom, and then weaken the G(*H) value. Meanwhile, for the Mo atom, an opposite tendency can be found, further confirming that the P doping is good for the HER process. Moreover, the work functions (W f ) based on ultravioletphotoelectron spectra (UPS) in Fig. 5f and Supplementary Figure 28 were calculated to be 4.8 and 4.53 eV for Mo 2 N and P-Mo 2 N, respectively. Note that a small W f of a catalyst will be helpful for electrons migrating through the sample and then assist with a better activity. Here a low W f further confirms that the P doping can promote electron transfer and enhance chemical activity.

Discussion
We describe here an air-stable Mo 2 N electrocatalyst prepared from a violent solid-state reaction of ethanediamine and MoO 3 by utilizing the strategy discussed above. Owing to the several reaction sites, the Mo 2 N was achieved with several exposed crystal facets, which lower the active Mo atom number. Moreover, the by-product of abundant crystal boundaries was beneficial to the HER performance. BP was used for P doping and the obtained P-Mo-N showed an enhanced HER activity compared with bare Mo 2 N. Experimental and computational results collectively confirmed more reaction sites provided by P doping, a lower carriers/ protons migration resistance across the catalyst/solution interface, as well as a suitable hydrogen adsorption free energy. Such a method for ensuring uniform multi-site solid-state reaction can be extended to synthesis of other materials, and thus enrich earthabundant elements based catalysts with satisfactory activity. Moreover, the P doping strategy can help us understand more about electrocatalytic reaction mechanisms such as HER presented here and also design more doped samples to enhance active site exposure and provide more adsorption and desorption characteristics.

Methods
Synthesis of Mo-N. In a typical synthesis, 0.5 g MoO 3 and a certain amount of ethanediamine was first mixed and grinded by a quartz mortar, after that, the mixture was moved to a ball mill and further treated thoroughly for 45 min. Then, the uniform mixture was transferred into one vacuum tube furnace and calcinated under 500°C for 4 h in a flowing N 2 -gas atmosphere. After naturally cooling to room temperature, the product was collected and stored for further use. Please note that without the treatment of the ball mill, an active Mo 2 N will be generated.
Synthesis of P-Mo-N. In a typical synthesis, 0.3 g of Mo-N was added into 100 mL BP solution, the mixture was then placed under vacuum and then stirred for 3 h to achieve complete adsorption. After that, the product was treated by centrifugation and vacuum drying. The obtained powder was mixed with ethanediamine by one quartz mortar and then moved into the vacuum tube furnace.
After 4 h 500°C calcination in a flowing N 2 -gas atmosphere, the product was collected, dried at 80°C and stored for further use.
Material characterization. X-ray diffraction (XRD) patterns of the samples were recorded on a Rigaku Smartlab-9 kW instrument using Cu Kα X-ray (λ = 1.54186 Å) radiation at a scanning rate of 4°/min in the region of 2θ = 10-80°. X-ray photoelectron spectra (XPS) were acquired on a Kratos Axis Ultra DLD spectrometer with Al Kα (hυ = 1486.6 eV) as the excitation source. Transmission electron microscopy (TEM) test was performed on a FEI Tecnai G2 F20 electron microscope at an acceleration voltage of 200 kV. Spherical aberration-corrected electron microscope was carried out on a JEM-ARM200F electron microscope (Beijing Zhongkebaice Technology Service Co., Ltd.). The field emission scanning electron microscopy (FESEM) were observed under a SU8020 electron microscopy. Brunauer-Emmett-Teller (BET) surface area analysis was carried out on a Quadrasorb SI-3 equipment.
HER electrode preparation. Typically, 10 mg sample was added into a 50 μL isopropanol, 50 μL nafion (Nafion 117 solution, Sigma-Aldrich) was dispersed into the above solution. The mixture was first fixated and then treated by an oscillator for about 30 min. After that, the mixture was further treated by ultrasound ca. Six hours to form the uniform ink, which was dropped onto glassy carbon electrodes (5.61 mm diameter, disk area: 0.2475 cm 2 ) with a mass loading of 0.2 mg cm −2 . The prepared electrodes were naturally dried in air.
Electrochemical characterization. All of the electrochemical tests were performed in the room temperature using a typical three-electrode setup. The graphite rod (5.0 mm, 99.997% metals basis) and Ag/AgCl electrodes were used as the counter and reference electrodes, respectively. The electrolyte was 0.5M H 2 SO 4 solution. A Zennium Zahner electrochemical workstation was used for the electrocatalytic measurements. Before the test, the as-prepared anodes were activated by a chronopotentiometry scan with the −30 mA cm −2 current density for 2 h. The Tafel slopes were obtained from the polarization curves by plotting overpotential against log(current density). The steady-state activity was evaluated by chronopotentiometry measurements under different potentials. Please note that the scan rate for the electrocatalytic measurements was 5 mV/s. The C dl was determined by measuring the capacitive current associated with double-layer charging from the scan rate CV-dependence. In our present work, the CV potential window was chosen to be −0.2 to 0 vs Ag/AgCl. The scan rates were 20, 40, 60, 80, and 100 mV s -1 . The double-layer capacitance (C dl ) was estimated Δj = (j charge − j off charge ) at −0.1 V vs Ag/AgCl against the scan rate. The linear slope was twice of the double-layer capacitance C dl . EIS measurements were carried out after the HER tests under the different potentials. The ZSim Demo software was used to fit the EIS results. The measured potentials vs. Ag/AgCl were converted to the reversible hydrogen electrode (RHE) scale according to the following Nernst equation: E RHE = E Ag/AgCl + 0.059 pH + E°A g/AgCl , where E RHE is the converted potential vs. RHE, E°A g/AgCl = 0.1976 at 25°C, and E Ag/AgCl is the experimentally measured potential against Ag/AgCl reference.
Computation model and details. The optimized lattice constants of the Mo 2 N unit cell with 4/mmm point group were a = b = 4.2627 Å and c = 8.0876 Å. Then four types of five-layered (112), (200), (312), and (220) surface with a 15 Å height were cleaved from their bulk phase and the top three layers were allowed to relax. The single-point energy and the geometry optimization were implemented in the Vinnea ab initio software package (VASP) within the framework of density functional theory 63 . The exchange-correlation functional adopted the generalized gradient approximation of the Perdew-Burke-Ernzerhof (GGA-PBE) formula 64 . The projector augmented wave method with a cutoff energy of 520 eV was set. The optimization convergence standards were 10-4 eV for energy and −0.05 eV/Å for force. The formation energy E f of different surface and P-doping configurations were calculated by: E f ¼ ðE total À n N μ N À n Mo μ Mo À n P μ P Þ=ðn N þ n Mo þ n P Þ; ð1Þ where the n N , n Mo , and n P were the total number of N, Mo, and P atoms and μ N , μ Mo and μ P were the chemical potential of the element N, Mo, and P defined from the nitrogen gas molecule, bulk Mo metal and black phosphorus phase. E total was the total energy of the configuration. For the adsorption free energy of hydrogen atom, it was calculated based on the formula: where the E *H , E * , and E H 2 were the total energy of the surface with and without hydrogen adsorption and the hydrogen gas molecule. Then 0.24 eV was added to correct the contribution of the zero-point energy and entropy based on previous report 65 .
Data availability