DFT Investigation of Ammonia Formation via a Langmuir–Hinshelwood Mechanism on Mo-Terminated δ-MoN(0001)

In this work, we employed density functional theory to elucidate the energetics associated with elementary steps along a Langmuir–Hinshelwood mechanism for the Haber–Bosch synthesis of ammonia from N2 and H2 on a hexagonal, Mo-terminated molybdenum nitride surface. Using nudged elastic band calculations, we determined the energy barriers involved in the reaction processes. An active site consisting of four nearest-neighbor Mo atoms, previously identified as an active site on similar surfaces, was chosen to investigate the reaction processes. Using this approach, we calculate a barrier of ∼0.5 eV for the dissociation of N2. The superior activity of the dissociation of the strong N2 bonds is rationalized based on the unique geometric and electronic configurations present at these active sites. Despite the favorable energetics for nitrogen dissociation, the energy cost for hydrogenation of NHx (0 ≤ x ≤ 2) species is shown to be energetically limiting for the formation of ammonia through the Langmuir–Hinshelwood mechanism at these sites, with elementary step activation barriers calculated to be as large as ∼2 eV. A comparison to Haber–Bosch results derived from a similar γ-Mo2N model system suggests the relative independence of surface chemistry and bulk stoichiometry for rhombic Mo4 active sites present on molybdenum nitrides.


INTRODUCTION
Ammonia is a very important compound used in the production of many useful goods, for example, fertilizers, other nitrogenous compounds, cleaning agents, and antiseptics, and for fermentation. Industrial synthesis of ammonia is carried out via the Haber−Bosch process, which consists of the direct reaction of nitrogen with hydrogen over iron-based catalysts. Efficiently achieving usable yields of ammonia requires ultra-ambient reaction conditions (P ≥ 100 atm. and T ≥ 400°C), which makes this process extremely energy intensive. For perspective, it is estimated that about 2% of the world's energy usage is consumed by the Haber−Bosch process, 1 with most of that supply being extracted from fossil fuels. Apart from energy utilization, this reaction also currently depends upon fossil fuels to provide hydrogen gas feedstocks. In an effort to produce relatively "green" ammonia, many researchers have devoted themselves to the search for better catalysts amenable to this process. Though Ru-based catalysts are found to be more efficient than iron ones, 2 their use is economically hindered by limited supplies of terrestrial ruthenium. Search for a more cost-effective catalyst is required to decrease the energy requirements for this process and bring down the cost of ammonia (in terms of energy consumption, materials costs, and environmental impact).
Thanks in large part to the Nobel Prize winning work of Ertl, 3 conventional Haber−Bosch reactions are known to follow a Langmuir−Hinshelwood mechanism by activating NN dissociation on fourfold Fe sites following electron transfer from the substrate. A similar mechanism is also predicted for Ru-based catalysts when modeling activity at fivefold-coordinated Ru steps. 4,5 For these active monometallic catalysts, adsorption and dissociation of triple-bonded N 2 is the most energy-intensive step in the reaction process. 6 More generally, Nørskov and his coworkers have shown that N 2 adsorption strength can be used as a proxy to predict the Haber−Bosch activity on different catalytic interfaces. 7 If nitrogen adsorbs too weakly on a material, then the reaction is limited by the inefficient dissociative adsorption of N. On the other hand, if N 2 adsorbs too strongly, then the reaction becomes limited by inefficient nitrogen hydrogenation and subsequent NH 3 product desorption. 8 Hence, optimal catalysts should present intermediate N 2 adsorption energies to best balance the needs of these opposing processes in accordance with the Sabatier principle.
In the past 2 decades, early transition-metal nitrides (TMNs) have been shown to activate ammonia formation. 9−12 Hybridization between orbitals of the metals and non-metals in these materials has been shown to lower the Fermi energy of these systems (relative to the pure metals), such that some typically less-active metals, like molybdenum, can be used to create binary complexes with electronic properties more closely approximating those exhibited by more catalytically active late transition metals. 13 Introduction of non-metal atoms (e.g., N) may also potentially change interatomic spacings in manners structurally favorable for improved adsorption and activation of N 2 molecules on such surfaces. Beyond Haber− Bosch, TMNs also show promising results as active sites for electrochemical conversion of dinitrogen to ammonia. 14,15 Within the TMN family, molybdenum nitrides have received considerable attention, with several recent studies suggesting enhanced catalytic properties important to the Haber−Bosch reaction process. 16−20 In addition to presenting promising catalytic properties, the molybdenum nitride family of materials also exhibit a host of other physical properties of broad scientific interest. For example, several crystalline phases are known to be low-temperature superconductors, 21,22 and hexagonal δ-MoN exhibits a hardness comparable to diamond and boron nitride. 23 Several groups have successfully demonstrated different approaches allowing for the controllable synthesis of various phase-pure molybdenum nitride structures. 24−26 Of particular interest to this paper is the recent work of Khaniya and Kaden, 24 which reports the epitaxial growth of MoN films on Ru(0001) via an ion-assisted and physical vapor deposition approach. The films produced in that work appear ripe for direct comparisons between surface-science experiments and density functional theory (DFT) predictions due to their single-crystalline and atomically planar terminations, which exhibit a long-range order and elemental stoichiometry compatible with δ-MoN(0001).
A computational work done by Zhao et al. investigated the catalytic properties of a γ-Mo 2 N(111) surface used to activate the Haber−Bosch reaction. 8 In that case, a rhombus containing four nearest-neighbor Mo atoms (Mo 4 ) was identified as an active site for the adsorption of dinitrogen, with NN bond dissociation requiring only 0.58 eV when binding in highly coordinated arrangements at these sites. By contrast, subsequent NH x hydrogenation at the same site was shown to be rate-limiting, with step-wise barriers reaching as high as 2.01 eV. Enhanced surface nitridation (likely to occur under Haber−Bosch conditions) was shown to generally lower the hydrogenation barriers, while impacts on nitrogen dissociation were shown to depend more strongly (and anisotropically) on N adatom locations relative to adjacent rhombic Mo 4 ensembles.
The first atomic layer in Mo-terminated γ-Mo 2 N(111) and δ-MoN(0001) is near-identical, suggesting that similar Mo 4 site behavior might be expected from the δ-phase interface (Mo-terminated surfaces of γ-Mo 2 N(111) and δ-MoN(0001) are shown in Figure 1 for comparison). Based on (i) previous predictions indicating advantageous conditions for Haber Bosch catalysis on rhombic configurations of the nearest neighbor Mo sites present on Mo-terminated γ-Mo 2 N(111), (ii) the existence of near-identical surface site configurations on Mo-terminated δ-MoN(0001), and (iii) the existence of an established recipe for controllable growth of single-crystalline, hexagonal MoN thin films, we aim here to compare and contrast predictions for the energetics associated with each of the elementary steps along a Langmuir−Hinshelwood mechanism for ammonia production on a δ-phase model to those previously reported on analogous γ-phase sites. Results from this comparison may help further disambiguate the degree to which surface versus bulk modification (i.e., N enrichment/depletion) might be expected to affect the resultant Haber−Bosch surfacechemistry in subsequent experimental investigations on similar model systems.

RESULTS AND DISCUSSION
2.1. Adsorption Characteristics. To inform our selection of initial state (IS) and final state (FS) configurations used in the nudged elastic band (NEB) calculations, we first considered several possible adsorption geometries for each of the sequentially produced atomic/molecular species required to form ammonia via a Langmuir−Hinshelwood mechanism (i.e., N 2 , N, H 2 , H, NH, NH 2 , and NH 3 ). A minimum of four distinct adsorption sites (bridge, fcc, hcp, and top) were considered in each case (see Figure 1). For clarity, the four sites can be distinguished as follows: top sites refer to adsorbates placed directly above surface Mo atoms, bridge sites refer to adsorbates placed directly between two neighboring surface Mo atoms, hcp sites refer to adsorbates placed in threefold hollows between three neighboring surface Mo atoms in positions directly above second-layer N atoms, and fcc sites refer to adsorbates placed in threefold hollows between three neighboring surface Mo atoms in positions directly above areas without second-layer N atoms. In order to clearly represent the adsorption configurations of the adsorbates, both "nearest bond lengths d" and "adsorbate heights [h]" are mentioned along with adsorption characteristics. Nearest bond lengths and adsorbate heights are calculated as the distance from the vertical position of the lowest atom within the adsorbed species to the nearest surface Mo atom and the average vertical position of top-layer Mo atoms, respectively.
2.1.1. Adsorption of N and N 2 . Of the four N/MoN trial configurations, only three proved stable. N atoms placed in bridge sites spontaneously migrate to hcp sites upon relaxation (see Table 1). Side and top views of the three stable binding motifs are shown in Figure S1. Of the stable configurations, individual N atoms were found to adsorb most strongly in hcp sites (E ads = 3.20 eV), indicating a preference for terminal adsorption in sites having continued registry with the remainder of the bulk structure. Top and fcc sites also provide stable, albeit less favorable configurations (E ads = 0.75 and 2.55 eV, respectively). Whereas weaker adsorption at top sites results in considerably less substrate-to-N electron transfer To accommodate the greater degrees of freedom associated with the molecular adsorption of dinitrogen, the number of trial configurations was increased from four to eight when searching for stable N 2 /MoN binding arrangements. A first set of molecules were placed in vertical configurations (with the N 2 bond aligned along the surface normal) at top, hcp, and fcc sites. A second set were placed in horizontal configurations with the N 2 molecule centered above the same series of positions in addition to the bridge site (oriented with the N N bond aligned along the Mo−Mo bridge bond). Within this second set, two azimuthally rotated N 2 configurations (NN surface orientations rotated by 30°with respect to one another) were separately considered at top sites to differentiate the molecules oriented parallel to the hexagonal arrangement of Mo rows from those bridging between threefold hollows on the opposite sides of the host atom. Optimization of these trial structures converged upon four stable adsorption arrangements depicted and characterized in Figure 2 and Table 2, respectively.
Most trial structures significantly distort to take on one of two new binding arrangements following relaxation (distinguished as "bridge" and "fcc-top"), suggesting a preference for non-vertical adsorption in more highly coordinated configurations. By contrast, the remaining arrangements (distinguished as "top-v" and "top-h") indicate the stable adsorption of molecules remaining essentially unchanged relative to their trial configurations following relaxation. Comparing across these two groups, we find greater adsorption energy, MoN-to-N 2 electron transfer, and NN bond distortion for the statistically more likely arrangements resulting from N 2 coordination to multiple Mo atoms. Consistent with these tendencies, the most strongly bound configuration (bridge; E ads = 2.54 eV) also features the greatest degree of N 2 −Mo coordination (both N atoms bind over three-fold hollow sites), the closest average molecular proximity to the surface (1.3 Å), the greatest N 2 bond elongation (0.20 Å growth vs gas phase geometry), and the largest amount of charge transferred from the surface to the molecule (1.31 e − evenly distributed across both atoms), such that the atoms in the molecule most closely resemble those in the two most favorable atomic adsorption configurations. Perhaps unsurprisingly, this arrangement bares a near-perfect agreement to that previously identified as the active configuration for N 2 dissociation on γ-Mo 2 N(111), 8 including a mild distortion of the Mo 4 site upon adsorption (0.4 Å contraction of the Mo atoms positioned at the vertices spanning the long axis of the rhombic ensemble). By comparison, top site configurations transfer far less charge onto the adsorbing molecule and leave the N 2 bond lengths unperturbed, suggesting poor activation toward dissociation. Adsorption in the fcc-top configuration results in intermediate effects, with the atom coordinating to the threefold hollow (designated as the "first N" in Table 2) exhibiting charge localization characteristics similar to atoms in the bridge configuration, and the atom coordinating to the top site taking on characteristics similar to atoms in the top-h configuration by Bader charge analysis. Given similar bond elongation characteristics (vs bridge site adsorption), this latter N 2 /MoN adsorption configuration appears to provide an alternative activation site also likely to be suitable for subsequent dissociation.
To allow for more direct comparison to previous literature, we reoptimized the pristine MoN, gas-phase N 2 , and N 2 /MoN structures to calculate a new bridge site E ads value using the Perdew−Burke−Ernzerhof (PBE) functional in a manner consistent with that used for the γ-Mo 2 N(111) model system. 8 Changing from PBE to optB88 causes E ads (bridge) to drop from 1.61 to 2.54 eV, indicating a 0.93 eV increase in this value presumed to predominantly derive from the inclusion of van der Waals interactions in the optB88 calculations. 27 The PBEcalculated E ads for an analogous (bridge) N 2 /γ-Mo 2 N(111) configuration was previously reported to be 2.27 eV, 8 indicating a 0.66 eV decrease in the strength of the attractive N 2 -Mo x N y interactions for this arrangement when transitioning from γ-Mo 2 N(111) to δ-MoN(0001). Whereas transition from the γto δ-phase model requires doubling the sub-surface concentration of N present in sixfold coordination sites between each subsequent Mo layer within the bulk, the introduction of a single N adatom into a threefold hollow adjacent to the bridge-bound nitrogen molecule has been shown to decrease the same configuration's (PBE-calculated) E ads (N 2 ) by as much as 1.23 eV on the γ-phase model system. 8 Consistent with this stark contrast between the impact from  Table 3 for NH/MoN, NH 2 /MoN, and NH 3 /MoN adsorption in different configurations. Trial configurations were, respectively, limited to linear, trigonal, and tetragonal molecular arrangements designed to maximize surface-H and H−H distances within the MoN-supported NH, NH 2 , and NH 3 structures. Optimization of the NH/MoN trial structures produces four stable adsorption configurations depicted in Figure S2. Adsorption energies associated with these structures range from 6.5 to 8.2 eV, with the most stable configuration resulting from adsorption at the hcp site. General trends in both surface-to-molecule charge transfer and molecular adsorption energy are similar to those observed for the adsorption of N atoms, with decreased surface coordination leading to decreased charge transfer and lower adsorption energies. Similarly, the fcc configuration again exhibits lower adsorption energy but slightly greater charge transfer to the adsorbed molecule relative to the hcp structure. A site-by-site Bader analysis comparison of analogous NH/ MoN and N/MoN structures implies net electron redistribution from the NH x adsorbate to the support upon N−H bond formation at all applicable N/MoN configurations. Unlike N, NH adsorption is also stable at bridge sites.
All NH 2 /MoN trial configurations converged into one of the only two stable adsorption structures (see Figure S3). In this case, similar characteristics are exhibited at both hcp and fcc sites, with a slight preference for fcc configurations indicated by a slightly greater adsorption energy (4.72 vs 4.53 eV) and the tendency for both top and bridge site migration to fcc rather than hcp sites upon relaxation.
H addition required for H−NH bond formation again leads to further NH x adsorbate-to-support electron redistribution when comparing the Bader analyses of analogous NH 2 /MoN and NH/MoN configurations. Unlike the first H, the addition of the second H significantly weakens the NH x adsorbate− surface interactions when comparing NH 2 and NH adsorption energies and adsorbate−surface separations at common sites.
Like NH 2 /MoN, we again find only two stable NH 3 /MoN adsorption structures after relaxing several trial configurations (see Figure S4). Unlike NH and NH 2 , which exhibit preferences for highly coordinated adsorption at threefold hollow sites, NH 3 instead exhibits a converse preference for less coordinated adsorption at top sites (E ads = 1.51 eV), with all but fcc trial configurations migrating to top sites upon relaxation. The general molecule surface bond weakening and electron-transfer trends continue when adding the third H, with Bader analysis showing NH 3 charge localization roughly equivalent to that in free molecules for both stable adsorption configurations. NH bond lengths (1.03 Å) and orientations (∼109°H−N−H bond angles) also closely match those for isolated NH 3 in both cases, suggesting favorable conditions for product desorption from either site once formed. NH bonds, within the adsorbed ammonia molecules, show an affinity for azimuthal alignment toward the nearest neighbor Mo atoms, while also exhibiting a modest repulsion away from nearby subsurface N atoms. These tendencies combine to result in a slight (∼7.2°) rotation of top site ammonia relative to the hexagonal pattern of Mo atoms, such that each NH bond orients with a slight preference toward adjacent fcc rather than hcp hollows.
2.1.3. Adsorption of H and H 2 . Adsorption energies, heights, and charge-transfer characteristics are given in Table  4 for H/MoN and H 2 /MoN adsorption in different configurations. Atomic H stably binds in four distinct adsorption structures following relaxation of each of our trial configurations (see Figure S5). Like N, H atoms also bind First N" refers to nitrogen atoms positioned closer to the surface plane in Figure 3.   Figure  S6, only two show evidence of significant H 2 −MoN interactions (top-h and top-h-R30). By contrast, all three vertical adsorption configurations show near-zero adsorption energies but nonetheless present trapping barriers sufficient to prevent relaxation into one of the horizontal arrangements upon optimization. In the more favorable top-h configuration, H 2 binds 1.85 Å above the surface with an adsorption energy of 0.82 eV and results in a 0.14 Å elongation of the H−H bond after a transfer of 0.13 electrons from the surface to the molecule. Azimuthally rotating the top-h H 2 molecule by 30°a bout the surface normal results in similar, albeit slightly weaker, adsorption characteristics (see top-h-R30 in Table 4/ Figure S6). While not discussed below, a representative barrier for the dissociation of top-h-R30 H 2 into two nearest neighbor fcc sites was found to be 0.05 eV, which we take to imply the facile dissociative adsorption of hydrogen on Mo-terminated δ-MoN(0001).
2.2. Ammonia Formation Mechanism. 2.2.1. N 2 Dissociation. NEB was used to track changes in total energy when forcing a transition from the most stable N 2 /MoN configuration (bridge) to a dissociated configuration consisting of two singly occupied threefold hollow sites immediately adjacent to that bridge site (hcp + fcc). One might intuitively expect this process follows a bond elongation deformation mechanism (see Figure 3). In the IS, a nitrogen molecule adsorbed in a bridge configuration is about 20 percent activated for dissociation (inferred by elongation in N−N bond length compared to gas-phase value). The highest energy image generated along the path connecting the IS and FS endpoints has been labeled TS and is calculated to lie 0.52 eV above the energy associated with the IS configuration. This value compares quite well with N 2 dissociation barriers calculated for both stepped Ru 4 and an analogous rhombic Mo 4 ensemble modeled as an active site on γ-Mo 2 N(111). 8 Consistent with individual hcp versus fcc N/MoN adatom adsorption characteristics, the N 2 molecule tilts such that N(1) repositions slightly closer to the surface than N(2) throughout the dissociation process. For clarity, Mo(1) bounds the hcp site and Mo(3) bounds the fcc site.
Unlike the Mo 2 N results, we note far smaller changes in the dihedral angle formed between the fcc and hcp threefold hollows contained within the Mo 4 active site following dissociation (∼5°here vs ∼20°in the previous work). Decreased torsional deformation could reflect a difference in active site rigidity resulting from changes in the bulk crystal structure. Despite near-negligible changes in a torsion angle, we do note significant lateral rearrangements within the Mo 4 ensemble during dissociation. Specifically, Mo(2) and Mo(4) contract ∼0.5 Å toward one another and Mo(1) and Mo(3) migrate ∼0.3 Å away from one another while the system traverses the IS−TS activation barrier. The reconstructed Mo 4 ensemble, which alleviates much of the surface strain introduced by the IS molecular adsorption configuration, is largely preserved upon relaxation from the TS into the fully dissociated FS structure.
2.2.2. Electronic Structure of N 2 during Dissociation. Bader charge differences tracking changes in electron accumulation at each of the four interacting molybdenums and both nitrogen atoms are provided in Table 5 for the IS, TS, and FS N 2 dissociation configurations described above. Transfer of electrons from the active site to the molecule is evident throughout the dissociation process. Beginning with the IS configuration, each N within the adsorbed molecule.
M carries more than 0.6 excess electrons (relative to isolated N 2 ), with the bulk of this excess charge balanced by ∼0.25 e − depletions from each of the four directly interacting Mo 4 atoms (relative to pristine MoN). As the adsorbed configuration distorts, surface-to-molecule charge transfer results in increased electron accumulation on both N atoms. While approaching the TS, increased N 2 electron accumulation is accompanied by anisotropically enhanced electron withdrawal from the Mo (2) and Mo(4) atoms, which both contract toward the mildly   (3) atoms, which the N adatoms migrate toward following dissociation. Figure 4a,b provides spatially resolved representations allowing for the direct visualization of charge accumulations and depletions described in the Bader analyses of the IS and FS configurations reported above. The overall increase in charge accumulation (depletion) surrounding the N (Mo) atoms upon adsorption and dissociation is readily apparent upon visual inspection. The general geometries and symmetries associated with the areas of enhanced/diminished charge density are consistent with those reported for analogous N 2 / Mo 4 configurations on γ-Mo 2 N(111). 8 As such, we presume qualitatively similar charge-transfer characteristics, whereby Mo atoms positioned along (across) the N−N bond axis interact with a single N via σ-bonds (simultaneously interact with both N 2 π-bonds). Interestingly, both structures (IS and FS) show the depletion of the electronic charge presumed to have been delocalized between Mo atoms prior to N 2 adsorption and dissociation. While also present, similar effects are less pronounced for analogous configurations present on the N 2 /γ-Mo 2 N(111) model system. 8 Figure 4c,d provides direct comparisons of the projected density of states (DOS) at each of the six directly interacting Mo and N atoms before and after dissociation by plotting DOS(FS) − DOS(IS) subtraction profiles for each atom common to both structures. For obvious reasons, dissociationinduced DOS changes are most apparent on the N atoms (note that the vertical axes in the Mo plots have been rescaled by a factor of ten to make the smaller changes more discernible for those sites). After dissociating, deep-lying states associated with N−N bonding interactions are completely lost (see negative density differences at energies below −6 eV) and replaced with new Mo-affiliated states located closer to the Fermi level. Increased N−Mo bonding interactions lead to mutually increased DOS near −5 eV at each of the six atoms after dissociation. Moreover, small hcp versus fcc related differences in the N(1) and N(2) plots correlate well with differences in the near-Fermi change in DOS at the correspondingly interacting Mo(1) and Mo(3) sites. General trends are again similar to those noted for Mo 4 -activated dissociation on γ-Mo 2 N(111); 8 however, N adatom-to-adatom asymmetries appear far less pronounced on our δ-MoN(0001) model, possibly reflecting a more homogeneous Mo-coordination to sub-surface N atoms in the δ-phase system.
In summary, N 2 dissociation on Mo-terminated δ-MoN(0001) can be viewed as a two-step process. The process is first fully activated following Mo 4 structural rearrangements associated with the donation of the electronic charge to the molecule via Mo(2,4) interactions across its π-bonds and concomitant N 2 bond distortion. Once fully activated, the elongated molecule readily dissociates to relax into its FS structure, within which the resulting adatoms accumulate more electronic charge through σ-bond interactions with Mo(1,3) atoms aligned along the parent molecule's primary axis. Despite the subtle differences in reaction-induced active-site reconstruction and qualitative electron redistribution characteristics, these results are otherwise consistent with those reported for N 2 dissociation on Mo-terminated γ-Mo 2 N(111). This includes reaction energies, activation barriers, and net adsorbate-support charge transfers all falling within ∼10% of those in the former work (with all values skewing smaller in our work). 8 In addition to the differences in the bulk structure and N-content between previous studies and this, discrepancies in reaction energies and energy barriers can arise due to a lack of inclusion of van der Waals interactions (PBE vs optB88), which are shown to be important for the selection of suitable catalysts. 28,35 2.2.3. NH x Hydrogenation. To most efficiently draw comparisons to past works on γ-Mo 2 N(111), we have only considered hydrogenation pathways resulting from the sequential H adatom attack of isolated N(H x ) adsorbates. Nonetheless, we expect similar anisotropic impacts to overall Haber−Bosch catalysis in the presence of co-adsorbed N for our δ-phase model due to its similar N 2 /Mo 4 active site geometry and interaction asymmetries. Using NEB, individual reaction barriers, energies, and transition states (TSs) were obtained for each of the three sequential hydrogenation iterations. IS configurations consist of NH x and H adsorbates positioned in their most proximate combination of individually most stable adsorption geometries. On this basis, all IS configurations consist of NH x species placed in threefold hollows (hcp for N and NH and fcc for NH 2 ) immediately adjacent to an fcc-bound H adatom. Using the same logic, the most stable NH x+1 configurations were chosen as FS products for each elementary hydrogenation step. This overall reaction pathway is consistent with the more direct of the two ammonia formation mechanisms explored in the previous γ-Mo 2 N(111) study. 8 Initial hydrogenation proceeds via simple H diffusion from an adjacent fcc-to top-site configuration (see Figure 5). The energy cost associated with this process (E a = 1.42 eV) is slightly larger than the difference in adsorption energy associated with top versus fcc H/MoN configurations on the pristine interface and again falls within ∼10% of the barrier As depicted in Figure 6, NH 2 formation proceeds in a manner analogous to the first N hydrogenation step, with the barrier to reaction again stemming from the energy required for H diffusion from a hollow onto an adjacent Mo (E a = 2.00 eV). In this case, the Mo−H TS bond tilts away from the adjacent NH, such that the H−NH distance only reduces from 3.30 → 2.70 Å when traversing the activation barrier. As the system relaxes out of the TS, H−H repulsion leads to the depicted NH 2 FS configuration (θ H−N−H = 101.6°), with the full hydrogenation step exhibiting an endothermic reaction energy of 1.68 eV. Unlike every prior step, the reaction barrier noted here is near-identical to (within 0.5%), while the NH−H TS geometry is significantly different from, that observed along an analogous reaction pathway in the γ-Mo 2 N(111) study. 8 Unlike N and NH hydrogenation, overcoming the final barrier to NH 3 formation (E a = 1.94 eV) includes the coordinated diffusion of both NH 2 and H adsorbates from the nearest neighbor fcc sites onto a mutually common adjacent top site (see Figure 7). During this process, the H and NH 2 moieties contract about 1 Å toward one another while the latter rotates to maximize H−H−H separations in the TS geometry. Once formed, the TS readily transforms into the FS NH 3 /MoN configuration, with the full hydrogenation step again proceeding endothermically, with a reaction energy of 1.03 eV. In this case, a similar TS configuration was again noted on γ-Mo 2 N(111), but the reaction pathway leading to this structure is different due to a preference for NH 2 adsorption at bridge sites within the previous model system. 8 The significantly larger barrier noted for this hydrogenation step on the γ-phase model (E a = 2.56 eV) most likely at least partially results from this IS discrepancy between the two models. Like the previous study, progressing from IS → TS results in significant adsorbate-to-Mo 4 electron transfer (not shown), with the diffusing H adatom contributing the majority of this charge while traversing the barrier to each sequential hydrogenation step.
A full summary of the energy profile associated with the most energetically favorable pathway for nitrogen dissociation and subsequent stepwise hydrogenation is provided in Figure  S7. While the focus of this work is specific to the Haber−Bosch reaction, results may also provide some insights into related nitrogen reduction reaction (NRR) mechanisms, which have become increasingly important within the electrochemical catalysis community. Impacts resulting from simultaneous interactions between multiple N and H atoms are likely to be important to elucidating the improved understandings of both reaction processes (Haber−Bosch and NRR), which we aim to address (at least in part) in a forthcoming publication.

CONCLUSIONS
Comparing our results for δ-MoN(0001)-activated N 2 dissociation and subsequent hydrogenation to those from γ-Mo 2 N(111) models exposing equivalent rhombic Mo 4 active sites suggests predominantly bulk-invariant surface chemistry at such sites. 8 Despite drastic changes to bulk stoichiometry and local (sub-surface) Mo−N coordination, Mo 4 sites present at both interfaces are predicted to exhibit extremely similar N 2 adsorption/dissociation and NH x hydrogenation characteristics, including adsorption configurations, electron redistributions, and overall reaction pathways. For example, when treated at the same level of accuracy, N 2 adsorption at the Mo 4 site only weakens by 30% when changing the structure from γ-Mo 2 N(111) to δ-MoN(0001); a change that increases the number of N atoms directly coordinating to the Mo 4 active site from five to eight. Notably, the placement of as few as only one N adatom above one of these three contrastingly populated hollows (omitting sites that become occupied following N 2 dissociation) is predicted to reduce E ads (N 2 ) by as much as ∼40% on the γ-phase model, with increased co-adsorption expected to only further weaken the interaction; a discrepancy highlighting the relatively inconsequential impact of bulk (vs surface) changes in the Mo/N stoichiometry/structure on Haber−Bosch catalysis. Similar to the γ-phase results, pristine Mo 4 sites are highly active toward N 2 adsorption and dissociation (E ads = 2.54 eV and E a = 0.52 eV), but rate limited by large barriers to subsequent hydrogenation [E a (NH + H → NH 2 ) = 2.00 eV]. Likewise, hydrogenation barrier height reductions will be expected when increasing N coverage around the δ-MoN-contained Mo 4 sites, with effects expected   to vary in a manner showing anisotropic dependence on the binding site locations of the co-adsorbates for reasons analogous to some of those previously described for the γ-Mo 2 N interface. 8 As such, improved Haber−Bosch activity on the δ-MoN-contained Mo 4 sites should be achievable by controlling the reaction conditions to beneficially tailor the [N] surface concentration. Moreover, results from this work further suggest that such adjustments can include radical changes to bulk Mo x N y stoichiometry without introducing significant changes to Mo 4 chemistry, so long as a hexagonal δ-MoN(0001)-like terminal structure is preserved.

COMPUTATIONAL DETAILS
Calculations were performed in the framework of the DFT 29 using plane wave pseudo-potentials implemented in VASP (Vienna Ab-initio Simulation Package 5.4.4). The generalized gradient approximation 30−32 was used with the optB88 package 33 to account for van der Waals interactions. This approach to modeling van der Waals effects has previously shown improved agreement between modeled and experimental parameters for establishing the activity of catalysts. 34,35 TSs were isolated using the NEB method to interrogate configurations and energies associated with three or five images found along pathways linking fixed IS and FS reaction configurations. Relaxed structures were optimized until Hellman−Feynman forces on all atoms were less than 0.01 eV/Å. An energy cutoff of 400 eV was used along with a Monkhorst−Pack grid 36 of 2 × 2 × 1. Lattice constants of bulk hexagonal MoN were determined using a Monkhorst−Pack grid of 18 × 18 × 18 with a cutoff energy of 800 eV. The calculated lattice constants (a = b = 2.88 Å and c = 2.86 Å) were found to be close to the experimentally determined values (a = b = 2.87 Å and c = 2.81 Å). 37 A molybdenum nitride surface was modeled using a 5 × 5 super cell with five layers. The bottom two layers of the slab were fixed to their bulk positions and a vacuum of more than 10 Å was maintained between periodic images of the slab in the z-direction (perpendicular to the surface). Adsorption energies are defined as with E sub representing the energy of the pristine substrate, E mol representing the energy of the molecule in the gas phase, and E mol−sub representing the energy of the surface after adsorbing the molecule. Using this convention, positive adsorption energies indicate thermodynamically favorable configurations (and vice-versa). For single atom adsorbates (N and H), half of the gas-phase energy associated with the corresponding diatomic species is used to calculate E mol . VASPKIT 38 and VESTA 39 were used for post-processing and visualization, respectively. A few key structures were selectively reoptimized using the PBE functional to provide direct comparisons to previous results reported at that level of accuracy for analogous γ-Mo 2 N(111)-supported configurations. 8 ■ ASSOCIATED CONTENT
Optimized configurations and reaction energy profiles (PDF)