Two-dimensional metal–organic TMTAP monolayers as a promising class of monoatomic electrocatalysts for oxygen reduction, oxygen evolution, and hydrogen evolution reactions

Abstract

Regenerative fuel cells, water splitting, and metal-air batteries all require high-performance electrocatalysts to drive the oxygen reduction reaction (ORR), oxygen evolution reaction (OER), and hydrogen evolution reaction (HER). Single-atom catalysts (SACs) have drawn considerable attention due to their superior catalytic properties. In this work, we evaluate 30 transition metals anchored to tetraazaporphyrin (TAP) monolayers as SACs for ORR, OER, and HER using density functional theory (DFT). Among the candidates operating via the 4e⁻ pathway, Zn-TAP (ORR overpotential: 0.43 V) and Rh-TAP (OER overpotential: 0.58 V) exhibit top-tier performance. Fe-TAP also shows good ORR activity (0.58 V). For the 2e⁻ pathway, Ag-TAP delivers overpotentials of 0.36 V for ORR, making it promising for H₂O₂ synthesis. In the 2e⁻ OER, the Cr, Mn, Fe, Co, and Zn systems show very low overpotentials (down to 0.17 V), albeit with poor 2e⁻ ORR selectivity. Regarding HER, Cd- and Os-TAP are identified as promising candidates, with near-ideal hydrogen adsorption free energies of −0.08 and −0.12 eV, respectively. Further mechanistic insights into O₂ activation are provided through analyses of spin density, charge density difference, and electronic structure. We note that this study is based on the computational hydrogen electrode (CHE) model under idealized acidic conditions (pH = 0) without explicit or implicit solvent effects. Therefore, the calculated overpotentials and free-energy profiles represent a thermodynamic screening framework, and extension to neutral/alkaline or experimental environments requires additional consideration of solvation and pH effects. Overall, this work highlights the potential of TMTAP SACs and offers guidance for the rational design of such catalysts.

---

1. Introduction

The growing global demand for energy and increasing environmental concerns have intensified the pursuit of green and sustainable energy technologies.¹ Although fossil fuels have historically driven industrialization, their widespread use leads to pollution, climate change, and resource depletion, underscoring the need for cleaner alternatives. Renewable sources such as solar, wind, and tidal energy hold promise, yet they face challenges including geographic limitations, intermittency, and instability, which restrict their ability to meet rising energy demands. Consequently, electrochemical energy conversion technologies, such as water electrolysis and fuel cells, have emerged as critical research frontiers. Fuel cells electrochemically convert the chemical energy of fuels like hydrogen into electricity, offering high efficiency, low emissions, and a compact design.^2,3 However, their commercialization is hindered by expensive catalysts, durability issues, and hydrogen storage challenges. Water electrolysis provides a clean route to high-purity hydrogen, supporting a future hydrogen economy with zero emissions. Both technologies rely on key electrochemical processes: the ORR, OER, and HER. The efficiency of these reactions largely dictates overall device performance.^4–6 Therefore, developing efficient, stable, and cost-effective electrocatalysts is essential for advancing hydrogen fuel cells, metal–air batteries, and electrolyzers, and is central to enabling the green energy transition.

Currently, Pt-based materials remain the benchmark electrocatalysts for ORR and HER, while Ir and Ru oxides are the standards for OER.^7–12 Nevertheless, their high cost, limited abundance, and tendency to dissolve and deactivate during prolonged operation lead to poor stability, hindering large-scale adoption. In renewable energy systems such as fuel cells and water electrolyzers, the cost, activity, and durability of electrocatalysts are critical to practical feasibility.^13–16 Hence, developing electrocatalysts with lower cost, higher efficiency, enhanced multifunctionality, and superior stability remains a key research goal. In line with the recent progress in polymetallic nanostructures, Yang et al.¹⁷ developed a library of binary to septenary polymetallic alloy nanotubes, demonstrating that the downshift of the d-band center in high-entropy alloys (e.g., Cu₃₀Ni₂₆Co₁₉Ru₁₄Ir₁₁) can significantly enhance both OER and HER catalytic performance in alkaline media. For acidic water splitting, Zhao et al.¹⁸ reported porous rare-earth-based high-entropy perovskite oxide nanosheets (RE-HEPN) with abundant oxygen vacancies and low-valence active sites, achieving ampere-level current density stability for over 1000 h in a PEMWE cell. Furthermore, He et al.¹⁹ systematically reviewed lanthanum-based OER catalysts, highlighting the critical roles of atomic doping, interfacial modulation, and defect engineering in optimizing the electronic structure and reaction kinetics. These pioneering studies clearly indicate that entropy stabilization, rare-earth incorporation, and defect chemistry are effective strategies to boost electrocatalytic water splitting. In recent years, SACs have attracted growing interest due to their exceptional atomic utilization efficiency, high selectivity, and tunable surface properties.^20–22 The remarkable catalytic performance of SACs largely originates from the coordination mismatch between the metal atoms at active sites and the supporting substrate.^22,23 Guided by this principle, researchers have designed various catalyst architectures and achieved promising progress by tuning the chemical environment of metal centers, including their local structure and coordination number, through conventional strategies such as doping, defect engineering, strain modulation, and chemical functionalization.^24–29 Numerous experimental and theoretical studies have confirmed that rationally designed two-dimensional (2D) single-atom catalysts can achieve superior electrocatalytic performance and optimal metal utilization through the targeted selection of metal centers and their coordination environments.^20,30,31

Inspired by these advances, the present study aims to evaluate the catalytic efficacy of 2D single-atom transition metal tetraazaporphyrin (TMTAP) materials for the ORR, OER, and HER. To this end, we anchored TM atoms from the 3d, 4d, and 5d series onto a chemically modified TAP substrate to construct TMTAP SACs. The TAP substrate has a porphyrinoid structure in which the four meso carbon atoms of the porphyrin ring are replaced by nitrogen atoms, resulting in a framework composed of four pyrrole rings and four nitrogen atoms. Experimental studies have demonstrated the feasibility of synthesizing single-atom catalysts derived from metallomacrocyclic compounds. For instance, Joo reported carbon-supported atomically dispersed M–N/C catalysts (M = Fe, Co, Ni, Pt) using metallomacrocycle precursors, achieving high activity and selectivity for ORR, CO₂RR, and CER.³² Similarly, Cao et al. employed a supramolecular macrocycle (pillar[4]arene[1]quinone) to regulate single-atom Co sites on MoS₂, yielding excellent OER performance.³³ Importantly, the synthesis of tetraazaporphyrin (TAP) macrocycles has been well-established since the pioneering work of Ficken et al. in 1958, with modern synthetic strategies enabling the isolation of stable hydrogenated TAP derivatives.^34,35 These experimental precedents collectively support the synthetic feasibility of the TMTAP catalysts proposed in our theoretical screening study. The potential of TMTAP SACs as efficient atomic-level catalysts for ORR, OER, and HER was systematically evaluated through DFT-based computational screening. Through this screening process, promising candidates were identified, and their reaction mechanisms and intrinsic material properties were investigated. This work introduces a novel TAP-based bifunctional oxygen electrocatalyst and provides a theoretical foundation for the design of multifunctional SACs on 2D platforms.

2. Computational methodologies

All spin-polarized calculations were performed using DFT as implemented in the Vienna Ab initio Simulation Package (VASP).³⁶ The projector-augmented wave (PAW) method³⁷ was used to describe the interactions between ion cores and valence electrons, with an energy cutoff of 450 eV. The exchange–correlation functional was treated within the generalized gradient approximation (GGA) using the Perdew–Burke–Ernzerhof (PBE) formulation.³⁸ A vacuum region of more than 15 Å was added along the z-axis to avoid artificial interlayer coupling arising from periodic boundary conditions. Energy and force convergence thresholds were set to 10⁻⁵ eV per atom and 0.02 eV Å⁻¹, respectively. For structural relaxation and electronic structure analysis, Γ-centered Monkhorst–Pack³⁹ k-point meshes with a resolution of 2π × 0.025 Å⁻¹ (a 5 × 5 × 1 mesh) were employed for Brillouin zone integration. K-point sampling and energy cutoff were both tested for convergence. Based on the convergence tests, we estimate the numerical uncertainty in our calculated free energies to be approximately ±0.03–0.05 eV. This translates to an uncertainty of roughly ±0.03–0.05 V in the derived overpotentials. Based on convergence tests for k-point sampling and energy cutoff, the numerical uncertainty in our calculated free energies is estimated to be ±0.05 eV, leading to an uncertainty of approximately ±0.05 V in the reported overpotentials. Van der Waals (vdW) interactions were included using Grimme's D3 semiempirical dispersion correction.⁴⁰ The thermal stability of the TMTAP catalysts was examined via ab initio molecular dynamics (AIMD) simulations performed on a 2 × 2 × 1 supercell at 500 K for 10 ps with a time step of 2 fs; temperature was controlled using a Nosé–Hoover thermostat.⁴¹ Bader charge analysis⁴² was carried out to gain insights into charge transfer and distribution characteristics.

The tendency of TM atoms to aggregate on the TAP substrate was evaluated by calculating cluster energies (E_clus = E_b + E_c). Here, E_b denotes the binding energy, defined as E_b = E_TMTAP − E_TAP − E_TM-single, and E_c represents the cohesive energy, given by E_c = E_TM-single − E_TM-bulk/N. In these expressions, E_TMTAP and E_TAP are the total energies of the TM-embedded TAP system and the pristine TAP substrate, respectively; E_TM-single is the energy of an isolated TM atom in vacuum; and E_bulk is the energy per TM atom in its stable bulk phase, with N indicating the number of atoms in the bulk model. To further assess system stability, dissolution potentials of the TM atoms were computed to evaluate their likelihood of leaching from the TAP substrate.^43,44 In this context, refers to the standard dissolution potential of the corresponding bulk metal, and n denotes the number of electrons transferred during the dissolution process.

The Gibbs free energy change (ΔG) for each elementary electrochemical step was calculated using ΔG = ΔE + ΔE_ZPE − TΔS + ∫C_pdT + ΔG_U + ΔG_pH within the CHE framework. In this expression, ΔE represents the DFT-derived reaction energy, T is the temperature (298.15 K), and ΔE_ZPE and ΔS represent the differences in zero-point energy and entropy, respectively, between the adsorbed species and their gas-phase counterparts; ∫C_pdT is the enthalpy correction. Values of ΔE_ZPE, ΔS, and ∫C_pdT for gaseous molecules were taken from the NIST database,⁴⁵ while those for adsorbed intermediates were extracted from vibrational frequency calculations. The contribution of an applied electrode potential U is given by ΔG_U = neU, with n being the number of electrons transferred. The pH-dependent free energy correction follows ΔG_pH = k_BT × ln 10 × pH, where k_B is the Boltzmann constant; pH is taken as zero in this study.

The elementary reaction step with the largest Gibbs free energy change (ΔG_max) in the proton–electron transfer process is identified as the potential-determining step (PDS). The corresponding limiting potential (U_L) or onset potential (U_onset), calculated as U_L = −ΔG_max/e, provides a measure of the ORR catalytic overpotential (η_ORR). Meanwhile, HER performance was assessed using the hydrogen adsorption free energy (ΔG_^*H), computed with the CHE model referenced to the reversible hydrogen electrode (RHE).⁴⁶ An optimal HER catalyst corresponds to ΔG_^*H = 0, with smaller absolute values |ΔG_^*H| indicating enhanced activity. The corresponding theoretical HER overpotential (η_HER) is given by η_HER = −|ΔG_^*H|/e.

3. Results and discussion

3.1. Structure and stability of catalysts

As shown in Fig. 1a, the TMTAP single-atom catalyst is based on a two-dimensional TAP monolayer, whose structural unit consists of four pyrrole rings forming a stable macrocyclic ligand. A key feature of this structure is that each pyrrole ring provides a pyrrolic nitrogen atom, which, together with its lone pair electrons, serves as a coordination site to anchor a transition metal ion. This multidentate coordination mode effectively immobilizes the metal center, preventing its migration and aggregation, and ensures that transition metals are uniformly dispersed as isolated single atoms, thereby providing an ideal platform for constructing high-performance SACs. By embedding 3d, 4d, and 5d transition metals into the TAP macrocyclic cavity, we systematically constructed a series of single-atom systems. Unlike the confined cavity of g-C₄N₃,⁴⁷ the cavity size of TAP allows most TM atoms to achieve planar coordination with the coordinating nitrogen atoms. However, Sc-, Ti-, Y-, Zr-, Nb-, Cd-, Lu-, Hf-, and Ta-TAP are clear exceptions, as their central metal atoms protrude noticeably out of the plane, resulting in puckered geometries (see Table S1). After TM incorporation, the lattice constants of TMTAP show only minor changes, indicating the excellent structural rigidity of the substrate framework. The TM–N bond lengths vary considerably with the metal species; 3d metals generally exhibit shorter bond lengths than their 4d/5d counterparts (e.g., Mn–N at 1.96 Å vs. Mo–N at 2.05 Å), which arises from atomic radius contraction due to differences in d-orbital electron occupancy. For ORR/OER, moderate TM–N bond lengths (e.g., 1.93–2.00 Å for Fe and Zn) facilitate optimized O₂ adsorption geometries, promoting O–O bond cleavage and favorable adsorption/desorption kinetics of intermediates (*OOH, *O, *OH).

Fig. 1 (a) Top and side views of SACs constructed by coordinating 30 TM atoms (TM = 3d, 4d, and 5d transition metal atoms) with the TMTAP substrate. (b) Binding energy E_b (purple), cohesive energy E_c (orange), and E_b + E_c (red) for TMTAP. (c) Dissolution potentials U_diss of metal atoms. A stable TMTAP system requires E_b + E_c < 0 and U_diss > 0.

The magnetic moment data reveal that, except for closed-shell metals such as Zn and Rh, most TMTAP catalysts retain unpaired d electrons (Table S1). The high magnetic moments in catalysts like CrTAP (3.84 µ_B) and WTAP (3.58 µ_B) originate from unquenched atomic magnetic moments under the influence of the coordination field, which may optimize the adsorption energies of oxygen-containing intermediates by enhancing spin polarization effects. Notably, the magnetic moments of catalysts such as VTAP (2.68 µ_B) and FeTAP (2.00 µ_B) are closely related to the charge transfer capability of oxygen intermediates, directly influencing the ORR/OER overpotentials.

Energy stability analysis (Fig. 1b) demonstrates that TMTAP catalysts possess excellent structural stability. E_b calculations yield negative values (ranging from −0.26 eV to −9.75 eV, see Table S2) for all catalysts, indicating strong chemical bonding between the TM atoms and the substrate, which effectively suppresses metal atom migration. According to the criterion E_b + E_c < 0 (where E_c is the cohesive energy of the bulk metal), most TMTAP catalysts satisfy the condition for single-atom dispersion, except for noble metals such as Re, Au, and Hg. For instance, FeTAP (E_b = −7.52 eV, E_c = 4.63 eV) achieves a total energy advantage of −2.89 eV, demonstrating its excellent anti-aggregation properties as an SAC, which is crucial for maintaining thermodynamic stability in electrocatalytic processes. Electrochemical stability was further evaluated by analyzing the U_diss (Fig. 1c). The results indicate that 23 catalysts are stable within the electrochemical window (U_diss > 0). Among them, RhTAP (U_diss = 1.77 V) exhibts strong resistance to dissolution, suggesting its potential applicability in acidic ORR environments. By simultaneously satisfying the criteria E_b + E_c < 0 and U_diss > 0, we identified 20 TMTAP catalysts (TM = Ti, V, Cr, Mn, Fe, Co, Ni, Cu, Zn, Mo, Tc, Ru, Rh, Pd, Ag, Cd, W, Os, Ir, and Pt) for further study.

To evaluate the thermal stability of the TMTAP SACs, we performed AIMD simulations at 500 K for 10 ps on the 20 pre-screened catalysts (see Fig. 2 and S1–S3). Throughout the simulations, the C–N framework within the unit cell showed no observable configurational distortion; all atoms exhibited localized thermal vibrations near their equilibrium positions, indicating excellent rigidity of the substrate structure. The TM atoms remained stably anchored at the center of the triazine cavities without any migration or aggregation, confirming the strong constraining effect of the TM–N coordination bonds. Both the total energy and temperature of the catalysts fluctuated within narrow ranges, confirming the preservation of structural integrity. These results demonstrate that the 20 TMTAP catalysts possess favorable thermal stability, providing theoretical support for their potential application under ambient conditions. Based on these findings, we proceeded to investigate the catalytic performance of these 20 kinetically stable catalysts for ORR/OER and HER, aiming to identify highly active SACs. It should be noted that the 10 ps AIMD simulation primarily probes resistance to thermal perturbation rather than electrochemical degradation mechanisms (e.g., dissolution, surface oxidation, or restructuring under applied potential). Therefore, these computational results should be interpreted as a preliminary stability filter rather than a definitive prediction of practical durability. The catalyst shows promising short-term thermal stability in AIMD simulations, serving as an initial computational screening indicator. However, its long-term electrochemical stability under operating conditions requires further experimental validation.

Fig. 2 Evolution of temperature and total energy during AIMD simulations for (a) FeTAP, (b) ZnTAP, (c) RhTAP, and (d) AgTAP at 500 K for 10 ps with a time step of 2 fs. The insets show snapshots of the initial and final frames of the AIMD simulations.

It is important to note that while our DFT calculations were performed under the CHE framework at pH = 0, many SACs are known to undergo demetallation in strongly acidic environments due to proton attack and metal leaching. The stability metrics we employed, including binding energies, anti-aggregation criteria, dissolution potentials, and short-term AIMD simulations, do not directly capture the kinetics or thermodynamics of demetallation under acidic operating conditions. For example, a low U_diss relative to the operating potential may indicate susceptibility to metal leaching. In this work, we have calculated U_diss for each TMTAP catalyst (see Table S2), and while some candidates show promising dissolution potentials, these values should be interpreted as thermodynamic trends rather than guarantees of acid stability. Experimental validation under acidic conditions is essential to confirm resistance to demetallation.

3.2. ORR and OER electrocatalytic activity

We systematically analyzed the thermodynamic pathway preferences and kinetic activity characteristics of the 20 screened SACs in 4e⁻ and 2e⁻ oxygen electrocatalytic reactions. Thermodynamic analysis shows that the electronic structure of the metal center imposes a selection on the reaction pathway via a ΔG_^*O threshold of 3.46 eV. The threshold is derived by equating the Gibbs free energies of the 2e⁻ pathway (O₂ → H₂O₂) and the 4e⁻ pathway (O₂ → H₂O). For a constant O intermediate, the thermodynamic overpotential for each pathway can be expressed as a function of ΔG_^*O. Solving ΔG_(2e⁻) = ΔG_(4e⁻) gives a critical ΔG_^*O value of 3.46 eV (ΔG_H2O2 − ΔG_H2O), below which the 4e⁻ pathway is thermodynamically favored. This derivation assumes all other steps scale linearly with ΔG_^*O and that the potential-determining steps are *OOH formation (2e⁻) and *O hydrogenation (4e⁻). It should be noted that the 3.46 eV threshold is purely thermodynamic. In practice, kinetic barriers, such as those associated with O–O bond cleavage or *OOH protonation, may alter the observed selectivity, especially near the threshold. Moreover, solvent effects, pH, and electrode potential can shift the effective free energies. Therefore, this value serves as a qualitative guide rather than a strict boundary. Specifically, Ni-, Cu-, Pd-, Ag-, and Pt-TAP exhibit ΔG_^*O values exceeding this threshold (Fig. S4); as a result, their *OOH intermediates preferentially follow the non-dissociative 2e⁻ pathway toward H₂O₂ production, owing to the precise modulation of adsorption strength through metal-support electronic interactions. Conversely, metal centers with ΔG_^*O below the threshold favor *OOH desorption to complete the 4e⁻ reduction pathway to H₂O.

Overpotential analysis reveals differences in reaction energy barriers among metal centers, as illustrated in Fig. 3. In the 4e⁻ pathway (Fig. 3a), ZnTAP (η_ORR = 0.43 V) and FeTAP (η_ORR = 0.58 V) exhibit highly efficient ORR activity by optimizing O₂ activation and proton-coupled electron transfer (PCET) steps. Their ORR overpotentials are comparable to or better than those of conventional noble-metal Pt (η_ORR = 0.45 V)⁴⁸ and recent reported SACs such as MnPor (η_ORR = 0.56 V),⁴⁹ Fe-Pc (η_ORR = 0.45 V),⁵⁰ Mn₃(HHTP)₂ (η_ORR = 0.44 V),⁵¹ Rh₃(HADQ)₂ (η_ORR = 0.31 V),⁵² and N-doped armchair graphene nanoribbon (η_ORR = 0.41 V),⁵³ and far superior to that of Cu@g-C₄N₃ (η_ORR = 0.70 V)⁵⁴ and N-doped graphene (η_ORR = 0.67 V),⁵⁵ indicating good ORR catalytic performance. Notably, ZnTAP combines a low ΔG_^*O (thermodynamic preference for the 4e⁻ pathway) with the lowest ORR overpotential (0.43 V), reflecting a synergistic optimization of electronic structure for both thermodynamic stability and reaction kinetics. RhTAP (η_OER = 0.58 V) demonstrates kinetic advantages in 4e⁻ OER due to the low energy barriers associated with *O desorption and O–O bond reconfiguration. Its OER overpotential is comparable to those of conventional metal oxide catalysts and 2D metal–organic SACs, such as RuO₂ (η_OER = 0.42 V),⁵⁶ IrO₂ (η_OER = 0.56 V),⁵⁶ Mn₃(HHTP)₂ (η_OER = 0.53 V),⁵¹ and Rh₃(HADQ)₂ (η_OER = 0.31 V),⁵² as shown in Table S3. For the 2e⁻ pathway (Fig. 3b), AgTAP exhibits ultralow ORR and OER overpotentials (0.36 V), indicating high intrinsic catalytic efficiency of the single-atom sites for H₂O₂ production via *OOH hydrogenation. Thus, AgTAP is a promising and efficient bifunctional ORR and OER catalyst. Additionally, although Cr-, Mn-, Fe-, Co-, and Zn-TAP show poor 2e⁻ ORR selectivity, they exhibit excellent 2e⁻ OER catalytic activity with low 0verpotentials of 0.49, 0.34, 0.39, 0.21, and 0.17 V, respectively. These results demonstrate that the oxygen catalytic performance of TMTAP catalysts is synergistically governed by the ΔG_^*O threshold and reaction overpotentials, and the metal-center specificity providing critical theoretical guidance for the rational design of high-efficiency oxygen reduction/evolution catalysts.

Fig. 3 (a) Overpotentials for the 4e⁻ ORR and OER for the 20 stable TM-TAP SACs. (b) Overpotentials for the 2e⁻ ORR for the 20 stable TM-TAP SACs.

The adsorption free energies of *OOH (ΔG_^*OOH), *O (ΔG_^*O), and *OH (ΔG_^*OH) on the TMTAP systems were calculated and are compiled in Table S4. To investigate the electrochemical activity trends of the TMTAP SACs, we constructed volcano plots by analyzing the linear relationships between ORR/OER overpotentials (η) and the adsorption free energies of reaction intermediates (ΔG_^*OH for η_ORR and ΔG_^*OOH for η_OER), as shown in Fig. 4. The results indicate that Zn- and Rh-TAP lie near the peaks of the ORR and OER volcano plots, respectively, while Fe- and Ni-TAP also approach the volcano summits. This distribution trend aligns with the overpotential calculations for the 4e⁻ pathway (Fig. 3), where Zn- and Rh-TAP show the lowest ORR (η_ORR = 0.43 V) and OER overpotentials (η_OER = 0.58 V), respectively. The study validates that the linear scaling relationships between ΔG_^*OH–η_ORR and ΔG_^*OOH–η_OER, captured by the volcano plots, serve as effective thermodynamic descriptors for quantitatively evaluating the intrinsic activity of TMTAP systems, thereby providing a theoretical paradigm for the rational design of SACs.

Fig. 4 Volcano plots of (a) η_ORR versus ΔG_^*OH and (b) η_OER versus ΔG_^*OOH.

To systematically elucidate the reaction mechanisms of the TMTAP SACs, we computed free energy diagrams for the 4e⁻ ORR/OER pathways on Fe-, Zn-, and Rh-TAP, as well as the 2e⁻ OER pathway on AgTAP, as depicted in Fig. 5. Taking the ORR process on ZnTAP as an example (Fig. 5b), the reaction follows a stepwise mechanism: * + O₂ → *OOH → *O → *OH → * + H₂O. In the initial stage, O₂ undergoes hydrogenation adsorption on ZnTAP, forming an *OOH intermediate (ΔG = −0.90 eV), accompanied by elongation of the O–O bond from 1.21 Å in free O₂ to 1.46 Å, while a Zn–O bond (1.96 Å) is simultaneously formed. Subsequently, under the synergistic action of proton–electron pairs (H⁺/e⁻), *OOH dissociates into *O and H₂O (ΔG = −0.96 eV), with the Zn–O bond in the *O species shortening to 1.88 Å. Further hydrogenation drives the conversion of *O to *OH (ΔG = −2.27 eV), and final protonation of *OH releases the second H₂O molecule (ΔG = −0.80 eV). According to the criterion that the PDS is the step with the largest ΔG governing the reaction rate, the *OH → H₂O transformation is identified as the PDS for ZnTAP, corresponding to an overpotential (η_ORR) of 0.43 V, which is superior to that of the Pt(111) catalyst (η = 0.45 V).⁴⁸ Comparative analysis of other systems reveals that the PDS for FeTAP is the *O → *OH hydrogenation step (η = 0.58 V), whereas RhTAP suffers a significant activity loss (η = 1.20 V) due to the high energy barrier (ΔG = 1.39 eV) in the *OOH formation step, underscoring the critical role of metal-center electronic structures in modulating intermediate adsorption energies.

Fig. 5 Free energy diagrams for the 4e⁻ ORR on (a) FeTAP and (b) ZnTAP at zero potential (pink), limiting potentials (orange), and equilibrium potential (blue). Free energy diagrams for the 4e⁻ OER on (c) RhTAP catalyst. Free energy diagrams for the 2e⁻ ORR on (d) AgTAP. The PDS and η are indicted in the diagrams.

To elucidate the influence of electrode potential on the 4e⁻ ORR process, we constructed free energy evolution profiles for the ZnTAP catalyst under open-circuit potential, limiting potential and equilibrium potential conditions (Fig. 5b). The analysis reveals that at open-circuit potential, all elementary ORR steps in the ZnTAP catalytic system exhibit a downward free energy trend (ΔG < 0), indicating an overall thermodynamic driving force for the reaction. However, localized energy barriers in certain intermediate transformation steps significantly narrow the operative potential window. The limiting potential for this reaction pathway is determined to be U_L = 0.80 V; when the applied potential exceeds this threshold, the free energy change for the *OH protonation step (the fourth step of ORR) becomes positive (ΔG > 0), transitioning from exothermic to endothermic. Under standard equilibrium potential (U = 1.23 V), the formation step of the second H₂O molecule (*OH → * + H₂O) requires overcoming an energy barrier of 0.43 eV, hindering ORR kinetics, while subsequent desorption processes remain thermodynamically spontaneous (ΔG < 0). This characteristic unambiguously identifies the *OH protonation step as the ORR PDS for the ZnTAP catalyst, corresponding to an overpotential of 0.43 eV.

As the reverse reaction of ORR, the OER exhibits opposite free energy variation characteristics. As illustrated in Fig. 5a–c, all elementary steps in the 4e⁻ OER pathways over Fe-, Zn-, and Rh-TAP shows an upward free energy trend, indicating that water splitting is endothermic and requires continuous external energy input. The calculation results demonstrate that for FeTAP, the PDS is the *O → *OOH transformation, with a free energy change of ΔG = 1.39 eV, corresponding to a theoretical overpotential (η_OER) of 1.16 V. For ZnTAP, the PDS occurs during the *OH → *O dehydrogenation step (ΔG = 2.27 eV), giving η_OER = 1.04 V. RhTAP exhibits the best OER activity; its PDS also involves the *O → *OOH conversion (ΔG = 1.81 eV), but owing to optimized adsorption of intermediates at the metal center, the overpotential is significantly reduced to η_OER = 0.58 V. These disparities highlight the differential modulating effects of distinct metal centers (Fe/Zn/Rh) on the TAP substrate through electronic structure regulation of OER intermediate adsorption/desorption behaviors, providing a theoretical foundation for the rational design of high-efficiency oxygen evolution catalysts.

Fig. 5d illustrates the free energy evolution characteristics of the AgTAP catalyst along the overall 2e⁻ ORR pathway. Similar to the thermodynamic response observed in the 4e⁻ ORR process, each elementary step of the 2e⁻ pathway (* + O₂ → *OOH → * + H₂O₂) exhibit a downward free energy trend (ΔG < 0), indicating an overall thermodynamic driving force for the reaction. Analysis of the PDS under open-circuit potential conditions reveals that the O₂ to form the *OOH intermediate (ΔG = −0.37 eV) is the rate-limiting step, corresponding to a theoretical overpotential (η_ORR) of 0.36 V. Under limiting potential conditions (U_L = 0.37 V), the free energy profile retains its downward trajectory, localized energy barriers in certain intermediate transformation steps significantly narrow the operative potential window. When the applied potential exceeds this threshold, the free energy change for the *OOH protonation step becomes positive (ΔG > 0), transitioning from exothermic to endothermic. Furthermore, under equilibrium potential conditions (U = 0.73 V), the free energy changes of all elementary steps become positive (ΔG > 0), indicating that this pathway is non-spontaneous under thermodynamic equilibrium. In summary, the ZnTAP catalyst exhibits optimal intrinsic activity in the 4e⁻ ORR with η_ORR = 0.43 V, while the RhTAP catalyst shows the lowest overpotential in the 4e⁻ OER (η_OER = 0.58 V), highlighting their excellent catalytic potential for ORR/OER. The AgTAP catalyst, owing to the thermodynamic advantage in 2e⁻ ORR pathway, shows promise for selective H₂O₂ production. This differential distribution of catalytic performance originates from the precise modulation of reaction intermediate adsorption/activation capabilities by distinct metal centers.

3.3. HER electrocatalytic performance

The HER pathway is typically described using a three-state diagram comprising the initial state (H⁺ + e⁻), the intermediate state (adsorbed *H), and the final product (H₂). The intrinsic HER activity of catalysts can be effectively evaluated using the thermodynamic descriptor (ΔG_^*H), the hydrogen adsorption free energy at the intermediate state. According to the Sabatier principle, optimal HER catalysts should exhibit ΔG_^*H values close to zero to balance hydrogen adsorption and desorption on the catalyst surface. Deviations from this ideal range (either excessively negative or positive ΔG_^*H) lead to diminished catalytic activity. In this study, we calculated ΔG_^*H for hydrogen adsorption on the 20 stable TMTAP catalysts, and the corresponding free energy diagrams are presented in Fig. 6a. Notably, the CdTAP (−0.08 eV) and OsTAP (−0.12 eV) systems show ΔG_^*H values closest to the theoretical optimum, indicating high potential for HER catalysis. Further analysis of the relationship between η_HER and ΔG_^*H (Fig. 6b) reveals a distinct volcano-type relationship for the TMTAP catalysts, with CdTAP and OsTAP located near the summit of the volcano plot. This positioning signifies their exceptional HER catalytic performance, which is superior or comparable to that of the Pt(111) surface catalyst (η_HER = 0.09 V),⁴⁶ Co@C₃N₂ (η_HER = 0.12 V)⁵⁷ and Mn₂C₁₂ (η_HER = 0.09 V) SACs.⁵⁸

Fig. 6 (a) Free energy diagrams for HER at equilibrium potential for the 20 TMTAP monolayers. (b) Relationship between HER overpotentials (η_HER) and ΔG_^*H for the 20 TMTAP monolayers. The η_HER of Pt(111) is shown in the diagram for comparison of HER catalytic activity.

3.4. Origin of catalytic activity

The adsorption behavior of O₂ on catalyst surfaces is a critical process in the initial steps of both ORR and OER. Fig. 7 systematically illustrates the spin density distributions and charge density differences for Fe-, Zn-, Rh-, and Ag-TAP SACs before and after O₂ adsorption, providing intuitive electronic-structure evidence to elucidate the adsorption mechanism. Before O₂ adsorption (Fig. 7a–d), systems such as Fe- and Rh-TAP exhibit pronounced spin polarization (red/green distributions) at the TM sites, indicating that unpaired d electrons occupy specific orbitals. After O₂ adsorption (Fig. 7a′–d′), the spin density distributions undergo reconfiguration: (1) spin-polarized regions extend from the isolated TM atoms toward the TM–O₂ coordination structure; (2) in certain systems (e.g., FeTAP), an extended spin density distribution appears over the O₂ molecule, suggesting involvement of the O₂ π orbitals in electron reorganization. This spin redistribution pattern aligns with the magnetic moment calculations (Table S5), revealing spin-dependent electron transfer between TM and O₂ during adsorption, where pairing effects between unpaired electrons in TM d orbitals and O₂ π orbitals significantly alter the system's spin state.

Fig. 7 Top and side views of spin density for (a) FeTAP, (b) ZnTAP, (c) RhTAP, and (d) AgTAP. (a′–d′) Top and side views of spin density after O₂ adsorption for Fe-, Zn-, Rh-, and Ag-TAP, respectively. Red and green regions denote spin-up and spin-down states, respectively. (e–h) Top and side views of charge density difference for Fe-, Zn-, Rh-, and AgTAP. Pink and light green regions represent charge accumulation and depletion zones, respectively.

Comparison of charge redistribution before and after adsorption (Fig. 7e–h) reveals a distinct pink electron accumulation zone at the TM–O bonding region, while the O–O bond region shows a light green electron depletion zone. This charge displacement pattern indicates: (1) the formation of strongly polarized TM–O bonds with electron transfer from TM to the O₂ antibonding orbitals; (2) reduced electron density at the O–O bond, which weakens bond order, confirming the activation mechanism of adsorbed O₂. This electron reorganization effect directly accounts for the elongation of the O–O bond length, creating favorable conditions for subsequent proton-coupled electron transfer steps.

Additionally, we selected Ti, Mo, and W metals with lower chemical reactivity for comparison with the more reactive Fe, Zn, Rh, and Ag. We calculated the Bader charges on the metal and O atoms, as well as the magnetic moment of O atoms after O₂ adsorption, for Fe-, Zn-, Rh-, Ag-, Ti-, Mo-, and W-TAP catalysts, as summarized in Table S6. The results indicate that, compared to Fe-, Zn-, Rh-, and Ag-TAP catalysts, O₂ molecules adsorbed on Ti-, Mo-, and W-TAP acquire more electrons from the metals. This leads to an increased number of paired electrons in the O-2p orbitals and a reduction in magnetism. Consequently, the interaction between oxygen and Ti, Mo, or W becomes stronger, occupying the active sites and hindering further participation of other reactants, ultimately resulting in lower catalytic activity.

To investigate the electronic structural features of the screened catalysts, we systematically analyzed the band structures and projected density of states (PDOS) of four systems (Fe-, Zn-, Rh-, and Ag-TAP) before and after O₂ adsorption. Band structure analysis (Fig. 8a–c) reveals that all four catalysts exhibit metallic characteristic, with continuous energy bands at the Fermi level for both spin-up and spin-down channels. This feature facilitates electron transport during catalytic processes. PDOS analysis (Fig. 8a′–a″, b′–b″, c′–c″ and d–d″) shows that before O₂ adsorption, significant hybridization occurs between the transition metal d orbitals and the 2p orbitals of coordinating pyrrolic nitrogen atoms near the Fermi level, indicating stable bonding interactions between the metal center and the nitrogen coordination environment. After O₂ adsorption, notable orbital overlap is observed between the 2p orbitals of adsorbed O₂ and the metal d orbitals near the Fermi level. This electronic structural characteristic provides an effective pathway for charge transfer between O₂ molecules and the active sites, thereby promoting the activation process of O₂ molecules and offering theoretical insights into their enhanced ORR catalytic efficiency.

Fig. 8 Spin-resolved band structures and PDOS before and after O₂ adsorption on the active TM atom for (a–a″) FeTAP, (b–b″) ZnTAP, (c–c″) RhTAP, and (d–d″) AgTAP. The Fermi level is set to 0 eV.

To gain further insight into the origins of OER/ORR/HER activity in TMTAP, we computed the d-band centers (ε_d) of the metal active sites. Fig. S5a and S5b reveal a volcano-shaped relationship between the ε_d values of certain metals and the ORR/OER overpotentials. This trend suggests that an optimal d-band center, positioned neither too close to nor too far from the Fermi level, signals a balanced interaction between the metal site and adsorbates, thereby promoting high activity in both ORR and OER. In the volcano plots, metals such as Fe (ε_d = −1.05 eV), Zn (ε_d = −1.39 eV), and Rh (ε_d = −1.29 eV) exemplify this favorable intermediate positioning. Additionally, the Bader charges of the metal active sites were calculated and are presented in Fig. S5c. Most metals exhibit a strong linear correlation with their Bader charge values. A high Bader charge corresponds to weak bonding between the metal and H atoms, which promotes hydrogen evolution. In contrast, a low Bader charge indicates stronger TM–H bonding, favoring hydrogen adsorption. Thus, for optimal hydrogen adsorption and desorption, the Bader charge should be intermediate (neither too high nor too low), placing it near the volcano peak, as illustrated by Cd (−0.90) and Os (−0.96) in the volcano plot.

4. Conclusions

In summary, we have performed systematic DFT calculations to explore the potential of single transition-metal atoms embedded in TMTAP monolayers for electrocatalyzing the ORR, OER, and HER. Among all candidates, Fe-, Zn-, and Rh-TAP emerge as promising ORR/OER catalysts via a 4e⁻ pathway, demonstrating high stability and low overpotentials. In particular, ZnTAP and RhTAP achieve the best catalytic performance with the lowest overpotentials (η_ORR = 0.43 V, η_OER = 0.58 V), which are comparable to or better than those of conventional noble-metal catalysts such as Pt (η_ORR = 0.45 V), RuO₂ (η_OER = 0.42 V), and IrO₂ (η_OER = 0.56 V). Meanwhile, AgTAP shows ultralow overpotentials for both ORR and OER (0.36 V), highlighting its high intrinsic activity toward H₂O₂ production and making it a promising bifunctional catalyst. Although Cr-, Mn-, Fe-, Co-, and Zn-TAP exhibit limited selectivity for the 2e⁻ ORR, they demonstrate excellent 2e⁻ OER performance with low overpotentials of 0.49, 0.34, 0.39, 0.21, and 0.17 V, respectively. Moreover, CdTAP (−0.08 eV) and OsTAP (−0.12 eV) display very small ΔG_^*H values, identifying them as highly active HER candidates. The effective O₂ activation on the three selected catalysts is further elucidated through spin density, charge density difference, and electronic structure analysis. It is important to emphasize the limitations of this theoretical study. All calculations were performed using the CHE model with pH = 0, assuming an idealized acidic environment. Explicit solvent molecules and implicit solvation effects were not included, despite their known influence on the adsorption free energies of key intermediates (e.g., *OOH, *OH, *O). Consequently, the reported overpotentials should be interpreted as intrinsic thermodynamic trends under the CHE framework rather than absolute predictions for experimental conditions. Future work should incorporate implicit solvation models and explicit water molecules to bridge the gap between theory and experiment. Overall, TAP-supported transition-metal monolayers offer a viable route to designing low-cost, promising single-atom electrocatalysts for the ORR, OER, and HER, with promising applications in clean and renewable energy technologies.

Conflicts of interest

The authors declare no competing financial interest.

Data availability

The data supporting this article have been included as part of the supplementary information (SI). Supplementary information: crystal structure information of 30 TMTAP monolayers; calculated binding energies, cohesive energies, cluster energies, and dissolution potentials of TM atoms; AIMD simulation results of other remaining TMTAP monolayers; calculated ΔG_^*O values of the 20 stable SACs; the bond length of adsorbed O₂ molecule, magnetic moment of TM atom before and after O₂ adsorption for Fe-, Zn-, Rh-, and Ag-TAP catalysts; calculated Bader charges on TM and O atoms, magnetic moment of O atom after O₂ adsorption for Fe-, Zn-, Rh-, Ag-, Ti-, Mo-, and W-TAP catalysts; volcano plot of the η_ORR as a function of the d-band center, volcano plot of the η_OER as a function of the d-band center, volcano plot of the η_HER as a function of the Bader charge. See DOI: https://doi.org/10.1039/d6ra03135h.

Acknowledgements

This work was supported by the National Natural Science Foundation of China (12504130), the Natural Science Foundation of Hubei Province (2024AFB364, 2024AFB424), and the Natural Science Foundation of Xiaogan Municipality (XGKJ2025010095).

References

1. C. F. Shih, T. Zhang, J. Li and C. Bai, Powering the Future with Liquid Sunshine, Joule, 2018, 2(10), 1925–1949.
2. M. S. Dresselhaus and I. L. Thomas, Alternative energy technologies, Nature, 2001, 414(6861), 332–337.
3. N. S. Lewis, Research opportunities to advance solar energy utilization, Science, 2016, 351(6271), aad1920.
4. J. Kibsgaard and I. Chorkendorff, Considerations for the scaling-up of water splitting catalysts, Nat. Energy, 2019, 4(6), 430–433.
5. S. Dresp, F. Luo, R. Schmack, S. Kühl, M. Gliech and P. Strasser, An efficient bifunctional two-component catalyst for oxygen reduction and oxygen evolution in reversible fuel cells, electrolyzers and rechargeable air electrodes, Energy Environ. Sci., 2016, 9(6), 2020–2024.
6. D. Zhao, Z. Zhuang, X. Cao, C. Zhang, Q. Peng, C. Chen and Y. Li, Atomic site electrocatalysts for water splitting, oxygen reduction and selective oxidation, Chem. Soc. Rev., 2020, 49(7), 2215–2264.
7. H. Zhu, Y. Feng, D. Zheng, X. Zhao, Y. Zhou, X. Fu, L. Zhao and X. Chen, Theoretical prediction of emerging high-performance trifunctional ORR/OER/HER single-atom catalysts: transition metals anchored into π–π conjugated graphitic carbon nitride (g-C₁₀N₃), Phys. Chem. Chem. Phys., 2023, 25(46), 31983–31994.
8. C. Fang, J. Zhou, L. Zhang, W. Wan, Y. Ding and X. Sun, Synergy of dual-atom catalysts deviated from the scaling relationship for oxygen evolution reaction, Nat. Commun., 2023, 14(1), 4449.
9. V. R. Stamenkovic, B. Fowler, B. S. Mun, G. Wang, P. N. Ross, C. A. Lucas and N. M. Marković, Improved Oxygen Reduction Activity on Pt₃Ni(111) via Increased Surface Site Availability, Science, 2007, 315(5811), 493–497.
10. M. Li, Z. Zhao, T. Cheng, A. Fortunelli, C.-Y. Chen, R. Yu, Q. Zhang, L. Gu, B. V. Merinov, Z. Lin, E. Zhu, T. Yu, Q. Jia, J. Guo, L. Zhang, W. A. Goddard, Y. Huang and X. Duan, Ultrafine jagged platinum nanowires enable ultrahigh mass activity for the oxygen reduction reaction, Science, 2016, 354(6318), 1414–1419.
11. L. Huang, S. Zaman, X. Tian, Z. Wang, W. Fang and B. Y. Xia, Advanced Platinum-Based Oxygen Reduction Electrocatalysts for Fuel Cells, Acc. Chem. Res., 2021, 54(2), 311–322.
12. P. Mukherjee, B. Kakade and A. Swami, Current Trends in Platinum-Based Ternary Alloys as Promising Electrocatalysts for the Oxygen Reduction Reaction: A Mini Review, Energy Fuels, 2022, 36(5), 2306–2322.
13. J. Liu, M. Jiao, L. Lu, H. M. Barkholtz, Y. Li, Y. Wang, L. Jiang, Z. Wu, D.-j. Liu, L. Zhuang, C. Ma, J. Zeng, B. Zhang, D. Su, P. Song, W. Xing, W. Xu, Y. Wang, Z. Jiang and G. Sun, High performance platinum single atom electrocatalyst for oxygen reduction reaction, Nat. Commun., 2017, 8(1), 15938.
14. W. Ye, S. Chen, Y. Lin, L. Yang, S. Chen, X. Zheng, Z. Qi, C. Wang, R. Long, M. Chen, J. Zhu, P. Gao, L. Song, J. Jiang and Y. Xiong, Precisely Tuning the Number of Fe Atoms in Clusters on N-Doped Carbon toward Acidic Oxygen Reduction Reaction, Chem, 2019, 5(11), 2865–2878.
15. M. Tong, F. Sun, Y. Xie, Y. Wang, Y. Yang, C. Tian, L. Wang and H. Fu, Operando Cooperated Catalytic Mechanism of Atomically Dispersed Cu−N4 and Zn−N4 for Promoting Oxygen Reduction Reaction, Angew. Chem., Int. Ed., 2021, 60(25), 14005–14012.
16. Y. Wang, R. Hu, Y. Li, F. Wang, J. Shang and J. Shui, High-throughput screening of carbon-supported single metal atom catalysts for oxygen reduction reaction, Nano Res., 2022, 15(2), 1054–1060.
17. C. Yang, J. Yao, S. Meng, P. Wang, M. He, P. Li, P. Xiao, J. Xiao, Y. Liu and Z. Li, A Library of Polymetallic Alloy Nanotubes: From Binary to Septenary, J. Am. Chem. Soc., 2025, 147(11), 9865–9878.
18. Y. Zhao, S. Meng, P. Wang, C. Yang, R. Tan, Y. Hou, X. Zhang, Y. Hou and Z. Li, Porous rare earth-based high-entropy perovskite oxide nanosheet for stable acidic water splitting at ampere-level current density, J. Energy Chem., 2026, 117, 578–590.
19. M. He, P. Wang, J. Yao, Y. Li, S. Meng and Z. Li, Lanthanum-based nanomaterials for oxygen evolution reaction, J. Rare Earths, 2025, 43(6), 1091–1099.
20. H. Tian, A. Song, P. Zhang, K. Sun, J. Wang, B. Sun, Q. Fan, G. Shao, C. Chen, H. Liu, Y. Li and G. Wang, High Durability of Fe–N–C Single-Atom Catalysts with Carbon Vacancies toward the Oxygen Reduction Reaction in Alkaline Media, Adv. Mater., 2023, 35(14), 2210714.
21. X. Zhao, Z. H. Levell, S. Yu and Y. Liu, Atomistic Understanding of Two-dimensional Electrocatalysts from First Principles, Chem. Rev., 2022, 122(12), 10675–10709.
22. L. Yan, P. Li, Q. Zhu, A. Kumar, K. Sun, S. Tian and X. Sun, Atomically precise electrocatalysts for oxygen reduction reaction, Chem, 2023, 9(2), 280–342.
23. L. Wang, C. Zhu, M. Xu, C. Zhao, J. Gu, L. Cao, X. Zhang, Z. Sun, S. Wei, W. Zhou, W.-X. Li and J. Lu, Boosting Activity and Stability of Metal Single-Atom Catalysts via Regulation of Coordination Number and Local Composition, J. Am. Chem. Soc., 2021, 143(45), 18854–18858.
24. R. Lang, W. Xi, J.-C. Liu, Y.-T. Cui, T. Li, A. F. Lee, F. Chen, Y. Chen, L. Li, L. Li, J. Lin, S. Miao, X. Liu, A.-Q. Wang, X. Wang, J. Luo, B. Qiao, J. Li and T. Zhang, Non defect-stabilized thermally stable single-atom catalyst, Nat. Commun., 2019, 10(1), 234.
25. R. Gusmão, M. Veselý and Z. Sofer, Recent Developments on the Single Atom Supported at 2D Materials Beyond Graphene as Catalysts, ACS Catal., 2020, 10(16), 9634–9648.
26. H.-Y. Zhuo, X. Zhang, J.-X. Liang, Q. Yu, H. Xiao and J. Li, Theoretical Understandings of Graphene-based Metal Single-Atom Catalysts: Stability and Catalytic Performance, Chem. Rev., 2020, 120(21), 12315–12341.
27. Y. Wang, D. Wang and Y. Li, Rational Design of Single-Atom Site Electrocatalysts: From Theoretical Understandings to Practical Applications, Adv. Mater., 2021, 33(34), 2008151.
28. K. Liu, J. Fu, Y. Lin, T. Luo, G. Ni, H. Li, Z. Lin and M. Liu, Insights into the activity of single-atom Fe-N-C catalysts for oxygen reduction reaction, Nat. Commun., 2022, 13(1), 2075.
29. Y. Zhao, D. P. Adiyeri Saseendran, C. Huang, C. A. Triana, W. R. Marks, H. Chen, H. Zhao and G. R. Patzke, Oxygen Evolution/Reduction Reaction Catalysts: From In Situ Monitoring and Reaction Mechanisms to Rational Design, Chem. Rev., 2023, 123(9), 6257–6358.
30. C. Wan, X. Duan and Y. Huang, Molecular Design of Single-Atom Catalysts for Oxygen Reduction Reaction, Adv. Energy Mater., 2020, 10(14), 1903815.
31. L. Peng, L. Shang, T. Zhang and G. I. N. Waterhouse, Recent Advances in the Development of Single-Atom Catalysts for Oxygen Electrocatalysis and Zinc–Air Batteries, Adv. Energy Mater., 2020, 10(48), 2003018.
32. S. H. Joo, Metallomacrocycle-Derived Atomically Dispersed Metal Electrocatalysts with Controlled Selectivity for Small Molecule Conversion Reactions, ECS Meet. Abstr., 2022,(14), 960.
33. S. Cao, W. Wu, C. Liu, L. Song, Q. Xu, H. Zhang and Y. Zhao, Supramolecular Macrocycle Regulated Single-Atom MoS2@Co Catalysts for Enhanced Oxygen Evolution Reaction, Energy Environ. Mater., 2024, 7(5), e12702.
34. G. E. Ficken, R. P. Linstead, E. Stephen and M. Whalley, Conjugated macrocycles. Part XXXI. Catalytic hydrogenation of tetrazaporphins, with a note on its stereochemical course, J. Chem. Soc., 1958, 3879–3886.
35. T. Fukuda and N. Kobayashi, Hydrogenated tetraazaporphyrins—old but new core-modified phthalocyanine analogues, Dalton Trans., 2008,(35), 4685–4704.
36. G. Kresse and J. Furthmüller, Efficient Iterative Schemes for Ab Initio Total-Energy Calculations Using a Plane-Wave Basis Set, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54(16), 11169.
37. P. E. Blöchl, Projector augmented-wave method, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50(24), 17953.
38. J. P. Perdew, K. Burke and M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett., 1996, 77(18), 3865.
39. H. J. Monkhorst and J. D. Pack, Special Points for Brillouin-Zone Integrations, Phys. Rev. B: Condens. Matter Mater. Phys., 1976, 13(12), 5188.
40. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu, J. Chem. Phys., 2010, 132(15), 154104.
41. S. Nosé, A unified formulation of the constant temperature molecular dynamics methods, J. Chem. Phys., 1984, 81(1), 511.
42. G. Henkelman, A. Arnaldsson and H. Jónsson, A fast and robust algorithm for Bader decomposition of charge density, Comput. Mater. Sci., 2006, 36(3), 354–360.
43. E. F. Holby, G. Wang and P. Zelenay, Acid stability and demetalation of PGM-free ORR electrocatalyst structures from density functional theory: a model for “single-atom catalyst” dissolution, ACS Catal., 2020, 10(24), 14527–14539.
44. X. Guo, S. Lin, J. Gu, S. Zhang, Z. Chen and S. Huang, Simultaneously achieving high activity and selectivity toward two-electron O₂ electroreduction: the power of single-atom catalysts, ACS Catal., 2019, 9(12), 11042–11054.
45. R. D. Johnson III NIST computational chemistry comparison and benchmark database, NIST standard reference database number 101, 2022. https://cccbdb.nist.gov/.
46. J. K. Nørskov, T. Bligaard, A. Logadottir, J. R. Kitchin, J. G. Chen, S. Pandelov and U. Stimming, Trends in the exchange current for hydrogen evolution, J. Electrochem. Soc., 2005, 152(3), J23.
47. S. Peng, X. Wu, W. Cheng, X. Li, Z. Cheng and D.-B. Long, Single metal atom anchored graphitic carbon nitride monolayers as promising multifunctional electrocatalysts for the oxygen reduction reaction, oxygen evolution reaction, and hydrogen evolution reaction, Appl. Surf. Sci., 2025, 695, 162840.
48. J. K. Nørskov, J. Rossmeisl, A. Logadottir, L. Lindqvist, J. R. Kitchin, T. Bligaard and H. Jónsson, Origin of the overpotential for oxygen reduction at a fuel-cell cathode, J. Phys. Chem. B, 2004, 108(46), 17886–17892.
49. Y. Yao, X.-T. Chen, X. Zhang, S. Jin, Z. Tian, G. Li and L.-M. Yang, Monoatomic metalloporphyrinoid catalysts for efficient oxygen reduction, Rare Met., 2025, 44(6), 3920–3933.
50. S. Wang, K. Meng, L. Qin, Y. Wu, Y. Wei, J. Rong, Y. Sui, X. Yu and Y. Yang, High-throughput screening of transition metal phthalocyanine electrocatalysts for oxygen reduction reactions, Int. J. Hydrogen Energy, 2024, 88, 850–857.
51. V. Hoyos-Sinchi, P. H. Souza and W. Orellana, Computational Insights into Two-Dimensional M₃(HHTP)₂ Metal–Organic Frameworks as ORR/OER Electrocatalysts, J. Phys. Chem. C, 2025, 129(31), 14002–14010.
52. Y. Jin, T. Liu, X. Liu, W. Deng, W. Xiao, Y. Su, X. Qi, G. Wang, D. Wang, M. Liu, Y. Wu, A. Aierken, X. Chen, X. Wang, C. Gao, H. Xiong, X. Wu, J. Ge and J. Bi, Efficient Bifunctional Electrocatalysts for Oxygen Evolution/Reduction Reactions in Two-Dimensional Metal–Organic Frameworks by a Constant Potential Method, Langmuir, 2025, 41(1), 745–754.
53. M. Li, L. Zhang, Q. Xu, J. Niu and Z. Xia, N-doped graphene as catalysts for oxygen reduction and oxygen evolution reactions: Theoretical considerations, J. Catal., 2014, 314, 66–72.
54. H. Deng, D. Deng, S. Jin, Z. Tian and L.-M. Yang, Unraveling the Activity and Mechanism of TM@g-C₄N₃ Electrocatalysts in the Oxygen Reduction Reaction, ACS Appl. Mater. Interfaces, 2024, 16(14), 17617–17625.
55. I.-C. Man, I. Trancă and S.-G. Soriga, First principle studies of oxygen reduction reaction on N doped graphene: Impact of N concentration, position and co-adsorbate effect, Appl. Surf. Sci., 2020, 510, 145470.
56. I. C. Man, H.-Y. Su, F. Calle-Vallejo, H. A. Hansen, J. I. Martínez, N. G. Inoglu, J. Kitchin, T. F. Jaramillo, J. K. Nørskov and J. Rossmeisl, Universality in Oxygen Evolution Electrocatalysis on Oxide Surfaces, ChemCatChem, 2011, 3(7), 1159–1165.
57. L.-C. Ma, H.-J. Wang, H. Li, P.-G. Yuan and J.-M. Zhang, Theoretical screening of highly efficient multifunctional single-atom catalysts supported by pc-C₃N₂ monolayers for the electrocatalytic HER, OER and ORR, New J. Chem., 2025, 49(5), 1672–1685.
58. S.-L. Li, G. Tian, Y. Chen, Y. Zhao, F. Cai and L. Qiao, Abundant active sites in M₂C₁₂ monolayer as single-atom catalysts with efficient catalytic activity for hydrogen evolution reaction, Int. J. Hydrogen Energy, 2025, 144, 1043–1050.

---