Ammonia cracking on single-atom catalysts: A mechanistic and microkinetic study

Ammonia cracking has been identified as a crucial step to unlocking a sustainable economy. Using the density functional theory, we modeled transition metal single-atom catalysts (SACs) supported on graphene and nitrogen-modified graphene to investigate the catalytic NH 3 cracking process. The results showed that (i) N-modified graphene secures the transition metal atoms (M) stronger than C-matrixes, and (ii) structures with three anchoring nitrogens (MN 3 ) are more reactive than MN 4 ones. On IrN 3 and RuN 3 SAC models, the N 2 evolution determines the total rate, while, on RhN 3 -SAC, it is the NH 3 dehydrogenation. Temperature-programmed desorption simulations on SACs showed variations compared to extended metal surfaces. Batch reactor simulations were employed to balance the sequence of elementary steps as a function of the temperature, revealing the overall NH 3 cracking activity. Results suggested IrN 3 and RhN 3 are strong candidates for NH 3 cracking at temperatures as low as 230 ◦ C.


Introduction
Hydrogen (H 2 ) is a viable and sustainable long-term energy vector that can address increasing global energy demands, and its utilization is compatible across vehicles, industries, and heating and cooling applications [1,2].However, H 2 is highly volatile and flammable, and its distribution using pipes presents challenges due to leakage [3].Instead, ammonia (NH 3 ) is a potential hydrogen vector with high hydrogen density (17.6 wt%), carbon-free, and suitable for storage at low pressure (~7.5 atm at 300 K) [4,5].
The endothermic NH 3 cracking process could provide a steady H 2 supply for clean energy technologies such as fuel cells [6].However, the cost and low energy efficiency of current NH 3 cracking catalysts hinder their implementation on a large scale [5,7,8].The literature on NH 3 cracking using heterogeneous catalysts is extensive; recent reviews can be found elsewhere [9][10][11][12][13][14][15].These traditional catalysts are commonly based on costly late transition metals.The initial capital of these catalysts can be significantly reduced using single-atom catalysts (SACs) [16].SACs present the highest level of atom efficiency as each atom participates in independent catalytic reactions.
SACs have recently attracted attention in the nitrogen reduction reaction (NRR), i.e., the reverse process of ammonia cracking [17][18][19][20][21]. Wang and co-workers found that Au SACs supported in g-C 3 N 4 have higher NRR electrocatalytic activity and selectivity than Au(211) surfaces [22].The Au/g-C 3 N 4 faradaic efficiency (FE) is 11.1 %, lower than for Ru SACs on N-doped carbon (FE = 29.6 %), which also benefits from reducing the hydrogen poisoning and the H 2 evolution reaction compared with Ru nanoparticles [23].A density functional theory (DFT) investigation found that the Ru atom in SAC structures is firmly anchored to the carbon matrix and presents a similar NRR limiting potential independently of the carbon-based materials, e.g., C 2 N, T-C 3 N 4 , and γ-graphene [24].Based on the fact that Fe is the industrial catalyst for the Haber-Bosch process and is present in the nitrogenase enzymes, [25].Li et al. modeled the highly spin-polarised FeN 3 center and found it to adsorb and activate N 2 [26].Lü et al. also reported a FeN 4 structure to be NRR active [27].Alternative transition metal atoms such as Sc, V, Mn, Co, Y, and Mo were embedded in graphene and tested as NRR catalysts, [26,28,29] the latest suggesting a Mars-van Krevelen reaction mechanism [29].
We considered the previous NRR observations and catalytic results on supported metal particles [7,8,10,13] to investigate the NH 3 cracking reaction on transition metal SACs supported in carbon and nitrogen-doped carbon structures using systematic DFT calculations.We went beyond the information in mechanistic energy profiles and implemented the derived information in microkinetic batch reactor models and temperature-programmed desorption (TPD) simulations, providing a multiscale approach to the computational results.

Computational details
We employed the Vienna Ab-initio Simulation Package (VASP) to model the NH 3 cracking mechanisms on transition metals single-atom catalysts (Fe, Co, Ni, Zn, Ru, Rh, and Ir) anchored on defective and pyridinic nitrogen-doped graphene layers [30,31].The spin-polarised revised Perdew-Burke-Ernzerhof (RPBE) method of the generalized gradient approximation (GGA) was adopted to describe the exchange and correlation energies with a plane-wave kinetic cutoff of 500 eV [32].Spin alignments were fully relaxed to the lowest energy configuration, avoiding any constraint to the system [33].Non-spherical contributions to atomic cores from the electron density gradient were represented by the projector augmented wave (PAW) [34].The zero-damping DFT-D3 method was used to describe long-range interactions [35].The optimized convergence thresholds of internal forces and electronic relaxation were set to 0.02 eV/Å and 10 − 5 eV, respectively.A 3 × 3 × 1 k-spacing Monkhorst-Pack grid sampled the Brillouin zone with a smearing broadening of 0.2 eV for bulk metals, 0.1 eV for SACs, and 0.02 eV for isolated molecules [36].Bader charges were calculated to assign partial charges to specific atoms [37].A negative value of the charge indicates the accumulation of electron density.
SACs' surfaces were represented by a graphite p(6 × 6) supercell slab model with a single atomic layer thickness.We added 15 Å of vacuum perpendicular to the slab to avoid spurious interaction with periodic images.Dipole correction perpendicular to the surface was applied upon the adsorption of any species.The strength of binding (E b ) a single metal atom on a defective carbon or nitrogen-doped carbon matrix was evaluated using Eq. ( 1).
Where E SAC , E g , and E m represent the energies of the supported SAC, the graphitic support (with or without N-doping), and the metal atom within the metal bulk in its most stable crystal structure.The last two energies define the binding reference energy.Table S1 of the Supplementary Information shows the optimized metal bulk lattice parameters.The molecular adsorption energies (E ads ) were defined by Eq. ( 2), and the relative energies (ΔE) along the energy profiles were calculated by Eq. ( 3).
Where E system is the total energy of the adsorbed structure, E SAC and E molecule denote the energy of the clean SAC and isolated molecular adsorbate, respectively.E NH3 and E H2 are the energies of isolated ammonia and hydrogen molecules.The last three terms are the reference energy maintained along all the reaction profiles.The half energy of a hydrogen molecule refers to the energy of one H atom, and n is the number of H dissociated from NH 3 .The reaction energy (E r ) of each step was given by the energy difference between the final state (FS) and the initial state (IS) (Eq.( 4)).An exothermic reaction step is associated with a negative E r .The transition states (TS) were determined using the climb-image nudged elastic band (ci-NEB) combined with the improved dimer method (IDM), ensuring a unique imaginary frequency along the reaction coordinate [38][39][40].Vibrational frequency calculations were carried out by constructing and diagonalizing the Hessian matrix, built from finite displacements of atomic positions of 0.05 Å in length.Only displacements of the adsorbates were regarded, i.e. the adsorbate vibrations were treated decoupled from the TM-SACs surface, since test calculations considering them yielded negligible changes in the gained vibrational frequencies at an exceedingly high computational cost.We defined the forward and reverse activation barriers (E a ) as the energy difference between TS and IS and between TS and FS, respectively (Eq. ( 5) and Eq. ( 6)).

Transition metal anchoring sites
We modeled single-and double-carbon vacancies in graphene and nitrogen-modified graphene (g-N 3 and g-N 4 ) supports [41].The vacant carbon and nitrogen sites are known to coordinate with metal atoms to saturate their dangling bonds [42][43][44].Fig. 1 represents the different structures according to the metal anchoring site.The significant surface properties of the SACs with the modeled transition metals are summarised in Table 1.
The binding energy of a transition metal atom within the nitrogenmodified graphene is stronger than in the non-doped defective carbon matrix, which explains the SACs' synthetic protocols [43,45].The binding energies in Table 1 show the SAC's formation of 4-coordinated structures (MX 4 ) to be thermodynamically more favourable than 3-coordinated ones (MX 3 ).Generally, the metal on MN 3 protrudes from the support plane driven by the N lone pair electrons (Fig. 1), affecting the M-N bonds and angles compared to the M-C bond.The resulting geometries in MN 4 and MC 4 are similar.Fe, Co, and Ni metal atoms have a stronger affinity with nitrogen than noble metals such as Ru, Rh, and Ir, meaning the former SACs are more stable than noble metal ones.A positive E b value indicates a thermodynamic driving force for the atoms to merge into a metallic bulk; our results on less durable SACs align with previous observations [42,44,46].

Ammonia cracking 3.2.1. Surface species
Based on the thermodynamically favourable metal interaction with the N-doped graphene, we investigated several non-equivalent NH x (x = 0-3) and H adsorption conformations on MN 3 and MN 4 SACs, and N and C sites neighbouring the M centres (M= Fe, Co, Ni, Zn, Ru, Rh, and Ir).A detailed description of adsorptions far from the metal centre can be found in Ref. [41].The most likely adsorption mode is on top of the transition metal atoms.The adsorption preference is explained by the partial occupation of hybridized sd metal orbitals interacting with the NH 3 lone pair electrons, supporting the electron back donation and activation of the N-H bond; see the MN x (x = 3 and 4) electronic details in Table 2 and density of states in Figs.S1-S7 [47].
A favourable NH 3 adsorption is essential to initiate catalytic molecular cracking.Table 3 summarises the adsorption energies of the adsorbents, i.e., NH 3 , N 2 , and H 2 .The NH 3 adsorption is generally more advantageous on MN 3 than on MN 4 .Another thermodynamic drive for ammonia cracking is the desorption of N 2 and H 2 , which should be more likely than their atoms remaining on the SAC.In this case, MN 4 structures predominate.The interaction strength and distance between the adsorbate and the adsorption site are commonly correlated on heterogeneous catalysts [48].Table S2 shows this is not the case with these SACs, as their interaction is confined to a single atom.
Another factor to consider in promising catalysts is a smooth energy profile along the reaction mechanism, i.e., intermediate species interacting neither too weak nor too strong [49].The focus on the first dehydrogenation is crucial as it has been identified previously as a rate-limiting step.Fig. 2  Considering the kinetic aspects is also crucial to identifying efficient catalysts; a generalized representation of the mechanism is shown in Fig. 3. From the shortlisted SACs, the activation energies on ZnN 3 are 2.13, 3.46, and 2.71 eV for the first, second, and third dehydrogenation, implying a kinetically hindered process.Overall, only three transition metal SAC candidates are not frustrated by NH 3 dehydrogenation, IrN 3 , RuN 3 , and RhN 3 , of which reaction energy profiles are shown in Fig. 4. The energy profiles (Fig. 4) include NH 3 adsorption, dehydrogenations, and nitrogen recombination and evolution, all linked through the transition state determining the kinetic feasibility of the overall endothermic NH 3 cracking; note that ΔE is not at standard conditions and includes the internal molecular energies leading to ΔE > ΔH of NH 3 formation (ΔH o f = − 46.1kJ⋅mol − 1 ).IrN 3 and RuN 3 profiles clearly display the highly endothermic nitrogen recombination, which is easier to achieve on RhN 3 .Another significant feature is the high activation energies for NH dehydrogenation on RuN 3 and RhN 3 .

Reaction Thermochemistry
The differences in Gibbs free energy (ΔG r ) and activation energy (ΔG a ) of each elementary surface reaction were determined to consider temperature effects on the NH 3 cracking process [48].We investigated the possibility of N 2 recombination through a Mars-van Krevelen mechanism but found energy barriers (>3 eV).Therefore, the surface NH 3 decomposition on 3 N-coordinated noble metals sites follows a Langmuir-Hinshelwood stepwise dehydrogenation (R1  On IrN 3 , the energy barriers for the stepwise dehydrogenation reactions are ~1 eV, more likely than on RuN 3 (ΔG a ~ 1.5 eV).Despite this, the endothermic N recombination reaction (R7) is the ratedetermining step with a barrier energy of ~2 eV (~1.5 eV for RuN 3 ), hindering the continuous generation of N 2 .Contrarily, on RhN 3 , the NH dehydrogenation (R5) controls the overall reaction rate with a ΔG a ~ 1.5 eV, although it is practically isoenergetic.The N 2 recombination (R7) on RhN 3 is just ΔG a ~ 0.6 eV.

Table 1
Geometric details and binding energies of modelled SACs.Label d(MX) (X = N, C) is the bond length between metal (M) and nitrogen (N) or carbon (C) atoms, ∠X 1 MX 2 indicates the angle between the MX bonds, and E b designates the binding energy.

Microkinetics
We calculated the partition functions of each gaseous and surface species involved in the process to derive their thermodynamic properties as functions of the temperature [50].For instance, the entropic terms include, for molecules, the electronic, translational, and rotational contributions besides the 3 N-6 vibrations in the vibrational partition function; for surface species, entropy includes frustrated translational and rotational contributions within the 3 N vibrations in the partition function.We implemented the temperature-dependent energies in the Eyring, Evans, and Polanyi approximation to compute the reaction constants of all elementary surface reactions [51][52][53].These constants were used to construct a kinetic model of the reversible NH 3 cracking reaction based on a microcanonical ensemble within the transition state theory framework, the details of which can be found elsewhere [48,54].The solutions for the system of ordinary differential equations, i.e. the variation of each species with time, were derived by propagating it to equilibrium (Fig. S23).

Temperature-programmed desorption (TPD)
The molecular desorptions from catalysts are crucial to ensure the availability of the active site and turnover; thus, we simulated the temperature-programmed desorption to identify trends and temperatures in which desorption occurs.We found significant differences between the N 2 and H 2 TPD on IrN 3 , RuN 3 , and RhN 3 SACs compared with their extended metals (Fig. 6).The H 2 desorption occurs on metallic Ir(111) at ~300 K, [34] whereas, on IrN 3 -SAC, it desorbs at a temperature below 200 K.Nevertheless, both catalysts present similar N 2 desorption curves expanding from 550-750 K.The N 2 and H 2 TPD curves on RuN 3 -SAC shifted to the lower temperature by ~150 K and   ~25 K with respect to the ones reported on the hcp-Ru(0001) surface.[34] These results prove the effective desorption of SACs compared to traditional metal catalysts as they minimize the lateral interaction with co-adsorbed species.The low desorption temperatures of H 2 (~300 K) and N 2 (~400 K) on the RhN 3 SAC suggest a significant NH 3 cracking conversion at mild temperatures upon NH dehydrogenation.

Batch reactor simulations
We have simulated the ratio between molecular species and active sites on IrN 3 , RuN 3 , and RhN 3 SACs as a function of the temperature and the reaction time.The equilibrium between SAC-adsorbed and gas NH 3 is established at low temperatures until the dehydrogenation and molecular recombination energy barriers are surmounted upon increasing temperatures.The steady-state ratio (χ) of NH 3 , N 2 , and H 2 as functions of temperature after 600 s of reaction is shown in Fig. 7.Note that 600 s is enough to reach equilibrium.
The H 2 evolution occurs between 425 and 550 K, resulting from an initial NH 3 decomposition.On IrN 3 , the content of H 2 in the gas phase rises at a low temperature (~425 K) and quickly reaches a plateau at around 500 K, indicating the saturation of active sites by NH x (x < 3) species.On RuN 3 , the H 2 does not reach an apparent plateau within the temperature range.It is related to the competing H 2 adsorption and the slight increase in the probability of desorbing with the increasing temperature.RhN 3 behaves similarly to IrN 3, although the χ(H 2 ) is ~1.5 lower due to the difficulty of overtaking the endothermic NH 2 dehydrogenation.Despite this, RhN 3 has an earlier (T~ 700 K) but slower H 2 desorption than IrN 3 (T~ 750 K).The main H 2 desorption is followed by N 2 evolution at T~ 725 K (~800 K on IrN 3 ), which competes with molecular adsorptions.N 2 desorption does not arise from RuN 3 , confirming the inability of this SAC to regenerate the catalytic site.Our results predict the NH 3 regeneration at ~500 K, which is in agreement with the NH 3 synthesis experiment on the Ru SACs [42].
The predominant species on the surfaces after 600 s of reaction within a temperature range of 200-1000 K are plotted in Fig. 8.At low temperatures (< 400 K), the adsorbed NH 3 molecule accommodates on the three catalysts.At 500 K, NH 2 becomes the dominant species on the IrN 3 surface, but it is quickly replaced by N as the temperature rises, displaying the slight endothermic path between dehydrogenation steps, see Fig. 4. The strong adsorption of N on IrN 3 hinders the reactivity of the site until the temperature reaches ~800 K, promoting the N 2 evolution.The dominant species on RuN 3 during the reaction are NH 2 and N at 460-520 K and > 520 K, respectively.The N ad-atoms remain on the metal site due to the high energy barrier for their recombination (R7).NH 2 is the predominant species on the RhN 3 -SAC at the temperature range of 500-720 K with a maximum coverage of 0.9 ML at 680 K, related to the endothermic R3 step (~0.75 eV).Once overtaken at ~700 K, molecular N 2 takes most of the metal centres due to its favourable re-adsorption on RhN 3 .It is worth mentioning that although R5 on RhN 3 has the highest activation energy of the energy profile, it is an apparent rate-limiting because NH does not accumulate on the surface at the temperature at which R3 is surmounted.

Conclusion
Single transition metal atoms supported on clean and nitrogenmodified graphene were simulated as MC 3 , MC 4 , MN 3 , and MN 4 (M= Fe, Co, Ni, Ir, Ru, Rh, and Zn).The results show that nitrogen-modified graphene has a more robust interaction with metal atoms than vacancy graphene, providing high resistance to metal sintering and leaching.We investigated the NH 3 cracking reaction mechanism on these SACs supported on N-modified graphene and found that the most likely adsorption site is directly on top of the metal atom.The relative energies suggested that MN 3 structures are more active than MN 4 as the former enhances the electron back-donation, activating the N-H bond.The results confirmed that, on SACs, the first NH 3 dehydrogenation and the N 2 evolution are the rate-determining steps for ammonia cracking.The most promising candidates were those formed by Ir, Ru, and Rh metals.The free energies of molecular and adsorbed species were calculated between 200 and 1000 K for a thermochemical analysis, revealing slight variations of reaction and activation free energies with the temperature.We implemented the free energies in temperature-programmed desorption simulations, finding that desorption temperatures from SACs are smaller than those on extended surfaces.Batch reaction simulations probed IrN 3 and RhN 3 SACs as promising catalysts for NH 3 cracking, being the latest active at lower temperatures than the former.

Declaration of Competing Interest
The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Alberto Roldan reports financial support was provided by Engineering and Physical Sciences Research Council.Xiuyuan Lu reports financial support was provided by China Scholarship Council.
illustrates the relative reaction energies (ΔE) starting from the clean SAC and isolated NH 3 molecule; adsorbed structures are schematically represented in Figs.S8-S21.The thermodynamic profiles show that only ZnN 3 , RuN 3 , RhN 3 , IrN 3 and ZnN 4 present favourable NH 3 adsorption and dehydrogenation.Most MN 3 SACs fail on the last dehydrogenation step, leading to adsorbed N. It differs from the MN 4 SACs, which, except ZnN 4 , are not able to stabilize NH 2 , meaning that the first NH 3 dehydrogenation thermodynamically hinders the reaction progress.The presence of NH 2 on a SAC is only likely on ZnN 3 , although these energies are close on NiN 3 , IrN 3 , RuN 3 , and RhN 3 cases.The NH 2 and NH dehydrogenations on FeN 3 , CoN 3 , and NiN 3 are thermodynamically impeded, explained by the small size of these atoms and the highly localized electrons at the Fermi energy Fig. S3.

Fig. 5 ;
for completeness, ΔG r and ΔG a of the other SACs are represented in Fig. S22.The dehydrogenation transition states (TS1-TS5) are shown in the schematic diagrams of Fig. S8-21.

Fig. 6 .Fig. 7 .
Fig. 6.Simulated N 2 and H 2 TPD on (a) IrN 3 , (b) RuN 3, and (c) RhN 3 at different initial atomic coverages (θ in ML) with the heating rate of 1 K/min.Fig. 7.The steady-state ratio (χ) of molecular NH 3 , N 2 , and H 2 as functions of temperature upon NH 3 cracking catalyzed by (a) IrN 3 , (b) RuN 3 , and (c) RhN 3 on batch reactor simulations.The initial simulated conditions are an NH 3 ratio of 5:1 with the surface and a reaction time of 600 s.

Table 2
Electronic details of modelled MN 3 and MN 4 SACs, including the metal Bader charge (q) and spin density (s), and the total surface magnetic moment (M s ).