Copper single-atoms embedded in 2D graphitic carbon nitride for the CO2 reduction

We report the study of two-dimensional graphitic carbon nitride (GCN) functionalized with copper single atoms as a catalyst for the reduction of CO2 (CO2RR). The correct GCN structure, as well as the adsorption sites and the coordination of the Cu atoms, was carefully determined by combining experimental techniques, such as X-ray diffraction, transmission electron microscopy, X-ray absorption, and X-ray photoemission spectroscopy, with DFT theoretical calculations. The CO2RR products in KHCO3 and phosphate buffer solutions were determined by rotating ring disk electrode measurements and confirmed by 1H-NMR and gas chromatography. Formate was the only liquid product obtained in bicarbonate solution, whereas only hydrogen was obtained in phosphate solution. Finally, we demonstrated that GCN is a promising substrate able to stabilize metal atoms, since the characterization of the Cu-GCN system after the electrochemical work did not show the aggregation of the copper atoms.


INTRODUCTION
The depletion of fossil fuels and the increase of the atmospheric CO 2 levels encourage the continuous search for alternative energy sources and routes to transform CO 2 into fuels or other valueadded chemicals, which would allow us to move the C-economy from a linear-to a cyclic-one 1,2 . Awkwardly, CO 2 is an extremely stable molecule and its reduction involves a multi-electron transfer process and a diversity of intermediates and products, such as CO, methane, methanol, formic acid, etc.
In the search for highly active, selective, and durable catalysts that will make possible the commercialization of this technology, the dispersion of single metal atoms on a substrate is taking much interest in the last years [3][4][5] . The most commonly explored singleatom catalysts (SACs) for CO2RR are transition metal elements due to their electronic configuration (vacant orbitals and d electrons) that can facilitate the M-CO 2 bonding as well as the desorption of the reduction products [6][7][8] . On one hand, the use of SACs allows maximizing the atomic efficiency of the metal catalysts (the wellknown atom economy concept of Green Chemistry), which often have shown unexpected activity [9][10][11] . On the other hand, it allows controlling the nature of the active sites, which could help to improve the products selectivity, since there is no presence of different phases. Unfortunately, these materials are often unstable and tend to aggregate into clusters/nanoparticles because of their high surface free energy, which enhances their mobility on the support surface 12,13 . In particular, under catalytic conditions, the interaction between the single atoms and the substrate decreases due to new interactions established between the single atoms and the reactant, favoring the migration of the single atoms on the support 6,14,15 . Therefore, finding a suitable and effective method to stabilize single atoms is a crucial requirement for the use of SACs.
The most used approaches in the literature for the stabilization of SACs are: (i) the establishment of strong metal−support interactions, and/or (ii) the confinement of the single atoms into a limited space, thus decreasing their mobility 16,17 . In this sense, graphitic carbon nitride (GCN) presents properties that make it an interesting support for SACs. It is a polymeric organic semiconductor that consists of 2D layered structures that can be easily synthesized from several precursors such as melamine, urea, thiourea, etc. 18,19 and shows a high physicochemical stability and non-toxicity. In its fullypolymerized configuration, GCN presents structural cavities (six-fold N cavities) where SACs could be confined and stabilized through strong metal-nitrogen interactions 20,21 , although the formation of less regular polymorphs has to be carefully evaluated in order to investigate the correct coordination of the metal atoms. Compared with nitrogen-doped graphene-related materials, which are frequently used for the dispersion of SACs 7,22-26 , GCN contains a higher amount and uniform type/arrangement of nitrogen species that offer abundant and more precisely defined coordination sites for single atoms [27][28][29] . No less important it is the fact that it has a relatively small band gap (2.7 eV; 30 ) therefore, it is able to absorb the light, reducing the external input necessary for the CO2RR 31 . Alone, GCN shows a high electron−hole recombination rate that limits the overall photocatalytic efficiency of the single-phase GCN 32,33 . However, combined with SACs, the separation of photogenerated charge carriers can be significantly enhanced since the metal atoms can act as electron traps as well as catalytic sites 28,[34][35][36] .
Among the different metals, copper is of particular interest since is the only one able to reduce CO 2 to relatively high-hydrocarbons (up to 6 C atoms) with a considerably high Faradaic efficiency (FE) [37][38][39][40][41][42][43] . Therefore, in order to convert CO 2 into C−C hydrocarbons, copper has to be included in the catalyst formulation. A few works about the study of Cu single atoms dispersed in GCN can be found in the literature 28,36,44 . For example, Quiao and coworkers reported a theoretical and experimental study on the electro-CO2RR in an aqueous medium catalyzed by Cu atoms dispersed in a GCN structure 44 . They obtained H 2 (>50%) and CO as main products, and low amounts of CH 4 , C 2 H 4 , C 2 H 6 , CH 3 OH, and C 2 H 5 OH. Wang et al. 1 have recently prepared a single-atom copper-modified GCN that has shown to be a competing photocatalyst for the photoconversion of CO 2 into CO, CH 3 OH, and CH 4 in an aqueous solution 36 . Moreover, Li and coworkers have reported a study of Cu atoms supported on crystalline GCN as photocatalyst for the selective reduction of CO 2 into CO 28 . However, the stability of this kind of materials has not been investigated.
In the study reported herein, we have modified the exfoliated GCN with copper single atoms in order to exploit the catalytic activity of the Cu-GCN hybrid system towards the electrochemical CO2RR, studying the effect of the visible light on the catalytic process. Both the substrate and the catalyst have been characterized by a large variety of techniques in order to get a comprehensive understanding of their physicochemical properties. The experimental observations have been combined with density functional theory (DFT) simulations of several models of 2D GCN, accounting for different coverages of Cu atoms, and of the XPS spectra to achieve an atomistic description of the systems under investigation. Regarding the electrochemical CO2RR, the dependence of the reaction products on the electrolyte has been studied by using quantitative techniques (such as gas chromatography and 1 H-NMR) coupled with rotating ring-disk electrode (RRDE) measurements, so evidencing an important role of the electrolyte in the products selectivity. In addition, the stability of the Cu-GCN hybrid system after the electrochemical work has been investigated to prove that GCN is able to stabilize the Cu atoms, thanks to the confinement in its structural cavities and to the strong copper−nitrogen interaction.

Physicochemical characterization
Two-dimensional (2D) GCN was obtained through an exfoliation treatment of the bulk material with H 2 SO 4 . A complete comparison of the bulk and exfoliated 2D GCN materials can be found in Supplementary Fig. 1. Subsequently, it was modified with copper atoms (20 wt.%) by using a well-established method for the introduction of single atoms (details in the "Experimental" section) 20,45,46 . The morphology of the GCN-based materials was investigated by scanning and transmission electron microscopies (SEM and TEM) (Fig. 1). The bulk material consists of stacked layers forming agglomerates of some micrometers, as seen in Fig. 1a, c, and in accordance with observations reported in the literature 47,48 . The exfoliated material shows the characteristic layered-like structure with nanosheets of some hundreds of nanometers, as seen in Fig. 1b, d. However, after the exfoliation process, a decrease in the size of the aggregates is observed. The change in the morphology can be attributed to the cleavage of the bridging -NHbetween heptazine units due to the hydrolysis/oxidation during the H 2 SO 4 exfoliation treatment 48,49 , as schematized in Fig. 1e. In addition, GCN materials usually contain structural defects (primary or secondary amines) due to the incomplete condensation of the heptazine units that may facilitate the breakage of the GCN sheets 48,50 . This hypothesis is further investigated below by using other characterization techniques. No copper clusters or nanoparticles are observed in Fig. 1f for the Cu-GCN sample with 20 wt.% of Cu, suggesting the successful formation of copper single atoms, but we cannot exclude the presence of small clusters not detected by the TEM microscope used for the analysis.
The X-ray diffraction (XRD) patterns of 2D GCN and Cu-GCNs show two diffraction peaks at around 27.3°and 12.2°, indexed to the (002) and (100) planes of GCN, respectively (Fig. 2a) 51 . The first peak is characteristic of the interlayer stacking of the conjugated aromatic systems; while the second one is associated with the inplane repeated units, i.e., the distance of the adjacent heptazine units. A shift of the (100) peak toward lower 2-theta values is observed for the exfoliated material (12.2°), compared to the bulk GCN (13.1°), in Supplementary Fig. 1a. This shift indicates an increase of the distance between the adjacent heptazine units, in agreement with the hypothesis of the cleavage of the bridging -NHbetween heptazine units and/or breakage of the GCN sheets due to the hydrolysis/oxidation of the structural defects. This modification of the exfoliated phase is also observed in the DFT simulated XRD spectra for the two models of 2D GCN considered (i.e., linear or alternated, Supplementary Fig. 3a, c), which are reported in Supplementary Fig. 1a: the feature occurring at~12.2°i s indeed generated by an in-plane diffraction, although not for (100) but for (200), (010), or (2-10) planes, depending on the arrangement of the monomers in Liebig's melon. No peaks related to copper nanoparticles are observed, suggesting the presence either of single Cu atoms or small Cu clusters embedded in the GCN structure (detection limit of XRD > 2 nm). In order to determine the maximum amount of single atoms that the GCN substrate can host, samples with 10 and 30 wt.% of Cu were also prepared ( Supplementary Fig. 2a). No peaks associated with copper nanoparticles are observed for Cu-GCN (10). However, when the Cu loading is increased up to 30 wt.%, new peaks attributed to the formation of copper nanoparticles appear. This result suggests that 20 wt.% is an upper limit of Cu single atoms that the GCN substrate can host without forming Cu clusters/ nanoparticles, which is in agreement with the value reported by Yang et al. 52 , but much higher than that usually reported in the literature for this kind of materials 34,[53][54][55] .
In order to investigate the structure and coordination of the Cu atoms, XAS measurements were performed on Cu-GCN (20). The analysis of the XANES region indicates that the Cu atoms have an average oxidation state of around +1.4 ( Supplementary Fig. 2), suggesting the presence of both Cu 2+ and Cu + , due to different Cu−N coordination arrangements. The presence of two different valence states has been already reported for Cu single atoms embedded in carbon nitride or N-doped graphene-related materials 6,53,56 . Fig. 2b shows the Cu K edge Fourier transformed (FT) extended X-ray absorption fine structure (EXAFS) spectra for Cu-GCN (20) and the samples used as references. The absence of an intense peak at 2.23 Å (not-phase shift corrected) suggests the absence of copper clusters (Cu-Cu bonds), confirming the dispersion of single Cu atoms in the GCN matrix 54 . Taking into account the characterization results for the Cu-GCN sample and the results reported in the literature for this kind of materials, the peak at 1.5 Å (not-phase shift corrected) can be related to Cu coordinated to N atoms in the first shell 28,52,53 , although the presence of Cu−O interactions cannot be excluded; while the one at 2.4 Å could be attributed to Cu−C interactions 6 . In order to fit the EXAFS data, we included Cu−N and Cu−C scattering paths, taking into account the most feasible structure of the Cu-GCN sample (Amino I and Amino II in  53 . Regarding the Cu−C and Cu−Cu shells, it is difficult to establish a comparison with the literature since usually only the fit of the first Cu−N shell is reported. The addition of the Cu−Cu shell improved the R f factor (compared to R f = 0.023 without including the Cu−Cu shell). However, it should be noted that the coordination number for this shell is very small. It could be attributed to (i) very small Cu clusters not detected by other techniques or, most probably, to (ii) other single Cu atoms coordinated to an adjacent melem unit, which would be in agreement with the 'high' Cu coverage regions determined by DFT calculations, as described below.
Another approach to get insights on the position and on the coordination sphere of Cu is to combine experimental XPS measurements with DFT simulations. According to the existing DFT calculations, single metal atoms should be mainly located in the empty cavities in between fully-polymerized heptazine units, as a consequence of the in-plane coordination with the nitrogen atoms that induces a significant stabilization 20,21,57 . Fig. 3a reports the Cu 2p XPS region for the Cu-GCN (20) sample. From the fit of the Cu 2p 3/2 region, it is evident that Cu is present in different oxidation states. The first component centred at 933.3 eV could be associated with both the presence of metallic Cu (933.2 eV) or Cu + (933.4 eV) 54,58,59 . Taking into account the XANES analysis, the presence of metallic Cu could be excluded; however, in order to confirm its absence, the Cu LMM Auger peak was analyzed as well (inset of Fig. 3a). The spectrum consists of a broad band between 915 and 920 eV with a maximum centred at 915.6 eV, which is characteristic of Cu + . On the other hand, the second component at 934.2 eV in the Cu 2p 3/2 peak, as well as the intense satellite peak at 940.7 eV, are related to the presence of Cu 2+60 , confirming the presence of Cu in two different oxidation states (Cu + and Cu 2+ ), in agreement with the XANES results. We take this result as a strong indication on the actual structure of GCN in the sample: a fully-polymerized GCN, in fact, is not expected to provide two different chemical environments to the adatoms, since Cu ions preferentially interact with the pyridinic N atoms and not with tri-coordinated N atoms (from our DFT calculations). Thus, we consider a partial polymerization, where the monomers retain some amino functional groups that can also involve Cu atoms in some coordinative bonds (see the melon models below), to be more consistent with the existence of two oxidation states for Cu.
The analysis of the C 1s and N 1s regions of Cu-GCN and GCN is reported in Fig. 3b, c and in Supplementary Tables 2-4. Regarding the C 1s region, both Cu-modified and pristine 2D GCN materials show the component at around 288 eV, which is attributed to C−(N) 3 and is characteristic of GCN 61,62 . The addition of copper atoms results in a 0.7 eV shift of this component toward higher binding energy (BE), which we can rationalize with the fact that copper atoms attract electrons from carbon atoms in Cu−N−C. The other three C 1s components observed for both Cu-modified and pristine 2D GCN systems are related to the adventitious carbon present in the sample holder. The adventitious carbon is visible due to the very thin film of the sample deposited in order to avoid charging effects.
Focusing on the N 1s region of 2D GCN, we deconvoluted the broad experimental band in four features at 398.0, 399.0, 400.7, and 401.7 eV as shown in Fig. 3c (bottom). These are in good agreement with the simulated peaks by DFT calculations reported in Fig. 3e (pristine) and also, with more details, in Supplementary Fig. 3e. We must note that the simulated spectrum is obtained by the overlap of the calculated lines for the alternated and linear models of melon (see Supplementary Fig. 3a, c). As reported in Fig. 3e, for pristine 2D GCN, we assign the four inequivalent nitrogen components to C−N=C, C−N−H 2 , (C) 2 −N−H, and N−(C) 3 species, respectively, in order of growing BE 61 . This theoretical assignment is in contrast with most of those present in the literature and based on chemical intuition, where graphitic and amino features are considered to be inverted 63,64 . However, in contrast with the experiments, our theoretical approach is unambiguous because each calculation specifically targets a single inequivalent atomic component. Moreover, our findings are in line with other previous experimental and theoretical works 61,65,66 .
We now consider what happens after the addition of copper atoms. Again, we fit the N 1s XPS region in Fig. 3c (top) with four components at 398.5, 399.3, 400.2, and 401.3 eV, respectively, thus at slightly different BE than before. Their intensities, especially for those at higher BEs, suggest the formation of a Cu-loaded melon phase. To rationalize the effect of the presence of the Cu atoms, we designed and optimized various melon models with different Cu concentration, reported in Supplementary Fig. 4. Then, we used such models to build two coverage regimes of Cu atoms, namely 'low' and 'high', corresponding to an average load of~10 wt.% and >50 wt.%, respectively. Since the DFT adsorption energy of Cu on the amino sites (in the range between −1.6 and −2.9 eV, see Supplementary Fig. 4) is larger than that on the top sites, where Cu atoms are only physisorbed, in the low coverage regime we have only single Cu atoms at different amino sites, whereas in the high coverage regimes we saturate all the amino sites and further add Cu atoms on top sites. If we only considered the models with lower coverage, we would match the nominal Cu concentration in the experimental samples, but we would not satisfactorily describe the experimental XPS spectrum, because the feature at 400.2 eV would not appear. Therefore, we must consider a mixture of all models, i.e., pristine melon together with the Cu-loaded systems, at both low and high coverage regimes, in the assumption that the sample is not homogeneous.
We built the simulated N 1s XPS spectrum by considering the calculated BEs and their corresponding chemical shift, as described in the computational methods, for all the N species in the investigated models. The experimental spectrum was then fitted by using this mixture of models with non-homogeneous Cu distribution, fixing the position of the peaks and their gaussian broadening in order to extract the optimal weight of the different coverage regimes, as reported in Fig. 3e and Supplementary Table  5. The resulting fit is satisfactory and leads to an extracted average Cu concentration of 14 wt.%, which is close to the total Cu loading determined experimentally (14.8 wt.%).
The presence of Cu atoms at different adsorption sites can explain the observation of two distinct Cu oxidation states in XAS and XPS spectra. However, since our computational setup is not suited to assess this specific aspect, we cannot reach a definitive conclusion. Qualitatively, the simulated N 1s XPS components indicate that Cu atoms directly interact with both pyridinic N and primary amines, whose XPS features are observed at increased BE because of some charge transfer to the metal atoms. Secondary amines and, to an even smaller extent, graphitic N are, instead, only indirectly affected by Cu atoms, as they experience the overall charge re-arrangement through the delocalized π system.

CO2RR activity and products analysis
Once Cu-GCN was chemically and structurally fully characterized, the electrochemical activity towards the CO 2 reduction was tested by linear scan voltammetry (LSV) in two different aqueous electrolytes, 0.5 M KHCO 3 (pH = 6.8) and 0.1 M phosphate buffer solution (pH = 7). The comparison of the LSV under N 2 and CO 2 atmosphere in 0.5 M KHCO 3 solution is reported in Fig. 4a, whereas the results in phosphate buffer solution are shown in Supplementary Fig. 5b. In both electrolytes, a catalytic behaviour for CO2RR is suggested since the onset of the reduction wave is positively shifted by 91 and 48 mV under the CO 2 atmosphere, respectively. In KHCO 3 , a crossing between the lines in N 2 and CO 2 is attained at around −1.5 V, suggesting a preferential production of H 2 from this applied potential. When the measurements are performed under visible irradiation, an increase of the current density and a shift of the onset potential toward less negative values are observed (Fig. 4a, and Supplementary Figs. 5b and 6), indicating a positive effect of the light attributed to the GCN substrate ( Supplementary Fig. 5a).
Several catalyst loadings were tested in order to optimize the current density at different applied potentials (Fig. 4b). Higher differences in the current density with the catalyst loading were observed when the applied potential became more negative. From these results, a catalyst loading of 66 µg cm −2 , was selected for the successive RRDE experiments. The RRDE studies were performed in order to qualitatively determine the CO 2 reduction products in both electrolytes 67,68 . With this aim, the catalystmodified disk was held at different potentials from −0.6 to −1.7 V vs SCE, while the Pt ring potential was scanned from −0.4 to +1 V vs SCE to detect (oxidize) the reduction products. In both electrolytic solutions, only the production of hydrogen was observed in the absence of CO 2 as expected ( Fig. 5a and Supplementary Fig. 7a), since only the wave related to the H 2 oxidation was observed at the Pt ring 69 . The current associated with this process is higher in the phosphate buffer solution, suggesting that the production of hydrogen is favoured in this media. In presence of CO 2 , a completely different CV is observed at the Pt ring in 0.5 M KHCO 3 : an oxidation current between −0.1 and +0.4 V, with peak potential at around +0.2 V, is observed (Fig.  5b), which can be attributed to HCOO − oxidation 67,68 . The attribution of this peak to the formation of HCOO − at the catalyst surface and subsequent oxidation at the Pt ring was confirmed by comparing the CVs in Fig. 5b with those performed at the Pt ring in the presence of added HCOOH (Supplementary Fig. 8). The peak position shifted towards lower potential values upon decreasing the cathodic potential at the disk, in agreement with a similar shift reported in the literature upon increasing the HCOOH concentration 67 . The attribution of this peak to HCOO − is also supported by the 1 H-NMR analysis of the liquid phase after controlled-potential electrolysis experiments (see below). RRDE experiments point out the importance of the electrolyte in the product selectivity and the beneficial effect of the bicarbonate salt with respect to the phosphate buffer towards CO2RR. The analysis of Supplementary Fig. 7 shows that H 2 is the main reduction product in 0.1 M phosphate buffer solution, both under N 2 and CO 2 atmosphere, since no evident re-oxidation peaks of carbon products were detected at the Pt ring. This observation is in strong agreement with previous reports in the literature, where the use of phosphate buffer led to a significant increase in the H 2 production 70,71 .
Measurements on pristine 2D GCN were carried out in both electrolytes at different rotation rates in order to confirm that the activity was due to the presence of Cu (Supplementary Figs. 9 and 10). No peaks related to the oxidation of HCOO − or other C-products were observed. Therefore, the catalytic production of HCOO − under CO 2 atmosphere in KHCO 3 is attributed to the presence of single Cu atoms in the carbon nitride layers (Fig. 6). Conversely, in phosphate buffer solution, the presence of the Cu atoms induces the promotion of the hydrogen evolution reaction, even under the CO 2 atmosphere (Supplementary Fig. 11).
Chronoamperometry tests were performed in 0.5 M KHCO 3 at −1.35 and −1.55 V vs SCE under a continuous flow of CO 2 , in order to analyze the products by 1 H-NMR and gas chromatography and compare them with the RRDE results. A 12% FE for HCOO − production was measured at both potential values for Cu-GCN (20), while the rest of the charge was attributed to H 2 evolution according to RRDE tests and chronoamperometry experiments carried out in a sealed cell under CO 2 . No other liquid products were detected, as shown in Supplementary Fig. 12. The same tests were also performed on pristine 2D GCN and Cu-GCN samples with 10 and 30 wt.% Cu loading. The current measured during the experiment for 2D GCN (0.04 mA/cm 2 ) was very low and almost negligible with respect to that obtained for Cu-GCN (3.6 mA/cm 2 ), as already observed in Supplementary Fig. 5, confirming that GCN does not reduce CO 2 . Regarding the samples with copper, we observed an increase of the current density with the amount of copper (Supplementary Fig. 13a). The current density for the Cu-GCN(10) sample during the electrolysis was very small and no liquid products were detected. The current density displayed by Cu-GCN(30) was higher than that obtained with Cu-GCN (20), as reported in Supplementary Fig. 13b; however, formate was obtained with a lower FE (9%, 40 µmol mg cat. −1 h −1 ) than with the sample with 20 wt% Cu (12%). This could be attributed to the presence of copper oxide nanoparticles. It should be highlighted that the results obtained with our Cu-GCN system are in line with those reported in the literature for similar systems 44 as the main CO 2 reduction product, together with a variety of different products such as CH 4 , CH 3 OH, and several C2 species with a very low FE (<4%). Even if the FE for some of the products was higher than that obtained with our system, the evolution of H 2 accounted for the largest charge consumption also in this case. Regarding the current density, similar results were obtained at the same applied potential values (5−10 mA cm −2 at −1.6 V vs SCE). However, the mixture of different liquids and gas products requires a further separation step. The catalytic activity, in terms of current density, of our Cu-GCN (14 mA/mg cat ) system is also comparable with the copper-doped cyanographene (10−20 mA/mg) reported by Pieta et al. 56 .

Cu-GCN stability
The stability of Cu-GCN under CO2RR conditions was investigated by subjecting it to a cycling treatment and then characterizing it by XRD and EXAFS. Results are reported in Fig. 7. The comparison of the XRD patterns for the fresh and cycled electrodes do not show the formation of copper or copper oxides nanoparticles after the cycling treatment under CO2RR conditions (Fig. 7a, b); however, the absence of peaks associated to copper species does not exclude the formation of clusters smaller than 2 nm (XRD detection limit). For this reason, the sample was analyzed by EXAFS. After the cycling treatment, the FT EXAFS spectrum does not differ significantly from the fresh one (Fig. 7c). The absence of an intense peak at 2.2 Å after the cycling treatment indicates that copper nanoparticles are not formed under catalytic conditions, differently to what reported in the literature for other Cu single atoms structures under CO2RR  conditions 6,14 . The structural parameters obtained from the fit, reported in Supplementary Table 1, indicate that the Cu-GCN sample is stable under CO2RR conditions, which could be attributed to the confinement of the Cu atoms in the structural cavities of GCN and to the strong copper−nitrogen interaction.

DISCUSSION
In this work, we synthesized and carefully characterized a 2D GCN functionalized with Cu single atoms by combining experimental and theoretical techniques. In this way, we demonstrated the formation of non-homogeneous samples with Cu adatoms at different adsorption sites and at different local concentrations, which explains the presence of two oxidation states. Subsequently, we tested Cu-GCN as electrocatalyst for the CO 2 reduction reaction in two different electrolytes in order to determine the effect on the products selectivity. In bicarbonate solution, formate was obtained as the only product in the liquid phase facilitating its separation; whereas, in phosphate solution, hydrogen was the only product detected. The catalytic results obtained are in line with those reported in the literature for similar systems. Even if a further improvement of the system will be necessary in order to increase the products selectivity, the results obtained suggest that GCN is a promising substrate able to stabilize metal atoms in its structural cavities thanks to the strong copper−nitrogen interaction established.

METHODS
Synthesis of carbon nitride and copper-modified carbon nitride 1 g of melamine was placed in an alumina crucible with a cover and introduced in a tubular furnace. The temperature was increased up to 100°C and maintained for 30 min in order to remove the humidity. Then, it was further increased up to 400°C and kept for 2 h in order to let the melamine condense and form heptazine units (melem). Finally, the temperature was increased up to 550°C and kept for 4 h to get the final polymerized product. The thermal process was carried out in Ar atmosphere (10 sccm) and using a temperature ramp of 2°C min −1 . The powder obtained had a yellow colour, characteristic of the bulk GCN (bulk GCN) 72 .
In order to obtain the 2D material, an exfoliation step was carried out. 1 g of bulk GCN was added into 30 mL of H 2 SO 4 (98 % w/w, Sigma-Aldrich) and kept under stirring overnight. Then, the suspension was cooled down in an ice/water bath and 100 mL of deionized water was added drop by drop. The diluted suspension was kept under sonication for 2 h and then under stirring for 24 h. The exfoliated GCN product was centrifuged at 9000 rpm for 10 min and washed with deionized water. The final product (GCN) was obtained by drying under vacuum. The powder showed the characteristic white colour of the 2D GCN 73 .
In order to modify 2D GCN with copper, 0.08 g of GCN and the necessary amount of Cu(NO 3 ) 2 ·2.5H 2 O (Sigma-Aldrich) to obtain a nominal metal loading of 10, 20, or 30 wt.% were dispersed in 100 mL of Milli-Q water and sonicated for 3 h under Ar atmosphere. Then, 20 mL of 0.5 M NaBH 4 solution were slowly added and let react for 1 h. Finally, the catalyst was separated by centrifugation, washed with Milli-Q water and ethanol, and

Characterization techniques
SEM images were acquired with a ZEISS Supra VP35 microscope, whereas the TEM images with a FEI Tecnai G12 microscope. The XRD characterization was performed with a Bruker D8 Advance, configured with a glancing angle geometry, operating with Cu K α radiation (λ = 0.15406 nm) generated at 40 kV and 40 mA. Raman spectra were collected using a ThermoFisher DXR Raman microscope using a laser with an excitation wavelength of 532 nm (1 mW), focused on the sample with a 50× objective (Olympus). Solid-state Fourier Transformed Infrared (FT-IR, KBr disk technique) absorption spectra were recorded with a Jasco FT-IR-4100spectrometer. X-ray photoelectron spectroscopy (XPS) measurements were acquired in a custom-made UHV system working at a base pressure of 10 −10 mbar, equipped with an Omicron EA125 electron analyzer and an Omicron DAR 400 X-ray source with a dual Al−Mg anode. Core-level photoemission spectra (C 1s, N 1s, O 1s, and Cu 2p) were collected at room temperature with a non-monochromatized Al K α X-ray source (1486.3 eV). Single spectra were acquired using 0.1 eV steps, 0.5 s collection time, and 20 eV pass energy. In order to analyze the single components of the C 1s and N 1s regions, the spectra were separated into chemically shifted components. For the C 1s region, an asymmetrical shape was used for the sp 2 component, whereas symmetrical Voight functions were used for the sp 3 component, the C−O functional groups, and C−(N) 3 species. For the N 1s region, symmetrical Voight functions were used. X-ray adsorption spectroscopy measurements were recorded on beamline B18 at Diamond Light Source (UK) with ring energy of 3 GeV and a current of 300 mA. The monochromator used was Si(311) crystals operating in Quick EXAFS (QEXAFS) mode. Pellets of Cu-GCN and the oxides used as references were measured in transmission mode at the Cu K (7709 eV) absorption edge at 298 K using a 36-element Ge detector. The Cu foil was measured simultaneously. In addition, a button electrode was prepared by painting a carbon paper electrode (1 × 1 cm) with the Cu-GCN sample, and subjected to an accelerated ageing treatment (1000 cycles between −0.6 and −1.5 V SCE at 50 mV s −1 in 0.5 M KHCO 3 ). The cycled electrode was also measured by ex situ XAS in fluorescence mode. Calibration of the monochromator was carried out using the Cu foil previously to the measurements. The acquired data were processed and analyzed using the Athena programme. Solid-state absorption spectra were recorded using a Cary 5000 UV−Vis-Spectrometer equipped with a Diffuse Reflectance Accessory consisting of an integrating sphere. The spectra were acquired and plotted as Kubelka−Munk function F(R) = A*C/S (A = absorbance, C = Concentration, S = scattering coefficients).

Electrochemical methods
Electrochemical measurements were conducted in a three-electrode cell using a saturated calomel electrode (SCE in KCl sat.) as a reference electrode and a graphite rod or a Pt wire as a counter electrode. All the experiments were performed at room temperature. The currents reported are normalized by the electrode geometric area and the potential values are referred to SCE. The working electrode consisted of a thin-film catalyst layer deposited on a glassy carbon (GC) electrode (3 mm diameter for the rotating-disk experiments, 1 or 2 cm 2 GC plate for the electrolysis measurements). Before the catalyst modification, the GC electrode was polished with 0.25 μm diamond paste and sonicated in isopropyl alcohol. The catalyst ink was prepared by mixing 2 mg of catalyst, 10 µl Nafion ® (5 wt.%, Sigma-Aldrich), 500 µl of Milli-Q water, and 500 µl of dimethylformamide. The ink was sonicated for 20 min in order to facilitate the dispersion of the catalyst. After ink deposition, the electrodes were dried in the air. Measurements were carried out in N 2 -and CO 2 -saturated 0.5 M KHCO 3 /0.1 M phosphate buffer, prepared from high purity reagents (Sigma-Aldrich). The RDE measurements were recorded in the absence and presence of light, using a white LED lamp (OSRAM-9 LED Star Oslon 150 white 1170 lm cm −2 @350 mA) as a light source, to investigate any activity enhancement by light. The ohmic drop compensation was performed during the data treatment, by using a value for the solution resistance measured by electrochemical impedance spectroscopy (generally 40 Ω for the KHCO 3 solution).
A RRDE with a GC disk (5 mm diameter) and the platinum ring was used for the characterization of the products obtained during the CO2RR 67,68 . The ring was cycled between −0.4 and +1.0 V vs SCE at 0.2 Vs −1 (3500 or 500 rpm as rotating speed), while the disk was held at different potential values between −0.6 and −1.7 V vs SCE.
A H-type cell was used for the chronoamperometry tests. Both the graphite counter electrode (or Pt wire) and the SCE were separated from the cathodic compartment by a ceramic frit. Two kinds of experiments were carried out: on one hand, the cell was completely sealed after 1-h purging with CO 2 and the solution was kept under magnetic stirring during the test. The gas and liquid products were analyzed by headspace gas chromatography and 1 H-NMR, respectively. Preliminary experiments showed that only H 2 and HCOO − were produced. On the other hand, the solution was continuously purged by CO 2 during the whole experiment and the liquid phase was analyzed by 1 H-NMR.

Analysis of the products
The liquid products were determined by 1 H-NMR analysis using a Bruker Avance (300 MHz) with the water pre-saturation method. Samples were prepared by mixing 400 µl of the electrolytic solution with deuterated water (D 2 O, 50 µl) and 50 µl of DMSO solution (prepared by dissolving 43.1 mg DMSO in 50 ml of Milli-Q water) as an internal standard (the DMSO peak was fixed at 2.71 ppm 74 ). The gaseous products (H 2 , CO, and CH 4 ) were analyzed by gas chromatography (Agilent Technologies 7890 A) equipped with a Molsieve 5 Å column, using a thermal conductivity detector and Ar as a carrier gas. The calibration curves for the different volatile products were obtained by injecting exact amounts of the pure gases. The measurements were carried out by sampling the headspace with a 500 μl gas-tight syringe.
The amount of HCOO − produced by a 0.26 mg cm −2 modified GC disc (1.1 cm 2 ) in a 4-h experiment under a continuous flow of CO 2 at −1.55 V vs SCE was 35 µmol. In more detail, the concentration of the product in the NMR tube was calculated by using Eq. 1: By keeping in mind that H DMSO ð Þ H HCOO À ð Þ is equal to 6, Area HCOO À peak Area DMSO peak = 0.4466, [DMSO] = 1.10 mM and the volume of the electrolyte in the CPE experiment was 9.5 ml, the amount of HCOO − produced resulted to be 35 µmol (corresponding to an average production rate of 31 µmol HCOO− mg cat −1 h −1 ). By considering that HCOO − is the result of a 2-electron reduction of CO 2 , the FE at a specific time t was calculated by using Eq. 2: FE% HCOO À ð Þ¼ 2 mol HCOO À t ð Þ F Charge ðtÞ 100 (2) with F (Faraday constant) = 96485 C mol −1 .

Computational methods
The calculations were carried out within the DFT framework, using the planewaves code QuantumESPRESSO [75][76][77] , through the Perdew−Burke−Erhenzoff exchange-correlation functional 78 , accounting for the Van der Waals dispersion forces in the form of the Grimme-D3 correction term 79 . Ultrasoft pseudopotentials 80 were employed from the most recent versions of the pslibray 81 , with wavefunction and charge density cutoffs of 46 and 326 Ry, respectively. Each model was constructed with a supercell containing several monomers, in order to allow the proper number of degrees of freedom for the structural optimization; for these structures, a converged grid of shifted Monkhorst−Pack 82 mesh of k-points was used. The N 1s XPS spectra were calculated through the ΔSCF scheme 83 , with the use of pseudopotentials modified as to contain a core-hole. In this way, given a chosen model, by running a separate calculation for each inequivalent coreexcited atom it is possible to reconstruct the shift of its binding energy (BE) from an artificial reference, called core level shift (CLS). In our setup, the reference BE was chosen as the average, different for each structure, of the total energies of the core-excited systems. The chemical shift between the XPS spectra of different models was calculated through the inclusion of an N 2 molecule far from the surfaces, acting as a common reference, which has already proven useful in the past 84 . In order to get a better comparison with experimental results, a Gaussian line shape was superimposed to each calculated CLS, with a half-width-half-maximum of 0.7 eV.
We constructed the pristine structure (2D GCN) as a 2D sheet of melon, whose simulated properties, namely XRD and XPS spectra, show a good agreement with the ones observed experimentally: the models and their results are shown in Supplementary Figs. 3 and 4. To fit experimental data, we take into account the formation of different melon structures, which we label 'alternated' or 'linear', differing by the relative arrangement of the monomers.
The Cu-GCN models were built by placing single Cu adatoms at three N sites in the monomer fragment, as shown in Fig. 3d or in Supplementary Fig. 4, which are labelled 'Amino I', 'Amino II' and 'Top', respectively. Two coverage regimes are considered, namely 'low' and 'high', as defined in Supplementary  Fig. 4 and detailed in the Supplementary Discussion.

DATA AVAILABILITY
All the data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Data related to this paper may be requested from the authors.

CODE AVAILABILITY
The code used in this work is available at www.quantum-espresso.org.