Morphological Attributes Govern CO2 Reduction on Mesoporous Carbon Nanosphere with Embedded Axial Co-N5 Sites

Although single-atom catalysts (SACs) have been widely employed in the CO2 reduction reaction (CO2RR), the understanding regarding the effect of morphological attributes on catalytic performance are still lacking, which prevents the rational design of high-performance catalysts for electrochemical CO2RR. Here, we developed a novel catalyst with axial Co-N5 sites embedded on controllable mesoporous carbon nanosphere with different graded pore structures. Benefiting from the precise control of porosity, the influence of morphological attributes on catalytic performance was well revealed. In situ characterization combined with density functional theory (DFT) calculations revealed that axial N-coordination induced local d-p orbitals coupling enhancement of Co with oxides and the optimal pore size of 27 nm promoted the interfacial bonding characteristics, which facilitate both the COOH* generation and CO desorption. Consequently, A superior selectivity of nearly 100% at -0.8 V vs. RHE and commercially relevant current densities of >150 mA cm could be achieved, and a strikingly high turnover frequency of 1.136*10 h−1 at -1.0 V has been obtained, superior to the most of Co-based catalysts.


Introduction
The ever-increasing energy demand accompanied with the excess emissions of greenhouse gases leads to the energy crisis and global warming [1][2][3] . Electrochemical reduction of CO2 provides a viable approach to convert greenhouse gas into valueadded products, utilizing renewable energy [4][5][6][7] . The key point of CO2RR is the construction of electrocatalyst [8][9][10][11] . Among various candidates, atomically dispersed metal catalyst has received extensive attention and get rapid development for their high atomic utilization, low-cost manufacturing, powerful metal-support interaction and modulating coordination environment, especially for single-site catalyst [12][13][14][15] . Despite extensive efforts by researchers, precise control of morphology still remains great challenging. The lacking understanding about the influence of morphological attributes on catalytic performance blocks the rational design of high-performance catalysts for CO2RR 16 .
For single-site catalyst, the selectivity of electrochemical CO2 reduction is closely related to the modulation of chemical environment in catalyst, the activity, meanwhile, is limited to the atomic utilization derived from the partial exposure of active atoms [17][18][19] . Thus, simultaneously engineering the electronic modification and morphological structure for tuning chemical environment and promoting exposure of active sites open new opportunities to outperform the well-established catalyst. Considering the unfilled 3d orbitals in transition metal-nitrogen co-decorated nanocarbon could enhance the charge transfer between M center and COOH* intermediate, it is a reasonable method to increase the d-states of M sites to reduce the energy barrier of COOH* formation, which is commonly believed to be an essential pathway for CO production [20][21][22][23] . And thus the selectivity of catalyst could be further optimized. Moreover, the increase of overall activity could be addressed by further improvement to achieve the most efficient utilization of active sites via morphological and structural control of catalyst. Hollow nanostructures combining functional shells and inner voids offer low density, large surface area. Meanwhile hierarchical hollow structures organized by secondary subunits with high-energy facets ensure highly dense exposure of active sites 24 . Thus the activity would be enhanced. Rationally increasing the complexity of chemical environment and morphological nanostructures will likely endow them with new functionalities and great significance for all kinds of catalytic reaction and other emerging technologies.
Herein, we reported a novel CO2RR electrocatalyst, in which atomically dispersed Co-N5 with axially tractive coordination was embedded in adjustable mesoporous carbon nanosphere. The axially tractive Co-N5 embedded on tortile N-doped carbon induced unique electronic structure, and the mesoporous foam structure promoted the contact between reactant and highly exposed active sites. Benefiting from the precise control of porosity, the influence of morphological attributes on catalytic performance was well revealed. As a result, the superior activity and selectivity have been reached. The rational designed catalyst showed a 96.1% CO faradaic efficiency (F.E.) at -0.8 V vs. RHE, and a remarkably high turnover frequency of 1.136*10 4 h −1 at -1.0 V has been obtained, which exceeded the most of Co-based catalysts. Commercially relevant current densities of >150 mA cm -2 also have been achieved in flow cell. The exceptional CO2RR performance originates from the active moiety with axial oxygen traction and precise control of porosity. Local structures and coordination environment was identified through X-ray absorption fine structure (XAFS) spectroscopy. DFT calculation revealed that axial traction of nitrogen atom produced a stronger d-p orbital coupling, leading to a reduced free energy of COOH*, further resulting in a promoted COOH* generation and CO desorption as well as suppressing the competitive hydrogen evolution reaction. Moreover, In-situ FIIR revealed the interfacial bonding characteristics, the optimal pore size of 27 nm has been obtained, which enhanced both the adsorption and conversion of intermediate for CO2RR. This work provides an efficient way to tailor and optimize the catalytic activity and selectivity by atomic-scale design and structural control.

Synthesis and characterizations of Co-MSPNC
As illustrated in Fig. 1a, silica colloid nanoparticles was used as a porosity-inducing template, polyaniline (PANI) was used as C and N sources and cobalt chloride was selected as metal precursors. Silica colloid nanoparticle and polyaniline self-assembled into nanospheres, where Co 2+ could be coordinated with the amine or imine nitrogen atoms (denoted as SiO2@Co-PANI). Then, after pyrolysis under an inert atmosphere at high temperature, the coated Co-PANI layers on the resultant SiO2@Co-PANI nanospheres carbonized into N-doped carbon and the coordinated Co 2+ ions were embedded on the walls(denoted as SiO2@M-N-C). Eventually, the annealed product was leached in NaOH ethanol aqueous solution to remove the SiO2 template, resulting in the formation of isolated single Co atoms immobilized on mesoporous N-doped carbon nanospher (denoted as Co-MSPNC). A uniform spherical structure containing silica colloid nanoparticles could be observed (Fig. 1b-c, Supplementary Fig. 1-3 The average pore size is about 27 nm and the thickness of wall is at the range of 4nm to 10nm, which depends on the measured positions ( Supplementary Fig. 6). No significant nanocrystals or nanoparticles (NPs) were observed throughout the porous nanocarbon, which agrees with the selected area electron diffraction (SAED) observations ( Supplementary Fig. 5). The enlarged HRTEM image (Fig. 1h) displayed lattice distortion defects characteristic, which might be attributed to the coordination of single cobalt atoms with nitrogen. In order to confirm the existence Co single sites, the aberration-corrected atomic-resolution high-angle annular dark field scanning transmission electron microscopy (HAADF-STEM) was carried out (Fig. 1i). Serried isolated bright dots evenly distributed on the carbon matrix, which demonstrates that the metal were atomically dispersed. The isolated Co atoms have been marked with green circles for clear observation, and no metal nanoparticles had formed. Notably, some Co atoms were anchored on curving support, which might lead to unique electronic modulation (inset of Fig. 1i). The elemental mapping acquired by the energydispersive X-ray spectroscopy (EDS) manifests homogeneous distribution of Co, N and C over the entire architectures (Fig. 1g, Supplementary Fig. 7), indicating that the single-site anchoring of Co on N doped carbon. The Co loadings are 0.21 wt% through the inductively coupled plasma mass spectrometry (ICP-MS) analyses. Additionally, Xray diffraction (XRD) was conducted, and two characteristic diffraction peaks at approximately 24° and 44° were displayed in the XRD pattern ( Supplementary Fig. 8) and were assigned to the graphite (002) and (100) planes, which further confirms the atomic dispersion of Co with non-formation of massive metallic and oxide crystals on carbon supports 25 . This result is according with the phenomenon observed in HAADF-STEM image. The relatively high D band-to-G band intensity ratio (ID/IG) was calculated to be 1.09 in Raman Spectrum ( Supplementary Fig. 9). This result implied its great degree of structural disorder within the catalyst, which may be attributed to the anchor of single cobalt atoms with nitrogen on carbon support.
In order to probe the role of the porosity in the CO2RR process, we established various Co-N single-site catalyst with different porosity. Different sizes of silica colloid nanoparticles was used as a porosity-inducing template. As shown in Fig. 2a  between Co metal foil and Co3O4 references, indicating that the Co average valence is between 0 and +2.6 ( Fig. 3b). The XANES spectra of Co-MSPNC contain similar characteristic features to those of CoPc, while differences lie mainly in the relative intensities and edge positions (labelled A-D in Fig. 3b). The inset of Fig. 3b highlights the pre-edge features of the XANES spectra, the integrated pre-edge peak (Peak A at approximately 7708 eV) intensity of CoPc is weaker than that of Co-MSPNC. The increases in the pre-edge peak intensity in Co-MSPNC indicates the distorted D4h symmetry 28 . Compared with that observed for CoPc, the transition of Co 1s to 4pz orbitals (peak B) as a fingerprint for square-planar M-N4 moieties was very weak in Co-MSPNC, which further confirms that Co-MSPNC was not a typical square-planar M-N4 structure 28,29 . Additionally, the peak C and peak D are attributed to the 1s to 4px,y transitions and multiple scattering processes, respectively 29 . Co-MSPNC showed larger relative intensity ratios (IC/ID) than CoPc, further implying that the Co-N4 structure was partially distorted due to an additional ligand 30 . Increased IC/ID has been found to correlate well with enhanced oxygen reduction reaction activity. Co-MSPNC is therefore expected to perform good electrocatalytic activities Furthermore, the Co K-edge Fourier transformed extended XAFS (FT-EXAFS) spectrum of Co-MSPNC was conducted. In Fig. 3c, the Co-Co scattering path vanishing at about 2.18 Å further suggests the atomic-level dispersion of the Co species rather than the crystalline metal structure in the Co-MSPNC. The major peak at 1.44 Å can be assigned to the backscattering of Co-N coordination 31 , further confirming the metalnitrogen coordination. Meanwhile, wavelet transform (WT)-EXAFS was conducted to probe the adjacent properties of metal atoms, allowing more intuitive comparisons on the radial distance k-space resolution. As shown in Fig. 3e, the WT signal related to Co-Co bond at 6.7 Å −1 in Co foil is not detected in Co-MSPNC, further verifying the absence of Co-containing nanoparticles. The observed intensity maxima at 3.3 Å −1 a High-resolution of N 1s, b Co K-edge XANES spectra, c Co K-edge k 3 -weighted EXAFS spectra, d Fitting curves of the FT-EXAFS spectra for Co-MSPNC, The greyish-white, blue and purple refer to C, N and Co atoms, respectively. e Co K-edge WT-EXAFS spectra. appeared in Co-MSPNC is respectively associated with the Co-N paths. Therefore, Co species was further confirmed as atomically dispersed and co-coordinated with N atoms.
To determine the local structure of Co-MSPNC, we established a theoretically optimized bimetal-nitrogen model based on density functional theory (DFT) and fitted FT-EXAFS curves. The EXAFS fitting curves at R-space and k-space are both consistent well with experiment spectra of Co-MSPNC (Fig. 3d). The main peak is the nearest neighbor shell, denoted as first Co-N coordination shell, in the radial distance of about 2.01 Å. The coordination numbers of Co-N are 5.3, which indicates the structure of Co-MSPNC is Co-N5. The calculated Co-N interaction path fits well with the experimental spectra and the structural parameters of coordination numbers (CN), bond distances (R), and Debye-waller factor (σ 2 ) are listed in Supplementary Table 1. Based on the distorted D4h symmetry and high fitting qualities (low R-factor) of the EXAFS, the isolated Co sites are determined to be inlaid in N-doped carbon with four square-planar Co-N bonds and one axial Co-N bond in a distorted structure.

Electrocatalytic CO 2 reduction performance
Inspired by the unique structural features, we tested the electrocatalytic CO2RR performance of Co-MSPNC in a two-compartment H-cell with a common three- cm -2 ) and Co-MSPNC-7 (-8.80 mA cm -2 ) (Fig. 4c), respectively. Meanwhile, a low H2 F.E. of < 10% was achieved on Co-MSPNC-27 at -0.5 V to -1.0 V, much lower than that of Co-MSPNC-43 and Co-MSPNC-7 (Fig. 4d), reinforcing the pivotal impact of axial Co-N coordination and porous morphology on enhance CRR performance. The enhanced activity of the catalyst is further disclosed by the calculation of the turnover frequency (TOF) based on the number of metal atoms involved in the reaction (Fig. 4e).
The Co-MSPNC-27 showed larger calculated TOFs than that of Co-MSPNC-43 and Co-MSPNC-7 across the entire potential window from -0.5 to -1.0 V. Notably, such a high TOFs of Co-MSPNC reported here was among the best relative to all of the previously reported Co-based CRR electrocatalysts (Fig. 4f, Supplementary Table 6) [31][32][33][34][35][36][37][38][39][40][41] . In view of the low solubility of CO2 aqueous electrolytes, the CO2RR performance of a gas diffusion electrode (GDE) with loaded Co-MSPNC-27 was further studied, where CO2 molecules diffused directly to the catalyst from the back of GDE, thus completely removing the mass-transport limitation of solution phase reaction. As shown in Fig 4g, the Co-MSPNC-27 reached a high current density of 150 mA cm −2 at −1.76 V with a CO FE of >90%. Co-MSPNC presents a low Tafel slope of 92.3 mV dec -1 for CO generation, which was significantly smaller than that of Co-MSPNC-43 (163.8 mV dec -1 ) and Co-MSPNC-7 (171.9 mV dec -1 ) (Fig. 4h). This reveals that the  Supplementary Fig. 14) than control samples. A stable current density at about -18.0 mA cm -2 was maintained during a 30h potentionstatic test of Co-MSPNC-27 (Fig. 4i). The structures of Co-MSPNC-27 after performance testing and stability test was well characterized by TEM, HADDF-TEM, XRD, Raman and XAFs ( Supplementary Fig. 15), which demonstrates the maintaining of mesoporous nanosphere and axial Co-N5 structure. Therefore, obvious enhanced CO2RR performance has been observed on mesoporous carbon nanosphere with the pore size of 27 nm and axial Co-N5 sites.

Intrinsic property of the Co-MSPNC catalysts.
To elucidate the key role axial Co-N played in the active Co-MSPNC for the CO2RR catalytic properties, we adopted the computational hydrogen electrode (CHE) model to calculate the free energy (G) of each intermediate state for the CO2RR. Besides, we calculated formation energies (ΔE) of Co-N5 and Co-N4 structures. The Co-N5 system has the lower energy than Co-N4, so the atoms tend to form Co-N5 (Supplementary Fig.   16). According to the computational hydrogen electrode (CHE) model, the stepwise atomic structures of the CO2RR process for the electrochemical reduction of CO2 to CO and the HER process for the reduction of H + ions to H2 are shown schematically in Fig. 5a. The free energy diagrams for the CO2RR to CO and HER at U = 0 V (vs RHE) are shown in Fig. 5b, presenting the free energy of each intermediate state on the metal atom sites in flat Co-N4-C and axial Co-N5 sites, and each step represents the transfer of one proton-electron pair. As shown, the free energy barriers (ΔG) for the whole reaction are 1.12 eV and 1.56 eV for Co-MSPNC and Co-N-C, respectively. The reason for this difference is that axial Co-N5 coordination can significantly decrease the free energy for the *COOH reaction, which is typically thought of as rate-determining step during CO2RR. Furthermore, the increase of free energy for *CO indicates that axial Co-N5 coordination facilitates the desorption of CO, leading to a high CO selectivity.
Additionally, the remarkably high positive free energy change required (0.44 eV) for the formation of adsorbed hydrogen intermediates during the dominant competitive HER (Fig. 5c) suggested the preferred high selectivity toward CO rather than H2.
Therefore, the introduction of axial Co-N5 coordination can not only promote the CO2RR but also maintain the inhibition of the competitive side reaction of the HER, which is in good agreement with the experiments Further we also calculated formation energies in axial N atoms to evaluate whether the exposed N atom is active for CO2RR. The free energy for the *COOH reaction in N atom is 7.376 eV, which is much higher than that in Co (0.670 eV). There is a great deal of difference between the two atoms. The almost 10 times of free energy more than Co atom reflect that the exposed N atom is almost non-active for CO2RR ( Supplementary   Fig. 17). In order to gain in-depth insights into the catalytic selectivity and the role of axial Co-N coordination, their electronic structures were further investigated. Firstly, the density of states (DOS), the crystal orbital Hamilton population (COHP) and their integrated-COHP (ICOHP) values for Co 3d orbitals and the oxide 2p orbitals. There is a negative correlation between the ICOHP value and the strength of d-p interaction, that is, a more negative ICOHP leads to a stronger d-p coupling 42 . For the *COOH step, the *COOH (2p) orbitals are strongly bonded with the Co (3d) orbitals on the Co site in Co-N5 because of their relatively bigger DOS overlap and more negative ICOHP value compared with those in Co-N4. Therefore, the interaction between the Co atom and C atom in COOH* is the stronger than Co-N4-C, consisting with its lower free energy. For the *CO step, the Co 3d electronic state in Co-N5 is further from EF than that of Co-N4 (Fig. 5d) and the ICOHP value of Co-MSPNC is significantly reduced to −1.466 eV (Fig. 5e), in agreement with its relatively low free energy to promote the CO desorption on axial Co-N5 sites. Additionally, to better visualize the electron transfer during the catalytic process, valence charge density difference (VCDD) and bader charge transfer number (Δe) were calculated. Fig. 5f-i presents axial traction of nitrogen atom produced an electron-pushing effect, leading to more significant charge polarization on Co-N5 sites. The larger charge transfer number demonstrates that axial Co-N5 sites can much more easily adsorb and bond with COOH than Co-N4-C, while the smaller bader charge for CO* step implies the promoting desorption of CO on axial Co-N5 sites. These results are consistent with DOS and COHP analysis.

Interfacial analysis of Co-MSPNC with different graded pore structures.
Besides the stronger d-p orbitals coupling activated by axially tractive Co-N5 sites, the porous structure also plays a significant role during the electrochemical reduction of CO2. The superior CO selectivity of Co-MSPNC-27 suggests the variations in the porosity could affect the eletrocatalytic performance, and 27 nm is the optimal pore size.
Considering the porosity directly affect the surface properties, the wetting properties of three samples were investigated by measuring the dynamic (advancing and receding) contact angles of the thin layers. The surface energies (Fig. 6a) were calculated from the contact angle hysteresis. As shown in Fig. 6b, the Co-MSPNC-27 catalyst was the most hydrophobic and presented the lowest hysteresis, resulting in the lowest surface energy (γgs,tot=~19.47 mJ m -2 ) 43,44 . Temperature-programmed desorption was employed to study the CO2 adsorption strength of the catalysts (Fig. 6c).While only one desorption peak was present in case of the Co-MSPNC, CO2 desorption from the porous samples happened in two distinct steps. We assume that the first step at a lower temperature (< 400 K) is connected to the adsorbed CO2 on the outer surface of the catalyst layers, and the second step is linked to the desorption of CO2 from the inner surface of the pores.
The temperature of the desorption from the pores was the highest (516 K) for Co-MSPNC-27, implying the highest CO2 adsorption strength. Also in the CO2 conversion, we assume a link between the CO2 adsorption and reduction activity. Considering that the studied Co-MSPNC catalysts are practically identical in terms of Co-N-containing moieties, we suggest that the distinct CO2 adsorption strengths arise from the differences in the pore structure. In this case, the strain induced on the CO2 adsorption sites by the differently curved surfaces is optimal for the 27 nm pores. In order to further understand the reaction pathway and verify the difference in adsorption properties towards reactive intermediates during the CO2 conversion, in situ Fourier transform infrared (FTIR) spectroscopy was performed. and Co-MSPNC-7, respectively, indicating a weaker adsorption and more difficult process to produce *CO. The adsorbed *COOH species, which are regarded as the crucial intermediates during the reduction of CO2 to CO, begin at a low potential of approximately -0.05 V vs. RHE for all three samples. As the conversion of *COOH to *CO is thermodynamically downhill, we speculate that Co-MSPNC-27 shows stronger *COOH adsorption than Co-MSPNC-43 and Co-MSPNC-7.

Discussion
In summary, we report the rational design of an axial Co-N5 sites embedded on controllable mesoporous carbon nanosphere with different pore sizes. Single site was confirmed by aberration-corrected electron microscopy and X-ray absorption spectroscopy. Benefiting from the precise control of porosity, the influence of morphological attributes on catalytic performance was well revealed. The optimized Co-MSPNC-27 possessed a high CO2RR catalytic selectivity and activity due to its unique electronic modification and morphological structure. A superior selectivity of nearly 100% at -0.8 V vs. RHE could be achieved, and a remarkably high turnover SiO2@M-PANI was obtained centrifugation, washing and dring. Next, sample was transferred into a quartz tube, which was then heated to 900 °C kept for 3 hours at a heating rate of 5℃/min in flowing argon atmosphere. The sample was collected and socked in 3 M NaOH H2O/EtOH (v: v = 1: 1) solution at 100 °C overnight, and then washed with copious deionized water and ethanol, dried at 70 °C under vacuum for overnight. The obtained material was directly used without further treatment.

Synthesis of Co-MSPNC-7:
The Co-MSPNC-7 was prepared with the same procedure as that of Co-MSPNC-27 except the using of silica colloid (The grain diameter is 7nm).

Synthesis of Co-MSPNC-43:
The Co-MSPNC-7 was prepared with the same procedure as that of Co-MSPNC-27 except the using of silica colloid (The grain diameter is at the range of 40-50nm).

Electrode
A three-electrode configuration was used for CO2RR measurements. To prepare the working electrode, 10 mg of catalyst was dispersed into 100 μL of Nafion solution (0.5 wt.%) and 900 mL of ethanol by ultrasonic for 2 h and stirring overnight to form a homogeneous ink, and then 100 μL of as-prepared ink was dropped onto the carbon paper with reactive area of 1 × 1 cm 2 , with a loading amount of 1.0 mg cm -2 . Platinum foil was served as counter electrode, and an Ag/AgCl(filled with saturated KCl solution) electrode was acted as reference electrode. All potentials were recorded against an Ag/AgCl reference electrode and converted to versus RHE with iR corrections.

Characterizations
The morphologies of samples was performed by field-emission scanning electron The X-ray absorption find structure spectra (Fe K-edge) were collected at 1W1B station in Beijing Synchrotron Radiation Facility (BSRF). The storage rings of BSRF was operated at 2.5 GeV with an average current of 250 mA. Using Si(111) double-crystal monochromator, the data collection were carried out in transmission/fluorescence mode using ionization chamber. All spectra were collected in ambient conditions. In situ FTIR spectra were acquired using a PerkinElmer spectrum100. All pectroscopic experiments are conducted at a 4 cm -1 spectral resolution and 16 scans.

XAFS Analysis and Results:
The acquired EXAFS data were processed according to the standard procedures using the ATHENA module implemented in the IFEFFIT software packages. The k 3weighted EXAFS spectra were obtained by subtracting the post-edge background from the overall absorption and then normalizing with respect to the edge-jump step.
Subsequently, k 3 -weighted χ(k) data of Fe K-edge were Fourier transformed to real (R) space using a hanning windows (dk=1.0 Å -1 ) to separate the EXAFS contributions from different coordination shells. To obtain the quantitative structural parameters around central atoms, least-squares curve parameter fitting was performed using the ARTEMIS module of IFEFFIT software packages.

Electrochemical measurements
Electrochemical measurements were conducted in a traditional H-type cell with a cation exchange membrane (Nafion 117) on electrochemical workstation (Ametek 1470E).
The working electrode was carbon paper (1 × 3 cm 2 ) with electrocatalyst. Platinum foil was served as counter electrode, and an Ag/AgCl electrode was acted as reference electrode. The CO2 or Ar saturated electrolyte was formed by bubbling CO2 (or Ar as a control) into electrolyte (0.5 M KHCO3) at least 30 min before CRR test. During the electrochemical measurements, the rate of gas maintained at 20 mL min -1 . The polarization curves were measured with a constant rate of 5 mV s -1 . The electrochemically active surface area (ECSA) was tested between -0.9V ~ -1.0 V (vs. Ag/AgCl). Electrochemical impedance spectroscopy (EIS) was performed in a frequency range from 10 5 Hz to 0.01 Hz. Gaseous and liquid products were quantified by gas chromatograph (Agilent GC7820) and 1 H nuclear magnetic resonance (BRUKER AVNCE III500M), respectively.

Density Functional Theory (DFT) Calculations
Spin-polarized DFT calculations were performed using the Vienna ab initio simulation package (VASP) code. Projector augmented wave (PAW) pseudopotential was used to describe the core electrons and a plane wave basis set with a kinetic energy cutoff of 520 eV was used to expand the wave functions. Electronic exchange and correlation were described within the framework of generalized gradient approximation (GGA) of the Perdew, Burke and Ernzernhof (PBE) exchange-correlation functional.
All electronic structure calculations were performed spin-polarized. The convergence of energy and forces were 1×10 -5 eV and 0.005 eV Å -1 , respectively. The dispersioncorrected DFT-D3 method was adopted to describe the van der Waals (vdW) interactions. The entropic correlation and zero-point energy have been considered in Gibbs free energy calculation. And the implicit solvation model (VASPsol) was used to examine the water-solvent effect on CO, COOH, H adsorbed on the N-doped graphene which adds implicit aqueous electrolyte into the system.
To simulate N-doped graphene supported Co-N4 and Co-N5 catalysts, a 7×7 graphene supercell with 98 carbon atoms was chosen. Periodic boundary conditions were employed along x and y directions with a vacuum region of 20 Å between graphene layers. The supercell is large enough to avoid the interaction between the adsorbate and its image. A 5×5×1 Monkhorst-Pack grid of k-points was used to sample the first Brillouin zones of the surfaces for structural optimizations, while a 6×6×1 grid of k-points was used for electronic states analysis. All atoms are fully relaxed during structural optimizations.

Computational Hydrogen Electrode (CHE) Model
The computational hydrogen electrode method developed by Nørskov was used to calculate the free energy of each intermediate state. Theoretically, the reaction pathway in the electrochemical reduction of CO2 to CO following the progress: CO2 (g) + H + + e -+ * → COOH* (1) COOH*+ H + + e -→ CO* + H2O (2) CO* → CO (g) + * in which the asterisk refers to the clean catalysis, and the asterisk-marked one refers to the species that adsorbed on the activity center. In the CHE, the chemical potential of a proton-electron pair, the free energy of (H + + -) equals to 1/2 H2(g) for the reversible hydrogen electrode (RHE), that G[H + + − ] =1/2G[H2] -eU, where U is the electrochemical potential relative to RHE. As the initial state (gaseous CO2 freely above an empty substrate) for the reference, the Gibbs free energy changes of interest can be given by, The above mentioned Gibbs free energy (G) of a species is calculated according to G = EDFT + EZPE -TS + ΔEsol + ΔHT (11) where EDFT is the DFT calculated energy change, EZPE is the zero-point energy, ΔEsol is the solvation energy, ΔHT is the temperature enthalpy correction from 0 to T K, The reversible potential of the overall half-cell reaction: CO2 + 2H + + 2e − → CO + H2O, which energy calculated from DFT is -0.09 eV vs RHE, in agreement with the experiment value -0.11V vs RHE.
The corrections and calculated free energy of gas-phase were listed in Supplementary  Table 3.

Data Availability
The data that support the plots within this paper and other findings of this study are available from the corresponding authors upon reasonable request.