Synergistic Effects of B-F/B-S and Nitrogen Vacancy Co-Doping on g-C3N4 and Photocatalytic CO2 Reduction Mechanisms: A DFT Study

Nonmetallic co-doping and surface hole construction are simple and efficient strategies for improving the photocatalytic activity and regulating the electronic structure of g-C3N4. Here, the g-C3N4 catalysts with B-F or B-S co-doping combined with nitrogen vacancies (Nv) are designed. Compared to the pristine g-C3N4, the direction of the excited electron orbit for the B-F-co-doped system is more matching (N2pz→C2pz), facilitating the separation of electrons and holes. Simultaneously, the introduced nitrogen vacancy can further reduce the bandgap by generating impurity states, thus improving the utilization rate of visible light. The doped S atoms can also narrow the bandgap of the B-S-Nv-co-doped g-C3N4, which originates from the p-orbital hybridization between C, N, and S atoms, and the impurity states are generated by the introduction of N vacancies. The doping of B-F-Nv and B-S-Nv exhibits a better CO2 reduction activity with a reduced barrier for the rate-determining step of around 0.2 eV compared to g-C3N4. By changing F to S, the origin of the rate-determining step varies from *CO2→*COOH to *HCHO→*OCH3, which eventually leads to different products of CH3OH and CH4, respectively.


Introduction
Fossil fuels are a vital source of energy for human beings [1]; however, the use of fossil fuels releases CO 2 to accumulate in the atmosphere, leading to a serious greenhouse effect [2,3]. Finding efficient ways to reduce CO 2 emissions or capture CO 2 and convert it into chemical commodities is the key problem to solve [4][5][6][7]. Lightweight materials with different nano-morphologies containing heteroatoms have been actively pursued for catalytic applications and as semiconductor materials because of their unusual properties [8]. Graphitic carbon nitride (g-C 3 N 4 ) is a star two-dimensional metal-free semiconductor material with a bandgap of 2.7 eV [9]. Unfortunately, the catalytic reduction of CO 2 using pristine g-C 3 N 4 is very inefficient due to the fast recombination of photogenerated holes and electrons and the low utilization efficiency of visible light [10,11]. To improve the photocatalytic efficiency of pristine g-C 3 N 4 , many studies have been carried out, including structure optimization, doping modification [12][13][14][15], structural and defect engineering [16], and composite heterojunction materials [17][18][19][20][21][22][23]. Experimental studies have shown that nonmetallic doping, such as B, O, S, P, and halogen doping, is an effective way to enhance the photocatalytic performance and maintains the metal-free photocatalytic properties of g-C 3 N 4 [24][25][26][27]. Yan found that boron-doped g-C 3 N 4 effectively utilized visible light, increasing the photocatalytic degradation rate of rhodamine B by 3.6 times [8], while Ohno et al. applied it in the photocatalytic CO 2 reduction process and generated methanol [24]. In the halogen-doped system, the F-doped g-C 3 N 4 was reported to have the smallest bandgap energy and the strongest light absorption ability, and the doped g-C 3 N 4 was corrugated to

Computational Methods
The electronic structural properties of the pristine and doped g-C 3 N 4 catalysts were calculated using the Vienna Ab-initio Simulation Package (VASP) [42], and the projected plus plane wave (PAW) pseudopotential [40] was used to describe the interaction between ionic real and valence electrons. The unit cell energy and lattice constant were calculated using the PBE functional of the generalized gradient approximation method (GGA) [43]. During the structural optimization calculations, a very precise accuracy of 1.0 × 10 −6 eV was used as the criterion for energy convergence. For the K-point setting of the system, the Mokhorst-Pack setting of 3 × 3 × 1 was used in the structural relaxation, and the K-point setting of 6 × 6 × 1 was used in the static calculation and performance calculation. Considering that the PBE method tended to underestimate the bandgap value, the band structures of the pristine and doped systems were calculated using the HSE06 hybrid functional [44,45]. In order to avoid mutual interference between adjacent systems, a vacuum space of 15 Å was introduced and the energy cutoff value of the system was set to 450 eV, and the effect of spin polarization was considered in the calculation. The PBE (D3) method with the Grimme van der Waals correction was employed because of the weak interactions between CO 2 Photocatalytic Reduction (CO 2 RR) Reaction species and catalysts [46].
Pure bulk g-C 3 N 4 is built on the hexagonal heptadiazine-based structure as it is the most stable graphitic phase [47], and it contains two or three different kinds of C and N atoms, which are marked as C1-C2 and N1-N3, respectively, as shown in Figure 1a,b. According to the experimental results, B-N bonds and N-F bonds existed in the XPS spectra of B-F-co-doped g-C 3 N 4 [25], indicating that B atoms preferentially replaced C atoms (C1), and F atoms tended to combine with N atoms (N2) to be doped into the interspace, as shown in Figure 1c. There are 5 possible sites where N v can be introduced, namely, N v 1-N v 5 in B-F-co-doped g-C 3 N 4 , as shown in Figure 1c. For B-S-co-doped g-C 3 N 4 , similar possible doping sites are provided, as given in Figure 1d, based on the experiment [26]. The optimal doping site can be determined by the formation energy (E form ), according to the following equation: point setting of 6 × 6 × 1 was used in the static calculation and performance calculation. Considering that the PBE method tended to underestimate the bandgap value, the band structures of the pristine and doped systems were calculated using the HSE06 hybrid functional [44,45]. In order to avoid mutual interference between adjacent systems, a vacuum space of 15 Å was introduced and the energy cutoff value of the system was set to 450 eV, and the effect of spin polarization was considered in the calculation. The PBE (D3) method with the Grimme van der Waals correction was employed because of the weak interactions between CO2 Photocatalytic Reduction (CO2RR) Reaction species and catalysts [46]. Pure bulk g-C3N4 is built on the hexagonal heptadiazine-based structure as it is the most stable graphitic phase [47], and it contains two or three different kinds of C and N atoms, which are marked as C1-C2 and N1-N3, respectively, as shown in Figures 1a,b. According to the experimental results, B-N bonds and N-F bonds existed in the XPS spectra of B-F-co-doped g-C3N4 [25], indicating that B atoms preferentially replaced C atoms (C1), and F atoms tended to combine with N atoms (N2) to be doped into the interspace, as shown in Figure 1c. There are 5 possible sites where Nv can be introduced, namely, Nv1-Nv5 in B-F-co-doped g-C3N4, as shown in Figure 1c. For B-S-co-doped g-C3N4, similar possible doping sites are provided, as given in Figure 1d, based on the experiment [26]. The optimal doping site can be determined by the formation energy (Eform), according to the following equation: Eundoped and Edoped are the electron energies of the catalysts before and after doping, respectively. μA and μB are the chemical potentials of the substituted atoms, i.e., μC, μN, μB, μF, and μS. Formation energy can also reflect the difficulty of doping; the lower the formation energy, the better the thermodynamic stability.

Geometric Structure and Thermodynamic Stability
The formation energy results are shown in Table S1 to determine the location of nitrogen vacancies. The Eform of the doped Nv4 site of B-F-Nv is lowest among the possible doping sites, indicating that it is the most stable structure, and the negative value of −0.55 eV suggests that the introduction of Nv is thermodynamically spontaneous [41] and easier to synthesize experimentally. In contrast, the lowest Eform of B-S-Nv is 1.07 eV, E undoped and E doped are the electron energies of the catalysts before and after doping, respectively. µA and µB are the chemical potentials of the substituted atoms, i.e., µC, µN, µB, µF, and µS. Formation energy can also reflect the difficulty of doping; the lower the formation energy, the better the thermodynamic stability.

Geometric Structure and Thermodynamic Stability
The formation energy results are shown in Table S1 to determine the location of nitrogen vacancies. The E form of the doped N v 4 site of B-F-N v is lowest among the possible doping sites, indicating that it is the most stable structure, and the negative value of −0.55 eV suggests that the introduction of N v is thermodynamically spontaneous [41] and easier to synthesize experimentally. In contrast, the lowest E form of B-S-N v is 1.07 eV, corresponding to the doped site of N v 2, suggesting that it is more difficult to introduce N v on the basis of B-S than B-F-doped g-C 3 N 4 .
To explain how the doped elements and vacancies affect the properties of pristine g-C 3 N 4 , the optimized bond length parameters for pure g-C 3 N 4 , B-F, and B-F-N v -doped g-C 3 N 4 are shown in Figure 2, Table S2 and compared with the experimental [48] and theoretical values [28]. It can be seen that the optimized geometric parameters agree well with the available results; the maximum deviations from theory and experiment are 0.02 and 0.13 Å, respectively, verifying the reliability of the calculation method and parameter settings. Figure 2a-d show the variation in the key bond lengths from g-C 3 N 4 , B-F, and B-F-N v to B-S-N v . Compared to the pristine g-C 3 N 4 , the introduction of elements and holes only changes the local geometry of the doping site and has little effect on other sites for the doped g-C 3 N 4 . The dramatic change lies in B (C2)-N1 bonds, which is elongated from 1.47 and 1.51 to 1.62 Å gradually from g-C 3 N 4 to the introduction of B and F atoms and then the nitrogen vacancies. This is because the radius of the B atom (0.88) is larger than that of the C atom (0.77) and the electronegativity of the F atom (3.98) is stronger. For B-F-co-doped g-C 3 N 4 , the introduction of the F atom with the stronger electronegativity can facilitate the electrons flowing to the F atom, thus weakening the B-N bond. With the further introduction of N v , the F atom shifts from N2 to the B atom, and importantly, the C atom changes from sp 2 hybridization to sp hybridization, leading to the enhanced electronegativity of the C atom. Thus, the interaction between the C and the N atom is further enhanced, and the interaction between the N and the B atom is weakened, leading to the further elongation of the B-N bond. Similarly, for B-S-N v , the introduction of B and S elements and N v extends the distance of the B-N1 bond to 1.53 Å. Due to the absence of N, the structure changes from a six-membered ring to a five-membered ring at the nitrogen vacancy ring, and local large deformation occurs. Generally, with the degree of local deformation of the geometric structure increased, the mutual repulsion between the lone pairs of electrons of the N atom in the system will reduce and the material will become more stable [49].
C3N4, the optimized bond length parameters for pure g-C3N4, B-F, and B C3N4 are shown in Figure 2, Table S2 and compared with the experimenta retical values [28]. It can be seen that the optimized geometric parameters the available results; the maximum deviations from theory and experime 0.13 Å, respectively, verifying the reliability of the calculation method and tings. Figures 2a-d show the variation in the key bond lengths from g-C3N Nv to B-S-Nv. Compared to the pristine g-C3N4, the introduction of elem only changes the local geometry of the doping site and has little effect on the doped g-C3N4. The dramatic change lies in B (C2)-N1 bonds, which is 1.47 and 1.51 to 1.62 Å gradually from g-C3N4 to the introduction of B an then the nitrogen vacancies. This is because the radius of the B atom (0.88 that of the C atom (0.77) and the electronegativity of the F atom (3.98) is s F-co-doped g-C3N4, the introduction of the F atom with the stronger electr facilitate the electrons flowing to the F atom, thus weakening the B-N bon ther introduction of Nv, the F atom shifts from N2 to the B atom, and imp atom changes from sp 2 hybridization to sp hybridization, leading to the en negativity of the C atom. Thus, the interaction between the C and the N enhanced, and the interaction between the N and the B atom is weakened further elongation of the B-N bond. Similarly, for B-S-Nv, the introduction ments and Nv extends the distance of the B-N1 bond to 1.53 Å. Due to th the structure changes from a six-membered ring to a five-membered ring vacancy ring, and local large deformation occurs. Generally, with the deg formation of the geometric structure increased, the mutual repulsion be pairs of electrons of the N atom in the system will reduce and the mater more stable [49].

Electronic Structure and Electronic Properties
The bandgaps of g-C3N4, and B-F-, B-F-Nv-, and B-S-Nv-co-doped g-C by different methods are listed in Table S3 together with experimental and ues [8,14,25]. The bandgap values calculated by the HSE06 hybrid function are in good agreement with those of previous studies [24,50], again

Electronic Structure and Electronic Properties
The bandgaps of g-C 3 N 4 , and B-F-, B-F-N v -, and B-S-N v -co-doped g-C 3 N 4 calculated by different methods are listed in Table S3 together with experimental and literature values [8,14,25]. The bandgap values calculated by the HSE06 hybrid functional in this paper are in good agreement with those of previous studies [24,50], again testifying the rationality of the structures. Due to the limitations of the GGA method, the multi-electron interactions cannot be fully described, and the bandgap value is often underestimated, but it does not affect the results of qualitative comparisons and the rules. The results show that B-F co-doping increases the bandgap value from 2.77 to 3.06 eV compared to g-C 3 N 4 , while the injection of N v significantly reduces the bandgap value from 3.06 to 2.67 eV, which facilitates the broadening of the visible light utilization.
The band structure and contribution of atoms and atomic orbitals to the density of states for g-C 3 N 4 , and B-F-N v -and B-S-N v -co-doped g-C 3 N 4 are presented in Figure 3. It can be seen that the valence band maximum (VBM) and conduction band minimum (CBM) for g-C 3 N 4 are located at point Gamma and point K, respectively. Due to the different positions of the VBM and CBM in k-space, g-C 3 N 4 is an indirect bandgap semiconductor with a bandgap value of 2.77 eV calculated by the HSE06 method. The valence band is mainly contributed by N (2px, 2py) orbitals of N atoms, while the conduction band is mainly C 2pz of C atoms and a small amount of N 2pz orbitals of N atoms. Due to the different orbital orientations, it is difficult to excite electrons from N (2px, 2py) to C 2pz or N 2pz orbitals, which is one of the reasons for the difficult transition of photogenerated electrons and low carrier mobility in g-C 3 N 4 [51]. rationality of the structures. Due to the limitations of the GGA method, the multi-electron interactions cannot be fully described, and the bandgap value is often underestimated, but it does not affect the results of qualitative comparisons and the rules. The results show that B-F co-doping increases the bandgap value from 2.77 to 3.06 eV compared to g-C3N4, while the injection of Nv significantly reduces the bandgap value from 3.06 to 2.67 eV, which facilitates the broadening of the visible light utilization. The band structure and contribution of atoms and atomic orbitals to the density of states for g-C3N4, and B-F-Nv-and B-S-Nv-co-doped g-C3N4 are presented in Figure 3. It can be seen that the valence band maximum (VBM) and conduction band minimum (CBM) for g-C3N4 are located at point Gamma and point K, respectively. Due to the different positions of the VBM and CBM in k-space, g-C3N4 is an indirect bandgap semiconductor with a bandgap value of 2.77 eV calculated by the HSE06 method. The valence band is mainly contributed by N(2px, 2py) orbitals of N atoms, while the conduction band is mainly C2pz of C atoms and a small amount of N2pz orbitals of N atoms. Due to the different orbital orientations, it is difficult to excite electrons from N(2px, 2py) to C2pz or N2pz orbitals, which is one of the reasons for the difficult transition of photogenerated electrons and low carrier mobility in g-C3N4 [51]. In contrast, B-F-Nv is a direct bandgap structure with a value of 2.67 eV calculated by the HSE06 method, smaller than g-C3N4, as both VBM and CBM are located at point Gamma. The direct transition of electrons reduces the energy loss caused by electron relaxation, which means a higher energy conversion efficiency [52][53][54]. The loss of nitrogen In contrast, B-F-N v is a direct bandgap structure with a value of 2.67 eV calculated by the HSE06 method, smaller than g-C 3 N 4 , as both VBM and CBM are located at point Gamma. The direct transition of electrons reduces the energy loss caused by electron relaxation, which means a higher energy conversion efficiency [52][53][54]. The loss of nitrogen atoms caused by vacancies destroys the π-conjugation of g-C 3 N 4 , leading to the reconstruction of the band structure and the formation of impurity states in the middle of the forbidden band. Due to the low doping concentration of B (1.73 at.%) and F (1.75 at.%) atoms, the contribution of both types of atoms to the band edges is not obvious [49]. The introduction of nitrogen vacancies generates additional electron and spin polarization in g-C 3 N 4 , and an impurity state is created in the middle of the forbidden band of the spin states. The valence band consists mainly of N 2pz and N (2px, 2py) orbitals and the conductive band is mainly composed of C 2pz orbitals. As the electron transition orbital direction is the same, it is easier to realize the electron transition and improve the carrier mobility. The impurity state consists of C 2pz and N 2pz orbitals, which belong to the hanging bonds of carbon atoms in nitrogen vacancies [27]. The composition of the impurity state is consistent with the composition at the bottom of the conduction band, which can receive the electrons excited at the top of the valence band and shift the light absorption wavelength to the long-wave direction. This facilitates continuous photoexcitation of the VBM to CBM in the visible light range.
Compared to g-C 3 N 4 and B-F-N v , B-S-N v is also a direct bandgap structure and the doped S atoms can further narrow the bandgap, as shown in Table S3. The PBE method is used to calculate the bandgap values for g-C 3 N 4 , B-F-N v , and B-S-N v as 1.21, 1.19, and 1.16 eV, respectively. By contrast, the impurity level is generated at near the VBM due to the introduction of nitrogen vacancies, which is mainly contributed by the N 2pz orbital. As shown in Figure 3c, the VBM is contributed by N atoms and a small amount of C and S atoms. The overlap of electronic states of S and C indicates the hybridization of S and C atomic orbitals, while the CBM is mainly contributed by C atoms and a small amount of N atoms. The B atoms have little contribution to the energy band edges and do not directly participate in the generation, separation, and migration of photogenerated carriers, but act as electron donors to enhance the photocatalytic performance of g-C 3 N 4 by affecting the charge distribution of surrounding atoms [8]. The formed impurity state at the top of the valence band of the spin-down state consists of hybridized C 2pz , N 2pz , and S 3pz orbitals. This is because the doped S atoms affect the distribution of C and N atoms in the lattice, resulting in p-orbital hybridization between C, N, and S atoms. Figure 4 presents the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO) of g-C 3 N 4 , B-F-N v , and B-S-N v . For g-C 3 N 4 , the HOMO covers all the corner N atoms, while the LUMO is mainly distributed on the C and N atoms. The bridge N atoms are not involved in either the valence or the conduction band, so they do not excite or accept electrons, hindering the transfer of charge carriers between the heptylhydrazine rings by bridging N atoms. Both photogenerated electrons and holes are concentrated on each heptylhydrazine ring, resulting in a high carrier recombination rate and poor photocatalytic activity due to the spatial overlap of HOMO and LUMO orbitals. In contrast, the charge density rearranges and electron-rich regions appear for both B-F-N v and B-S-N v . The HOMO orbitals are mainly distributed on the N v -introduced building blocks, while the LUMO orbitals are mainly on the remaining unmodified heptyltriazine rings, which differs from the spatial overlap of the HOMO and LUMO orbitals of g-C 3 N 4 , as shown in Figure 4b,c. The HOMO and LUMO orbitals of B-F-N v and B-S-N v are completely spatially separated, which will effectively improve the separation efficiency of photogenerated electrons and holes. In addition, the HOMO and LUMO orbitals are also distributed on the bridged N atoms, which facilitates the carrier migration between the heptatriazine ring building blocks and increases the carrier mobility [55,56].

Adsorption Properties of CO2
The effects of modification methods on the electronic and optical properties of g-C3N have been focused on extensively; however, the theoretical studies of their applications on specific photocatalytic reactions are still limited. In this work, the effects of B-F-Nv and B-S-Nv modification on the photocatalytic CO2 reduction of g-C3N4 are investigated. The first step of the photocatalytic reaction is the adsorption of CO2. Figure S1 presents the possible adsorption sites for CO2 on B-F-Nv and B-S-Nv, together with adsorption energy to determine the optimal adsorption site. The adsorption energy of CO2 placed above the nitrogen vacancy in parallel (Figures S1a,b) and vertically ( Figure S1c) is −0.16 eV, while that placed above the interstitial in parallel ( Figure S1d) is −0.26 eV, indicating that the optimal adsorption site is above the interstitial. This may be because after modification and under the visible light excitation, holes appear around the nitrogen vacancies to gen erate oxidation sites and exhibit oxidation activity, while the charge density increases from −1.12 e to −1.39 e on the corner N atoms around the interstitial, as shown in Figure  S3 (red cycle). Therefore, the optimal adsorption site of CO2 on B-F-Nv is not the nitrogen vacancy with strong oxidation activity, but the interstitial with strong reduction activity.
For B-S-Nv, the adsorption energy is −1.99, −2.68, and −1.97 eV for CO2 adsorbed above the nitrogen vacancies, above the interstitial, and on top of B and S atoms, respec tively, as shown in Figures S1e-g. Thus, the optimal adsorption site of CO2 on B-S-Nv is also above the interstitial, which is similar to B-F-Nv, but the adsorption of CO2 on B-S-N is stronger than on B-F-Nv. This can be explained by the charge transfer determined by the differential charge density (CDD) map. Figure 5 is a CDD plot of CO2 adsorption on g-C3N4, B-F-Nv, and B-S-Nv. The yellow and blue color represents electron enrichment and depletion, respectively, so charge flows from the blue area to the yellow area. In the g C3N4 adsorption system, there are less charge distributions on the interstitial edge corne N atoms and bridging N atoms. In the B-F-Nv adsorption system, the charge transfers from the B atom to the adjacent N atom active site, and the charge around the nitrogen vacancy flows to the adsorption active site, forming a charge-rich environment at the edge N atom site. F atoms serve as bridges to transfer electrons to CO2. However, the charge density on the bridging F atom is too high, and the electron repulsion between the bridging F and the O atom of CO2 becomes larger, weakening the adsorption effect. In contrast, for the B-S Nv adsorption system, the charge concentration flows to the adsorption active site, and the charge density enrichment degree on the corner N atoms at the gap edge increases thereby enhancing the adsorption of CO2.

Adsorption Properties of CO 2
The effects of modification methods on the electronic and optical properties of g-C 3 N 4 have been focused on extensively; however, the theoretical studies of their applications on specific photocatalytic reactions are still limited. In this work, the effects of B-F-N v and B-S-N v modification on the photocatalytic CO 2 reduction of g-C 3 N 4 are investigated. The first step of the photocatalytic reaction is the adsorption of CO 2 . Figure S1 presents the possible adsorption sites for CO 2 on B-F-N v and B-S-N v , together with adsorption energy, to determine the optimal adsorption site. The adsorption energy of CO 2 placed above the nitrogen vacancy in parallel ( Figure S1a,b) and vertically ( Figure S1c) is −0.16 eV, while that placed above the interstitial in parallel ( Figure S1d) is −0.26 eV, indicating that the optimal adsorption site is above the interstitial. This may be because after modification and under the visible light excitation, holes appear around the nitrogen vacancies to generate oxidation sites and exhibit oxidation activity, while the charge density increases from −1.12 e to −1.39 e on the corner N atoms around the interstitial, as shown in Figure S3 (red cycle). Therefore, the optimal adsorption site of CO 2 on B-F-N v is not the nitrogen vacancy with strong oxidation activity, but the interstitial with strong reduction activity.
For B-S-N v , the adsorption energy is −1.99, −2.68, and −1.97 eV for CO 2 adsorbed above the nitrogen vacancies, above the interstitial, and on top of B and S atoms, respectively, as shown in Figure S1e-g. Thus, the optimal adsorption site of CO 2 on B-S-N v is also above the interstitial, which is similar to B-F-N v , but the adsorption of CO 2 on B-S-N v is stronger than on B-F-N v . This can be explained by the charge transfer determined by the differential charge density (CDD) map. Figure 5 is a CDD plot of CO 2 adsorption on g-C 3 N 4 , B-F-N v , and B-S-N v . The yellow and blue color represents electron enrichment and depletion, respectively, so charge flows from the blue area to the yellow area. In the g-C 3 N 4 adsorption system, there are less charge distributions on the interstitial edge corner N atoms and bridging N atoms. In the B-F-N v adsorption system, the charge transfers from the B atom to the adjacent N atom active site, and the charge around the nitrogen vacancy flows to the adsorption active site, forming a charge-rich environment at the edge N atom site. F atoms serve as bridges to transfer electrons to CO 2 . However, the charge density on the bridging F atom is too high, and the electron repulsion between the bridging F and the O atom of CO 2 becomes larger, weakening the adsorption effect. In contrast, for the B-S-N v adsorption system, the charge concentration flows to the adsorption active site, and the charge density enrichment degree on the corner N atoms at the gap edge increases, thereby enhancing the adsorption of CO 2 .

Catalytic Performance for the Photocatalytic CO2 Reduction
Carbon dioxide is a very stable linear molecule, and the photocatalysts with a suitable conduction band potential to convert carbon dioxide into hydrocarbon fuel are highly desirable. The selectivity of g-C3N4 to the reduction product is still controversial and the yield is low. In addition, how the doping modification affects the catalytic performance and product selectivity of g-C3N4 should also be investigated. The Gibbs free energy curve and adsorption energies of products in different stages for g-C3N4, B-F-Nv, and B-S-Nv are presented in Figure 6, as well as the representative geometries of the stable points in the specific reduction process of CO2 on g-C3N4. The photocatalytic CO2 reduction pathways are similar for all three catalysts, i.e., *CO2→*COOH→*CO→*CHO→*HCHO→*OCH3→ *CH3OH→*CH3→*CH4 as proposed by Azofra et al. [49]. The speculated probable products are CO, HCHO, CH3OH, and CH4, which are 2, 4, 6, and 8 electron processes, respectively. Further, Yu et al. [57]. Experimentally synthesized a g-C3N4 catalyst and applied it to the reduction of CO2 under UV-Vis radiation, confirming that the main product of CO2 reduction was CH3OH and that it reached a yield of 0.81 μmol·g −1 . In the Gibbs free energy diagram calculated in this paper, it can be seen that △G = −0.48 eV when CO2 is reduced to CH3OH, and 0.53 eV of energy needs to be absorbed in the subsequent reaction. Therefore, it can be seen that methanol exists stably, and methanol is the final reduction product, which is consistent with the experiment [57].

Catalytic Performance for the Photocatalytic CO 2 Reduction
Carbon dioxide is a very stable linear molecule, and the photocatalysts with a suitable conduction band potential to convert carbon dioxide into hydrocarbon fuel are highly desirable. The selectivity of g-C 3 N 4 to the reduction product is still controversial and the yield is low. In addition, how the doping modification affects the catalytic performance and product selectivity of g-C 3 N 4 should also be investigated. The Gibbs free energy curve and adsorption energies of products in different stages for g-C 3 N 4 , B-F-N v , and B-S-N v are presented in Figure 6, as well as the representative geometries of the stable points in the specific reduction process of CO 2 on g-C 3 N 4 . The photocatalytic CO 2 reduction pathways are similar for all three catalysts, i.e., *CO 2 →*COOH→*CO→*CHO→*HCHO→*OCH 3 →*CH 3 OH→*CH 3 →*CH 4 as proposed by Azofra et al. [49]. The speculated probable products are CO, HCHO, CH 3 OH, and CH 4 , which are 2, 4, 6, and 8 electron processes, respectively. Further, Yu et al. [57]. Experimentally synthesized a g-C 3 N 4 catalyst and applied it to the reduction of CO 2 under UV-Vis radiation, confirming that the main product of CO 2 reduction was CH 3 OH and that it reached a yield of 0.81 µmol·g −1 . In the Gibbs free energy diagram calculated in this paper, it can be seen that G = −0.48 eV when CO 2 is reduced to CH 3 OH, and 0.53 eV of energy needs to be absorbed in the subsequent reaction. Therefore, it can be seen that methanol exists stably, and methanol is the final reduction product, which is consistent with the experiment [57].

Catalytic Performance for the Photocatalytic CO2 Reduction
Carbon dioxide is a very stable linear molecule, and the photocatalysts with a suitable conduction band potential to convert carbon dioxide into hydrocarbon fuel are highly desirable. The selectivity of g-C3N4 to the reduction product is still controversial and the yield is low. In addition, how the doping modification affects the catalytic performance and product selectivity of g-C3N4 should also be investigated. The Gibbs free energy curve and adsorption energies of products in different stages for g-C3N4, B-F-Nv, and B-S-Nv are presented in Figure 6, as well as the representative geometries of the stable points in the specific reduction process of CO2 on g-C3N4. The photocatalytic CO2 reduction pathways are similar for all three catalysts, i.e., *CO2→*COOH→*CO→*CHO→*HCHO→*OCH3→ *CH3OH→*CH3→*CH4 as proposed by Azofra et al. [49]. The speculated probable products are CO, HCHO, CH3OH, and CH4, which are 2, 4, 6, and 8 electron processes, respectively. Further, Yu et al. [57]. Experimentally synthesized a g-C3N4 catalyst and applied it to the reduction of CO2 under UV-Vis radiation, confirming that the main product of CO2 reduction was CH3OH and that it reached a yield of 0.81 μmol·g −1 . In the Gibbs free energy diagram calculated in this paper, it can be seen that △G = −0.48 eV when CO2 is reduced to CH3OH, and 0.53 eV of energy needs to be absorbed in the subsequent reaction. Therefore, it can be seen that methanol exists stably, and methanol is the final reduction product, which is consistent with the experiment [57].  For B-F-N v , the formation of free radicals is endothermic because breaking carbonoxygen single or double bonds requires external input energy, i.e., *CO 2 →*COOH 1 , *CO→*CHO 3 , and *CH 3 OH→*CH 3 7 , which is similar to g-C 3 N 4 . However, the *HCHO→*OCH 3 5 process for B-F-N v releases 0.46 eV of energy, which is contrary to g-C 3 N 4 . As shown in Figure S3 (blue cycle), the charge on the adsorption site of the C atom for B-F-N v is +1.60, which is larger than those of g-C 3 N 4 and B-S-N v with the values of +1.45 and +1.49, respectively. The increased charge density facilitates the formation of OCH 3 to become an exothermic process, and thus is conducive to the subsequent reactions. The first hydrogenation process for B-F-N v is the decisive step of the whole reduction reaction with the barrier of 1.05 eV, which is lower than that of g-C 3 N 4 (1.28 eV), as shown in Figure 6a, facilitating the formation of COOH free radicals. In addition, the reaction releases energy (∆G < 0) when a neutral product is formed, such as *·COOH→*CO 2 , *CHO→*HCHO 4 , *·OCH 3 →*CH 3 OH 6 , *·CH 3 →*CH 4 8 , and there is no chemical bond with the substrate in product of stage. However, the formed radicals such as ·COOH, ·CHO, ·OCH 3 , and ·CH 3 can be bonded with the catalyst, and the adsorption sites are obviously different. For the species of C-terminal free radicals, they are easier to bond with the negatively charged N atom, and the species of O-terminated radicals are more likely to form chemical bonds with C atoms in the substrate, as shown in Figure 6c, which affects the elementary process of the reduction reaction to some extent. The energy of 0.55 eV is released when ·OCH 3 transforms to CH 3 OH 6 , and the adsorption energy is also favorable for desorption of methanol with the value of 0.82 eV, as shown in Figure 6b, so the final product is also CH 3 OH, which is similar to g-C 3 N 4 and consistent with the previous experimental results [57].
Unlike B-F-N v , in the process of CO 2 reduction on B-S-N v , there are two consecutive energy absorption processes 1 and 2 . The reaction barrier of the first hydrogenation with the value of 0.22 eV is much lower than those of both g-C 3 N 4 and B-F-N v , which is favorable for the initial reaction. The reduction of HCHO to CH 3 O 5 is the rate-determining step by the absorption energy of 1.03 eV, lower than that of g-C 3 N 4 . Once crossing the reaction rate-determining step, it undergoes three exothermic processes until reduction to the final product CH 4 with the small adsorption energy.

Conclusions
In this paper, the photocatalytic performance and photocatalytic CO 2 reduction of the nonmetallic co-doping combined with nitrogen vacancy systems, B-F-N v and B-S-N v , are investigated and compared to g-C 3 N 4 to explore the doping effects. It is found that B-F-N v and B-S-N v can combine the advantages of the three modification methods to exert a synergistic effect. B-doping can compensate for the defect of the conduction band drop caused by F-doping and N v , and maintain the high reduction potential of CO 2 . With the doped F, B-F-N v has a corrugated configuration, which reduces the repulsion between the lone pair electrons of the N atom, and improves the stability of the g-C 3 N 4 structure. At the same time, F can act as a bridge for the transfer of charges from the doping system to CO 2 . N v and B atoms form a charge-rich environment at the reactive site, which makes up for the charge loss caused by F doping and S doping.
The B-F-N v -and B-S-N v -co-doped g-C 3 N 4 have a direct bandgap structure compared to the indirect one for g-C 3 N 4 . The introduction of heteroatoms and nitrogen vacancies affects the charge distribution and improves the separation of the HOMO and LUMO, thereby enhancing the separation efficiency of photogenerated electrons and holes. The synergistic effect of B-S-N v makes the charge flow to the adsorption active site, which increases the charge enrichment of corner N atoms at the gap edge in g-C 3 N 4 , resulting in enhanced CO 2 adsorption. B-F-N v -and B-S-N v -co-doped g-C 3 N 4 have a lower rate-determining step than g-C 3 N 4 for the rate-determining step photocatalytic CO 2 reduction. B-S-N v has the strongest adsorption capacity for CO 2 , and the desorption capacity of the reduction products is stronger than that of g-C 3 N 4 , but weaker than that of B-F-N v . The suggested main reduction product is CH 3 OH for CO 2 reduction catalyzed by g-C 3 N 4 and B-F-N v , and is CH 4 catalyzed by B-S-N v . Therefore, the nonmetallic co-doping combined with N v in this work provides a new idea for improving the catalytic performance of photocatalysts.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/molecules27217611/s1, Figure S1: The possible adsorption sites for CO 2 on B-F-N v (a-d) and B-S-N v (e-g), and corresponding adsorption energy (eV); Figure S2: The reaction pathways for B-F-N v (a) and B-S-N v (b); Figure S3: The bader charge (unit e) of g-C 3 N 4 (left), B-F-N v (middle), and B-S-N v (right). Negative and positive values mean gaining and losing electrons. The elements are colored in pink for B, gray for C, blue for N, red for O, cyan for F, and yellow for S; Figure S4: Calculated band structure (left), density of states (middle), and partial density of states (right) diagrams of g-C 3 N 4 (a), B-F-N v (b) and B-S-N v (c), the energy range is from −10 eV to 10 eV; Table S1: Formation energy (eV) of B-F-N v and B-S-N v co-doped g-C 3 N 4 ; Table S2: The bond length parameters of intrinsic and doped g-C 3 N 4 [28,48]; Table S3: The band gap energy (eV) for the intrinsic and doped g-C 3 N 4 [28,31,50].