Computational investigation of CO₂ hydrogenation and electrochemical reduction to formic acid and methanol using chalcogen-doped nitrogen-graphene nanoflake as a metal-free catalyst

Abstract

Clean fuel production and pollutant removal are critical industrial challenges. This study presents a model for CO₂ adsorption, focusing on non-metallic, biocompatible catalysts. We demonstrate the potential of chalcogen-doped graphene nanoflake modified with pyridinic nitrogen atoms as biocompatible catalysts for hydrogenating CO₂ to formic acid and methanol. Using dispersion-corrected density functional theory, the catalytic hydrogenation and electrochemical reduction of CO₂ over single chalcogen atoms (Se, Te) embedded in nitrogen-doped graphene nanoflake were analyzed. High hybridization between Se/Te and N states near the Fermi level stabilizes the chalcogen atoms on the graphene nanoflake. Se-doped nitrogen-containing graphene nanoflake exhibits lower energy barriers and higher catalytic performance than Te, facilitating CO₂ conversion to formic acid and methanol with improved stability. Migration barriers, electronic structures, and adsorption energies highlight Se-doped nitrogen-containing graphene nanoflake as an efficient catalyst for CO₂ hydrogenation and electrochemical reduction at room temperature. In particular, Se doping not only stabilizes the material but also improves catalytic performance for the selective production of CO/CH₃OH. This work introduces Se-doped nitrogen-containing graphene nanoflake as a promising non-metal biocompatible catalyst.

Introduction

Anthropogenic greenhouse gas emissions are driving climate change, with CO₂ from industrial operations accounting for a significant portion of these emissions. Gas flaring—widely practiced in oil and gas extraction and chemical manufacturing to burn off surplus gases—releases large volumes of CO₂, exacerbating global warming and posing serious public health risks due to byproduct pollutants and radiation exposure near flare sites^1,2,3. Mitigating flare gas emissions through capture and utilization strategies is therefore essential not only for meeting climate targets but also for promoting environmental and economic sustainability⁴. Current carbon-capture efforts focus largely on large point sources such as power plants and manufacturing facilities, employing functional materials, including amine‐based solvents, porous adsorbents, and selective membranes, to sequester CO₂⁵. However, the high energy penalty associated with regenerating CO₂ from its oxidized state remains a critical hurdle⁶. Integrating capture with on‐site conversion, rather than sequential capture then reaction, can optimize energy use and process intensification, especially for low‐concentration, intermittent streams like flare gases^7,8.

Among advanced sorbents and catalysts, metal–organic frameworks (MOFs) and nanostructured metal catalysts have attracted attention for CO₂ capture and electrochemical reduction. For example, MgO-nanoparticle‐embedded MOF/polymer hybrids achieve high CO₂ uptake from flue gas⁹, while Pd‐decorated WO₂ nanowires enable selective photocatalytic CO₂-to-CO conversion¹⁰. Density functional theory (DFT) and climbing image nudged elastic band (CI-NEB) studies of Fe, Co, and Ni clusters on carbon nanotubes have highlighted Ni’s superior catalytic activity toward CO₂ reduction¹¹.

Organic chalcogen compounds, particularly those containing selenium and tellurium, represent a versatile class of materials with growing importance in chemistry and materials science. Their incorporation into organic and inorganic frameworks imparts unique structural, electronic, and catalytic properties that underpin a wide range of applications in technology, medicine, agriculture, and industry. Selenium-based reagents have proven to be highly efficient and recyclable catalysts in various organic transformations, offering significant advantages for industrial applications. Meanwhile, tellurium derivatives exhibit remarkable optical, thermal, and electrical characteristics, supporting their integration into advanced materials such as quantum dots, semiconductors, electronic devices, and biomedical systems. The distinctive dual nucleophilic and electrophilic reactivity of selenium and tellurium further emphasizes their value in modern synthetic methodologies and highlights their expanding role in innovative materials development^12,13,14,15. More recently, benzyl-ditelluride exhibited a near-thermoneutral reversible CO₂ capture–release cycle, outperforming other chalcogens in our earlier work¹⁶.

Formate and methanol, two key C₁-products, have also been explored on metal surfaces (e.g., Sn, Cu) and doped carbon nanotubes, revealing promising electrocatalysts for selective HCOO⁻ and CH₂OH formation with low overpotentials^17,18. Among the possible CO₂ reduction products (CO, CH₄, C₂H₄, CH₃OH, HCOOH, etc.), we focus on formic acid (HCOOH) and methanol (CH₃OH) because they require fewer proton–electron transfer steps. This makes them thermodynamically more accessible and easier to achieve with high selectivity compared to multi-carbon products that demand complex C–C coupling. Both HCOOH and CH₂OH are liquid at ambient conditions, simplifying storage, handling, and transport. Importantly, they are also high-value products: HCOOH is used as a chemical feedstock and as a precursor for liquid organic hydrogen carriers (LOHCs), while CH₂OH serves as a key industrial solvent, a precursor to formaldehyde and acetic acid, and a sustainable liquid fuel. These advantages explain why HCOOH and CH₂OH are prioritized in catalytic CO₂ hydrogenation and electrochemical reduction studies^19,20.

Graphene-based materials offer a versatile platform for electrocatalytic CO₂ reduction (ECR) due to their high conductivity, tunable surface chemistry, and robust mechanical properties²¹. Heteroatom doping (N, S, Se, Te) and generation of defects create active sites that enhance CO₂ adsorption and activation, while single-atom catalysts on doped graphene further improve selectivity and stability for hydrogenation to formic acid^22,23,24. Se- and Te-doped graphene and carbon nanotubes have particularly shown favorable Gibbs free energies for methanol production in the gas phase²⁵. Yet, despite these advances, a detailed atomistic understanding of non-metal active sites, especially chalcogen-doped, nitrogen-modified graphene nanoflake, remains underexplored.

The design of transition metal-free electrocatalysts for energy conversion processes, such as the oxygen reduction reaction (ORR) and oxygen evolution reaction (OER), holds great promise for advancing renewable energy technologies, particularly fuel cells and metal–air batteries²⁶. Electrocatalytic activity of MN₄-doped (M=Cu, Ag, and Zn) single-walled carbon Nnanotubes in ORR was studied by Kuzmin and Shainyan. Due to the obtained results, the relative catalytic activity toward ORR was increased in the order Zn < Ag ≪ Cu. The combined catalytic activity of the two sites (AgN₄/C₂ sites) was predicted for the AgN₄-CNT catalyst²⁷.

In the recent study on the role of single and double atom dopants on the catalytic activity of carbon-based materials, density functional theory calculations demonstrate that γ-graphyne doped with various heteroatoms exhibits distinct nitrogen reduction reaction (NRR) activities, with B(sp²)-C(sp)-B(sp) configuration showing superior performance. Boron-based doping of graphyne and g-C₂N monolayers is effective for promoting NRR activity²⁸. In 2024, a DFT study of defect engineering in hexagonal boron nitride (h-BN) nanosheets declared that vacancy sites such as VN and VB₂N selectively promote CO₂ electroreduction to methanol and methane, while suppressing hydrogen evolution²⁹. In 2025, DFT calculations demonstrate that transition-metal-doped BN (TM@BN, TM = Sc, V, Mn, Fe, Ni) presents distinct CO₂ reduction activities, with V@BN emerging as a promising catalyst. V@BN demonstrates favorable selectivity toward CH₄, HCOOH, and CH₂OH formation with low limiting potentials³⁰. Point defects in bilayer borophene (BL-B) were investigated for CO₂ activation and electrochemical reduction to CH₂OH, showing how defect-induced charge modulation enhances catalytic performance³¹.

In this work, we employ dispersion interactions and van der Waals effects (DFT + D3) calculations to investigate graphene nanoflake sheet co-modified with pyridinic nitrogen and doped with selenium (Se) or tellurium (Te). Building on evidence of Se and Te’s catalytic efficacy for CO₂ electrochemical reduction to methanol and formic acid^25,32,33, the stability of graphene nanoflake sheet doped with Se and Te and modified with pyridinic nitrogen was first evaluated. Finally, the hydrogenation of CO₂ to formic acid (HCOOH) and its electrochemical reduction to methanol (CH₂OH) on these biocompatible catalysts are systematically analyzed. HCOOH is used directly as a chemical feedstock and as a liquid organic hydrogen carrier (LOHC) for hydrogen storage and transport, while CH₂OH is a key industrial chemical (solvent, precursor to formaldehyde, acetic acid), and an excellent liquid fuel for sustainable energy systems^19,20.

The findings highlight the potential of chalcogen-doped graphene nanoflake sheet modified with pyridinic nitrogen atoms as biocompatible and efficient catalysts. This theoretical study provides valuable insights into CO₂ adsorption, reduction, and hydrogenation mechanisms, paving the way for advancing CO₂ conversion technologies and improving catalytic efficiency.

Results and discussion

Geometry, electronic Structure, and thermodynamic stability of SeN₄Gr and TeN₄Gr

Prior theoretical and experimental studies have shown that introducing pyridinic N at graphene vacancy defects perturbs the electronic structure near the frontier region, increases local charge density, and enhances catalytic activity^34,35,36,37. In our nanoflake model, substituting the four carbons around a di-vacancy with N atoms eliminates dangling bonds and stabilizes the lattice, although the formation energy of N₄Gr is 3.97 eV, confirming that the pristine N-decorated vacancy is not thermodynamically favorable. NBO analysis assigns a net charge of − 0.22 e to each pyridinic N, consistent with electronic redistribution at the defect.

SeN₄Gr and TeN₄Gr were constructed by embedding a single Se or Te atom into the pyridinic-nitrogen–modified di-vacancy (Fig. 1). In both cases, the chalcogen binds to the four surrounding N atoms, giving a quasi-planar N₄ coordination. Upon relaxation, the Se–N bond length is 2.09 Å, whereas Te–N averages 2.19 Å; the out-of-plane displacement of Te (~ 1.16 Å) exceeds that of Se, consistent with Te’s larger covalent radius. Natural Bond Orbital (NBO) analysis confirms net positive charges on the chalcogens (+ 0.47 e on Se; + 0.71 e on Te), indicative of partial electron donation to the N₄ framework. Thermodynamically, Se incorporation is more favorable: the formation energies are − 2.26 eV (SeN₄Gr) and − 1.97 eV (TeN₄Gr). The adsorption energy (E_ads) of the chalcogen also follows this trend, with E_ads(Se) = − 3.52 eV and E_ads(Te) = − 1.47 eV (Table 1). The latter aligns with reported Te doped graphene values (E_ads ≈ − 1.04 eV for Te-Gr)²⁵.

Table 1 Calculated adsorption energy (E_ads, eV), formation energy (E_formation, eV), dipole moment (µ_B) (Debaye), and the atomic charge on the chalcogen atom (Se, Te) atom (Q_M, e).

The electronic density of states, DOS, (Fig. 1a–c) provides further insight. N₄Gr shows a finite gap when referenced to mid-gap (0 eV), with only modest N contributions spread across the valence and conduction bands, confirming that pyridinic N introduces localized states without generating dominant frontier features. Upon anchoring Se or Te, additional states emerge near the mid-gap region, increasing the TDOS intensity at the band edges. SeN₄Gr shows stronger redistribution in this region than TeN₄Gr, consistent with its shorter bond length, more favorable embedding energy, and lower dipole moment (2.22 D vs. 2.50 D). In contrast, TeN₄Gr displays weaker intensity near the frontier and more spectral weight deeper in the valence band, reflecting weaker stabilization of the active site.

Electrostatic potential (ESP) maps (Fig. 2) further confirm this trend, the regions around the chalcogen atoms (Se, Te) and their neighboring nitrogen atoms are slightly positively charged (green), while the areas surrounding the nitrogen atoms in N₄Gr are negatively charged (red). This positive charge indicates a strong affinity for interactions with negatively charged species, suggesting that SeN₄Gr and TeN₄Gr exhibit high adsorption potential for gas molecule. This polarization pattern suggests a high chemisorption affinity for electron-rich species at the chalcogen site. Localization function (ELF) and localized orbital locator (LOL) plots reinforce this picture: SeN₄Gr exhibits continuous bonding basins along the N–Se bonds and a clear ELF maximum at Se, whereas TeN₄Gr shows comparatively weaker localization at Te. Overall, the covalent contribution follows Se–N > Te–N, consistent with shorter Se–N bonds and the more favorable Se embedding energy. The more continuous LOL features in SeN₄Gr also indicate greater electronic delocalization and stability relative to TeN₄Gr. In N₄Gr, by contrast, lower ELF/LOL values around the N atoms highlight the primarily electrostatic character of the bare N-decorated vacancy.

Taken together, the structural, thermodynamic, DOS, ESP/ELF/LOL, and dipole moment data consistently identify SeN₄Gr as the more stable and electronically responsive platform. This superior balance of stability and reactivity suggests that SeN₄Gr will provide stronger adsorption and activation of reactants in the subsequent hydrogenation and electrochemical reduction pathways.

Adsorption of H₂, CO₂, and H₂–CO₂ on SeN₄Gr and TeN₄Gr

Efficient adsorption and activation of reactants at the catalyst surface are prerequisites for electrochemical CO₂ hydrogenation^38,39,40. To evaluate this, we examined the adsorption of H₂, CO₂, and their coadsorption on SeN₄Gr and TeN₄Gr. Optimized geometries are shown in Fig. 3, and E_ads and charge-transfer values are summarized in Table 2.

Table 2 Adsorption energy (E_ads (eV)) and electronic charge changes (Q |e|) on H₂, CO₂, and H₂−CO₂ adsorbed on the surface of SeN₄Gr and TeN₄Gr.

On SeN₄Gr, H₂ binds with E_ads = − 0.33 eV, compared to − 0.29 eV on TeN₄Gr. In both cases, the H–H bond elongates from 0.73 Å in free H₂ to 0.79 Å (SeN₄Gr) and 0.81 Å (TeN₄Gr), indicating significant activation toward dissociation. NBO analysis shows charge depletion of + 0.67 e around Se and − 1.21 e around Te upon H₂ adsorption, consistent with back-donation from the chalcogen–nitrogen framework into the H₂ σ* orbital.

CO₂ adsorption is comparatively weaker than H₂, with E_ads = − 0.28 eV on SeN₄Gr and − 0.15 eV on TeN₄Gr. Adsorption induces bending of the O–C–O angle from 180° to 136° (SeN₄Gr) and 148° (TeN₄Gr), together with elongation of one C–O bond by ~ 0.030 Å. NBO analysis indicates partial charge transfer of − 0.12 e (SeN₄Gr) and − 0.07 e (TeN₄Gr), confirming modest activation. These results reinforce that both catalysts preferentially adsorb and activate H₂ over CO₂, a favorable trait for initiating hydrogenation.

When H₂ and CO₂ are coadsorbed, adsorption energies become more exergonic (Table 2). In these complexes, H₂ remains chemisorbed to the chalcogen center, while CO₂ physisorbs via its O atoms onto adjacent pyridinic N sites (Fig. 3). Pre-adsorbed H₂ induces further CO₂ bending (O–C–O ≈ 130° on SeN₄Gr), reflecting synergistic charge redistribution that may lower the barrier for the first hydrogenation step.

Overall, SeN₄Gr consistently exhibits stronger adsorption and greater activation of both H₂ and CO₂ compared to TeN₄Gr, as evidenced by more negative adsorption energies, larger bond distortions, and higher charge transfer. These results establish SeN₄Gr as a more effective non-metal catalyst for CO₂ hydrogenation, combining moderate binding strength with pronounced activation, key features for efficient conversion to formic acid or methanol.

Reaction mechanism for CO₂ hydrogenation to formic acid on SeN4Gr and TeN4Gr

The catalytic hydrogenation of CO₂ on chalcogen-doped, nitrogen-modified graphene nanoflake (SeN₄Gr and TeN₄Gr) was investigated by considering two plausible pathways: (i) coadsorption of H₂ and CO₂, leading to a carboxyl intermediate, and (ii) initial H₂ adsorption followed by CO₂ insertion to form a formate intermediate. Previous theoretical studies indicate that the coadsorption pathway typically involves prohibitively high activation energies for H₂ and CO₂ dissociation on doped carbon systems, rendering it less feasible^41,42. Our results are consistent with this, showing that the formate pathway is energetically preferred, particularly on SeN₄Gr, which exhibits lower reaction barriers than TeN₄Gr.

The minimum energy pathway (MEP) and optimized structures are presented in Fig. 4. The sequence begins with H₂ adsorption and dissociation on the catalyst surface. At the first transition state (TS1), the H–H bond elongates from 0.74 Å in free H₂ to 1.45 Å on SeN₄Gr and 1.64 Å on TeN₄Gr. During this step, one H atom remains coordinated to the chalcogen center, forming H–Se and H–Te bonds of 1.60 Å and 1.93 Å, respectively, while the second H atom migrates toward the graphene nanoflake support.

Formation of the first intermediate (MS₁) requires activation energies of 0.09 eV (SeN₄Gr) and 1.11 eV (TeN₄Gr) (Table 3). At this stage, the H–H bond is completely cleaved, yielding a new H–N bond with lengths of 1.09 Å on SeN₄Gr and 1.14 Å on TeN₄Gr, confirming successful H₂ dissociation. CO₂ then interacts with MS₁ and proceeds through the second transition state (TS₂). The associated barriers are 0.50 eV for SeN₄Gr and 1.46 eV for TeN₄Gr, demonstrating the catalytic advantage of Se. At TS₂, the carbon atom of CO₂ engages one hydrogen atom bound to the chalcogen, while one oxygen atom coordinates with the chalcogen site. This reduces the Se–O and Te–O bond lengths to 1.90 Å and 2.06 Å, respectively, while forming a new C–H bond of 1.12 Å, signaling the onset of formate (HCOO^-) intermediate generation.

Table 3 Calculated Gibbs energies (ΔG) (eV) and activation energies (ΔG^#) (eV) for the studied steps of the CO₂ hydrogenation to formic acid on the proposed catalytic cycle of SeN₄Gr and TeN₄Gr.

The reaction then advances through the third transition state (TS₃), culminating in HCOOH formation. The activation energies for this step are 0.69 eV on SeN₄Gr and 0.48 eV on TeN₄Gr. However, the overall maximum barrier for the pathway is markedly lower for SeN₄Gr (0.69 eV) compared to TeN₄Gr (1.46 eV), underscoring the superior catalytic efficiency of Se in driving CO₂ hydrogenation. Following product formation, desorption occurs via MS₃, with Se–O and Te–O bond lengths extending to 3.08 Å and 3.98 Å, respectively. Thermodynamically, desorption is more favorable for SeN₄Gr (ΔG = − 0.79 eV) than for TeN₄Gr (ΔG = − 0.22 eV), confirming easier product release and higher turnover potential on Se-based catalysts.

Overall, both systems demonstrate catalytic competence, but SeN₄Gr consistently outperforms TeN₄Gr (Table 3). SeN₄Gr features a much lower initial activation barrier (0.09 eV) and a more favorable desorption step, reinforcing its superior activity and thermodynamic favorability. The rate-determining step (RDS) for formic acid formation is most likely the conversion of the HCOO intermediate into HCOOH. Separation of HCOOH from the SeN₄Gr surface is thermodynamically accessible, with a Gibbs free energy of − 0.36 eV, confirming efficient product release. This performance aligns well with previous theoretical studies^42,43,44.

When benchmarked against alternative catalysts—including Ni-, Pt-, or Cu-doped graphene and noble-metal catalysts⁴⁵—the activation energies of SeN₄Gr are substantially lower, highlighting its promise as a cost-effective and sustainable alternative. The intrinsic role of defect sites is critical: doping with chalcogens at pyridinic N vacancies modulates adsorption energies, stabilizes intermediates, and enhances hydrogenation efficiency. Calculated Gibbs free energies (ΔG) for hydrogen adsorption on SeN₄Gr and TeN₄Gr are smaller than those for bare carbon and graphitic N ⁴⁶, confirming that these doped systems are thermodynamically better suited for hydrogen adsorption/desorption.

In summary, SeN₄Gr demonstrates pronounced advantages in CO₂ hydrogenation due to lower energy barriers, stronger stabilization of intermediates, and more favorable desorption. These attributes underscore the unique role of Se in enhancing the catalytic function of doped graphene and establish chalcogen-doped nitrogen–containing graphene nanoflake as a next-generation non-metal electrocatalyst for efficient CO₂conversion.

Electrochemical CO₂ reduction to methanol on SeN4Gr and TeN4Gr

The primary objective of CO₂ conversion technologies is to generate hydrocarbon molecules that can serve as “green fuels,” thereby reducing greenhouse gas emissions. Considerable advances have been achieved with metal-based catalysts for the photoelectrochemical conversion of CO₂ into hydrocarbons. Nevertheless, their high cost, scarcity, and limited scalability have motivated the search for alternatives.

In CO₂ electrochemical reduction, methane formation is typically less favorable than methanol production. This is because CH₄ generation requires more proton–electron transfers, complete C–O bond cleavage, higher overpotentials, and stabilization of difficult intermediates. In contrast, methanol, as a partial reduction product, involves fewer kinetic and thermodynamic barriers, making its formation comparatively more efficient^47,48,49.

In this context, chalcogen-doped graphene structures have emerged as particularly promising due to their ease of preparation, relatively low cost, and potential for large-scale application⁵⁰. Although challenges remain—including ongoing debates over product selectivity (CO vs. CH₂OH or C₂H₅OH) and the generally low yields reported to date—recent studies highlight the capacity of chalcogen-doped systems to act as effective CO₂ conversion catalysts. These developments establish chalcogen-doped materials as a viable and cost-efficient alternative to traditional metal-based catalysts in the pursuit of sustainable energy solutions.

Mid-gap–referenced DOS/PDOS overlays (Fig. 5a,b), shown here as examples, reveal that CO₂ adsorption introduces additional CO₂-derived states near the HOMO–LUMO region (grey traces). However, the substrate PDOS (C and N) remains essentially unchanged, with no clear peak shifts or broadening within ± 1 eV of mid-gap. This lack of significant redistribution indicates limited hybridization between CO₂ frontier orbitals and the chalcogen–N₄ site, consistent with predominantly physisorptive binding. SeN₄Gr shows slightly stronger intensity than TeN₄Gr, but without notable spectral reshaping.

NBO analysis further reveals strong charge redistribution during this step: the carbon atom of CO₂ acquires − 0.35 e, while Se carries + 0.46 e. On TeN₄Gr, CO₂ activation is also observed, with C carrying − 0.37 e and Te + 0.72 e. This charge separation creates a strong electrostatic attraction between the negatively charged carbon atom and the positively charged chalcogen center, promoting bond formation, CO₂ activation, and subsequent reduction. The effects of electronic charge redistribution and intermolecular orbital interactions on the reaction pathway were examined using the NBO approach. Donor–acceptor interaction energies are calculated at the CO₂ adsorption in the presence of SeN₄Gr and TeN₄Gr. These calculations indicate the effects of SeN₄Gr and TeN₄Gr on the stabilization energies of CO₂ adsorption, using second-order perturbation energies, E⁽²⁾. The stabilization energies are related to the lp_Se → σ*_C−O (2.51 kcal mol^− 1) and lp_Te → σ*_C−O (1.68 kcal mol^− 1) in the presence of SeN₄Gr and TeN₄Gr, respectively, which are based on bond formation. The simultaneous presence of a negatively charged carbon atom in CO₂ and a positively charged chalcogen atom promotes a electrostatic attraction, ultimately facilitating chemical bond formation.

To assess this potential, we investigated CO₂RR on SeN₄Gr and TeN₄Gr, focusing on how the chalcogen dopant modulates CO₂ adsorption and activation. The reaction pathway follows the formaldehyde route after CO formation, with SeN₄Gr showing pronounced selectivity for CH₂OH, as depicted in Fig. 6. According to Table 4, CO₂ physisorption on SeN₄Gr begins spontaneously with a Gibbs free energy change of − 0.18 eV at 298 K.

Table 4 Gibbs reaction energies (eV) for the catalytic performance for the electrochemical reduction mechanism of CO₂ to CH₃OH.

The first hydrogenation step, involving transfer of a single H⁺/e⁻ pair to form the HOCO° intermediate, requires a Gibbs free energy of 0.14 eV. Since this initial step is widely regarded as the rate-determining step in CO₂ reduction, the low barrier highlights SeN₄Gr’s efficiency as a catalyst. Once HOCO° forms, the reaction profile becomes strongly exothermic, with subsequent transformations yielding neutral products such as CO (with H₂O release), H₂CO, and finally CH₂OH. Along this route, radical intermediates HCO° and CH₂O° are generated. As summarized in Table 4, their formation is highly favorable, releasing − 2.47 eV and − 2.51 eV of energy, respectively. These radicals interact strongly with the Se site in nitrogen-doped graphene nanoflake, underscoring the critical role of selenium in stabilizing key intermediates. Interestingly, while HOCO° and HCO° coordinate to the Se atom through their carbon centers, DFT-D3 calculations reveal no evidence of O–N bonding during the H₂CO or CH₂O° steps. Instead, as shown in Fig. 6, the Se atom facilitates the chemisorption of the CH₂O° radical, confirming its importance in the catalytic network that drives methanol production. Overall, this pathway demonstrates the potential of SeN₄Gr as a highly effective and selective non-metal catalyst for CO₂ electrochemical reduction to methanol.

Mechanistic DFT-D3 (PBEPBE-D3/6-31G+(d, p)) analysis of TeN₄Gr reveals distinct behavior compared with SeN₄Gr (Fig. 6). The first hydrogenation step, involving transfer of one H⁺/e⁻ pair to form the HOCO° radical, requires a Gibbs free energy of 0.24 eV—noticeably higher than the 0.14 eV required on SeN₄Gr. This indicates that TeN₄Gr provides weaker stabilization of the intermediate and makes the initial hydrogenation step the rate-determining step (RDS). In later stages, differences between the two catalysts become more pronounced. For example, during the third H⁺/e⁻ pair addition, the energy release on TeN₄Gr decreases to − 2.04 eV, compared with − 2.47 eV on SeN₄Gr after HCO° formation. This reduced exothermicity reflects the weaker interaction of Te with key intermediates. Furthermore, the lower O-philicity of the TeN₄Gr surface diminishes its ability to stabilize oxygen-containing species, further hindering the reaction pathway. Collectively, these results show that while TeN₄Gr exhibits catalytic activity and can promote CO₂ reduction, it is consistently outperformed by SeN₄Gr, which offers superior stabilization of intermediates and more favorable energetics along the methanol formation pathway.

In the final step of electrochemical CO₂ reduction, methanol (CH₂OH) formation was analyzed. On SeN₄Gr, this step is markedly more thermodynamically favorable than on TeN₄Gr, as indicated by the more negative Gibbs free energy (ΔG). The formation of CH₂OH proceeds through protonation of the O-moiety of the CH₂O° intermediate, which remains coordinated to the chalcogen center. A direct comparison of computed Gibbs free energies (Table 4) underscores the superior catalytic efficiency of SeN₄Gr: its lower ΔG reflects stronger stabilization of intermediates and more favorable energetics for CH₂OH production. This suggests that SeN₄Gr not only accelerates the electrochemical reduction of CO₂ but also enhances the yield of methanol relative to TeN₄Gr. Overall, the more negative Gibbs free energy profile and faster conversion rates observed on Se-doped, nitrogen-modified graphene nanoflake highlight its superior catalytic performance, reinforcing SeN₄Gr as a more effective and efficient platform than its Te-based counterpart.

Recent efforts to design efficient and stable electrocatalysts for the CO₂ reduction reaction (CO₂RR) have leveraged dispersion-corrected periodic DFT to evaluate transition-metal-free homo- and hetero-biatomic systems composed of Al, Be, B, and Si anchored on a 7,7,8,8-tetracyanoquinodimethane (TCNQ) monolayer. Free-energy analyses identified Al₂-TCNQ, AlBe-TCNQ, and BeSi-TCNQ as promising candidates for selective CO₂ reduction to methanol, with BeSi-TCNQ in particular achieving methanol production at a remarkably low free energy of − 0.29 V, demonstrating exceptional catalytic efficiency⁵¹. In comparison, SeN₄Gr and TeN₄Gr exhibit an even more negative Gibbs free energy profile (Table 5), further reinforcing their promise as cost-effective, non-metallic catalysts. These findings underscore the potential of chalcogen-doped, nitrogen-modified graphene nanoflake as a next-generation platform for selective electrochemical CO₂ reduction, while also providing design principles for the development of sustainable and scalable catalytic systems.

Table 5 Comparison of transition-metal-free homo-/hetero-biatomic, single atom (SA)-functionalized graphitic carbon nitride (g-C₂N), and chalcogen-doped nitrogen-containing graphene nanoflake for CO₂ hydrogenation and electrochemical reduction.

Recent investigations on single-atom (SA) functionalized g-C₂N monolayers have demonstrated their strong potential in CO₂ activation and reduction, with Al-SACs exhibiting remarkable activity toward HCOOH formation and B-SACs enabling methanol production through favorable desorption energetics⁵². Building on these advances, our results show that non-metal chalcogen-doped nitrogen-containing graphene nanoflake provides an even more advantageous thermodynamic landscape. Specifically, SeN₄Gr and TeN₄Gr display a more negative Gibbs free energy profile for CO₂ conversion to both HCOOH and CH₂OH (Table 5), underscoring their potential as efficient, sustainable alternatives to conventional metal-based electrocatalysts. The Table 5, compares catalyst type, host structure, selectivity toward HCOOH/CH₂OH, and their associated ΔG values and overpotentials, which serve as key performance descriptors for CO₂ hydrogenation and electrochemical reduction. Our results show that SeN₄-Gr exhibits the most favorable free energy (− 0.79 eV) and overpotential (− 2.06 eV) among the non-metal chalcogen-doped graphene systems, outperforming not only TeN₄-Gr (− 0.22 eV/−1.71 eV) but also several state-of-the-art SAC/g-C₂N and TCNQ-based biatomic catalysts. This strongly supports our conclusion that Se-doped, nitrogen-modified graphene is a highly efficient and cost-effective non-metal biocompatible catalyst for CO₂ hydrogenation and electrochemical reduction.

Conclusion

DFT-D3 (PBEPBE-D3/6–31 + G(d, p)) calculations were performed to investigate the catalytic pathways of SeN₄Gr and TeN₄Gr for CO₂ hydrogenation and electrochemical reduction. Both systems were found to be thermodynamically stable and catalytically active, with SeN₄Gr consistently outperforming TeN₄Gr. The activation of H₂ molecules was facilitated by hybridization between H₂–σ states and vacant Se-4p/Te-5p orbitals. For CO₂ hydrogenation, the formate (HCOO) intermediate was confirmed as the key step, with maximum activation barriers of 0.69 eV (SeN₄Gr) and 1.46 eV (TeN₄Gr). Product desorption was also more favorable on SeN₄Gr (ΔG = − 0.79 eV) compared to TeN₄Gr (ΔG = − 0.22 eV), further underscoring Se’s catalytic advantage. In electrochemical reduction, SeN₄Gr again exhibited lower energy barriers: the initial hydrogenation step required only 0.14 eV versus 0.24 eV for TeN₄Gr. Subsequent radical and product-release steps proceeded exothermically, with Se stabilizing key intermediates such as CH₂O·, thereby enabling efficient methanol formation. These results confirm that SeN₄Gr offers a more favorable thermodynamic landscape and higher catalytic efficiency than its Te analog. Together, the findings establish chalcogen-doped, nitrogen-modified graphene nanoflake as a promising class of non-metal, biocompatible catalysts for CO₂ conversion.

By selectively generating C₁ products such as formic acid and methanol—both high-value, liquid fuels and industrial feedstocks—these systems combine catalytic efficiency with practical handling advantages. Beyond their current performance, the insights presented here provide a conceptual framework for tailoring defect-engineered, chalcogen-based graphene materials as sustainable, cost-effective alternatives to traditional metal catalysts in future CO₂ conversion technologies.

Computational details

. Quantum chemical calculations and reaction mechanisms were performed using density functional theory (DFT) as implemented in the Gaussian 09 software package^27,53, employing the gradient-corrected Perdew-Burke-Ernzernhof (PBE) exchange-correlation functional⁵⁴. The Gaussian 6-31G+(d, p) basis set for H, C, N, O, and Se atoms^55,56, while LANL2DZ57(Los Alamos National Laboratory 2-double-ζ) was used for Te .To account for dispersion interactions and van der Waals effects, the DFT + D3 method with Grimme’s scheme was applied. This approach was crucial for accurate structure determination, geometry optimization, and energy calculations^57,58,59.

The catalyst models were constructed from a finite graphene nanoflake (C₂₀H₁₄ cluster), in which edge carbon atoms were hydrogen-passivated to remove dangling bonds. Di-vacancy defective graphene was generated by removing two carbon atoms, and nitrogen doping was introduced by replacing the four carbons surrounding the vacancy with pyridinic N atoms (N₄Gr). Embedding Se or Te into this site yielded the SeN₄Gr and TeN₄Gr structures (Fig. 1).

Vibrational frequency calculations were used to confirm local minima and transition states, and to evaluate thermodynamic parameters, including enthalpy (ΔH) and Gibbs free energy (ΔG). All values were obtained at 298.15 K and 1 atm, with zero-point energy (ZPE) and thermal corrections included.

The DOS and PDOS were obtained directly from the formatted checkpoint (.fchk) files generated by Gaussian 09. These were post-processed with a custom Python script to extract orbital contributions and plot mid-gap–referenced DOS/PDOS profiles (HOMO–LUMO midpoint set to 0 eV), using a Gaussian broadening of 0.20 eV and an energy step of 0.01 eV. Since all investigated systems converge to closed-shell ground states, the DOS/PDOS are unpolarized.

To validate the structural stability of defective nitrogen-doped graphene nanoflake, the formation energy (E_form) was calculated according to Eq. (1):

$$E_{{form}} = E_{{N4Gr}} - m\mu _{C} - n\mu _{N}$$

(1)

Here, E_N4Gr is the total energy of the defective nitrogen-doped graphene nanoflake, while m and n denote the number of carbon and nitrogen atoms, respectively. The chemical potential of nitrogen (µ_N) was defined as half the total energy of an isolated N₂ molecule, and the chemical potential of carbon (µ_C) was derived from pristine graphene.

For the chalcogen-doped systems (XN₄Gr, where X = Se or Te), the formation energy was determined using Eq. (2):

$$E_{{form,XN4Gr}} = E_{{XN4G}} - E_{{N4Gr}} - \mu _{X}$$

(2)

where E_XN4G is the total energy of the XN₄Gr structure and µ_X represents the chemical potential of the chalcogen atom. For Se and Te, µ_X was obtained from their respective bulk crystalline phases, a standard reference in formation energy calculations.

The adsorption energy (E_ads) was calculated to evaluate the interaction strength between graphene nanoflake and the reactants. It is defined as:

$$E_{{ads}} = E_{{G + reactants}} - \left( {E_{{reactants}} + E_{G} } \right)$$

(3)

In Eq. (3), E_G+reactants is the total energy of N₄GrX with adsorbed reactants, E_G is the energy of N₄Gr, and E_reactant represents the energy of free reactants (e.g., Se, Te, H₂, CO₂, H₂−CO₂). A negative E_ads indicates an energetically favorable and exothermic adsorption process. Additionally, natural bond orbital (NBO) analysis was performed to evaluate atomic charges and charge-transfer values, providing deeper insights into the interactions and stability of the system.

Fig. 1

The optimized structures and their corresponding TDOS/PDOS of (a)N₄Gr, (b) SeN₄Gr, and (c) TeN₄Gr, referenced to mid-gap (0 eV). Black linesrepresent TDOS, orange lines the carbon nanoflake contribution, blue lines thenitrogen contribution, and purple/green lines the selenium and telluriumcontributions, respectively.

Fig. 2

Electrostatic potential (ESP) maps, The ELF plots, and the LOL plots of (a,d,g) N₄Gr, (b,e,h) SeN₄Gr, and (c,f,i) TeN₄Gr. The iso surface value of ESP maps is set to be 0.004 a.u.

Fig. 3

Top and side views of (a,b) H₂ adsorption, (c,d) CO₂ adsorption, and (e,f) coadsorption of H₂ and CO₂ on Se and Te supported on N₄Gr, respectively.

Fig. 4

Gibbs energy diagram (thermodynamics, in eV) and structures corresponding to the reaction path followed by the CO₂ hydrogenation to formic acid in the presence of (a) SeN₄Gr and (b) TeN₄Gr.

Fig. 5

Total and projected density of states (DOS/PDOS) for (a) SeN₄Gr before (top) and after (bottom) CO₂ adsorption, and (b) TeN₄Gr before (top) and after (bottom) CO₂ adsorption. Black lines represent TDOS, orange lines the carbon nanoflake contribution, gray lines represent the CO₂ molecule contribution, blue lines the nitrogen contribution, and purple/green lines the selenium and tellurium contributions, respectively.

Fig. 6

Gibbs free energy diagram (thermodynamics, in eV) and structures in the electrochemical CO₂ reduction mechanism: (a) SeN₄Gr; and (b) TeN₄G.