Modulation of Photocatalytic CO2 Reduction by n–p Codoping Engineering of Single-Atom Catalysts

Transition metal (TM) single-atom catalysts (SACs) have been widely applied in photocatalytic CO2 reduction. In this work, n–p codoping engineering is introduced to account for the modulation of photocatalytic CO2 reduction on a two-dimensional (2D) bismuth-oxyhalide-based cathode by using first-principles calculation. n–p codoping is established via the Coulomb interactions between the negatively charged TM SACs and the positively charged Cl vacancy (VCl) in the dopant–defect pairs. Based on the formation energy of charged defects, neutral dopant–defect pairs for the Fe, Co, and Ni SACs (PTM0) and the −1e charge state of the Cu SAC-based pair (PCu−1) are stable. The electrostatic attraction of the n–p codoping strengthens the stability and solubility of TM SACs by neutralizing the oppositely charged VCl defect and TM dopant. The n–p codoping stabilizes the electron accumulation around the TM SACs. Accumulated electrons modify the d-orbital alignment and shift the d-band center toward the Fermi level, enhancing the reducing capacity of TM SACs based on the d-band theory. Besides the electrostatic attraction of the n–p codoping, the PCu−1 also accumulates additional electrons surrounding Cu SACs and forms a half-occupied dx2−y2 state, which further upshifts the d-band center and improves photocatalytic CO2 reduction. The metastability of Cl multivacancies limits the concentration of the n–p pairs with Cl multivacancies (PTM@nCl (n > 1)). Positively charged centers around the PTM@nCl (n > 1) hinders the CO2 reduction by shielding the charge transfer to the CO2 molecule.


Introduction
Artificial photocatalytic CO 2 reduction is an intriguing research area aiming to reduce fossil consumption and mitigate the greenhouse effect [1,2].Photocatalysts are the crucial materials for readily converting CO 2 into environmentally friendly fuels via photolysis using abundant solar energy [3][4][5].With maximal atom-utilization efficiency, transition metal (TM) single-atom photocatalysts (SACs) exhibit excellent catalytic performance comparable to precious metals [6][7][8][9][10][11]. Due to its sufficient optical absorption and adequate activation centers, two-dimensional (2D) bismuth oxyhalide BiOX (X = F, Cl, Br, I) is an outstanding host to anchor the TM SACs photocatalyst and accelerate the CO 2 reduction [12][13][14][15][16].For TM-doped bismuth oxyhalide, the dopant of TM SACs can form impurity levels in the forbidden band, promoting the generation and separation of photogenerated carriers [17,18].However, the most outstanding advantage is that doped TM SACs regulate the surface state, which enhances the CO 2 adsorption by promoting the charge transfer between the activation centers and absorbed CO 2 molecules [17 -19].For example, isolated Cu SACs on BiOBr (Cu@BiOBr) establish a strong built-in electric field which serves as an electron trap to facilitate charge transfer and stabilize charge carriers.As a result, 0.5% Cu@BiOBr has a higher CO 2 absorption uptake (2.7 cm 3 g −1 ) than BiOBr (2.3 cm 3 g −1 ) [20].
In Co-SAC-doped Bi 3 O 4 Br, the isolation of Co 2+ by replacing Bi 3+ enables Co@Bi 3 O 4 Br layers to be more negatively charged, improving the CO 2 adsorption and stabilizing the *COOH intermediate on the surface [17].Charge localization induced by Fe SAC doping strengthens CO 2 bonding strengths and improves the CO 2 absorptive capacity in both porous Bi 5 O 7 I Micro-flower and B i4 O 5 I 2 [21,22].
Some issues remain to be addressed in understanding the CO 2 photolysis by the TM SACs on BiOX.The defect configurations around the TM SACs are ambiguous, especially whether halogen atoms are present or not near the TM impurity sites [17][18][19][20][21].Both experimental and theoretical works confirm that the [Bi 2 O 2 ] 2+ layer is sandwiched by two X − layers in BiOX [23][24][25][26][27][28][29].However, the unlocked bismuth (Bi) surface without a halogen covering is usually used to simulate photocatalytic CO 2 reduction around the TM SACs [17][18][19][20].This divergence highlights the necessity of assessing the effects of halogen vacancy (V halogen ) on the microstructures and reduction performance around the TM SACs.V halogen usually acts as an n-type defect, whereas TM SACs exhibit p-type characteristics and form negatively charged centers to promote the reduction reaction [23,24].
The n-p codoping concept is established based on the electrostatic attraction between the nand p-type dopants (or defects) with opposite charge states, which affects the ionization, solubility, and charge transfer of the doped semiconductors [30][31][32][33].For example, the electrostatic attraction within the n-p pair enhances both thermodynamic and kinetic solubilities, creating tunable intermediate bands to effectively narrow the band gap and enhance the visible-light photoactivity of TiO 2 [32].The n-p pairs limit the applications of Ga 2 O 3 by affecting the dopant ionization [33].On BiOX, n-p codoping can be established via the electrostatic interactions in the combination of p-type TM dopants and n-type V halogen defects.The effects of n-p codoping on the stability and reducing capacity of TM SACs must be determined.For TM SACs, the d-band theory is universally recognized as a means of evaluating the reducing capacity [34,35].Consequently, the dependence of the orbital alignment and d-band center of TM_d on n-p codoping deserves systematic study regarding the CO 2 reduction on TM SACs on BiOX.
In this paper, we investigate the effects of n-p codoping on the stability and reducing capacity of TM SACs for photocatalytic CO 2 reduction reaction (CO 2 RR) on BiOCl by using first-principles calculation.Previous experimental works found that Fe, Co, Ni, and Cu are all effective SACs for accelerating CO 2 RR on BiOX [17][18][19][20].Consequently, we constructed a dopant-defect combination by pairing these TM SACs and their surrounding Cl vacancies.Thermodynamically, the balance between nand p-defects is determined by the defect equilibrium, which can be solved by calculating the formation energies of charged defects (∆H f ) [36][37][38].Based on the ∆H f , we find that the Cl monovacancy (V Cl ) is more stable than the Cl multivacancies.n-p codoping is established by the coulomb attraction between the negatively charged p-type TM SACs and positively charged n-type V Cl .The electrostatic attraction of n-p codoping enhances the stability of TM SACs.Neutral pairs are formed for the Fe, Co, and Ni SACs (P TM 0 ), while the Cu SAC-based pair is dominated by the q = −1 charge state (P Cu −1 ).In the n-p pairs, electrostatic interaction settles the electron accumulation around TM SACs.Accumulated electron occupies localized TM_d orbital and upshifts the d-band center toward the Fermi level, enhancing the reducing capacity of TM SACs based on the d-band theory.As a result, the CO 2 absorption is improved along with the enhancement of charge transfer and the decrease in Gibbs free energies.The P Cu −1 also locates additional electrons surrounding the Cu SACs and forms a half-occupied d x 2 −y 2 state, further unshifting the d-band center and improving the reducing capacity of TM SACs.The metastability of Cl multivacancies limits the concentration of the n-p pairs with Cl multivacancies (P TM@nCl (n > 1)).Positively charged centers around P TM@nCl (n > 1) shield the charge transfer between the P TM@nCl (n > 1) and the CO 2 molecule, hindering the CO 2 reduction.

Computational Details
Density functional theory (DFT) [39] calculations were performed in the Vienna Ab initio Simulation Package (VASP) [40].The interaction of ions with electrons was determined by projector augmented wave potentials (PAWs) [41].The electron exchange-correlation energy was treated based on local-density approximation (LDA) [42] because the lattice constants and band structure of BiOCl calculated by LDA are in good agreement with experimental values [43].The results were also checked based on the Perdew-Burke-Ernzerhof (PBE) functional [44].The kinetic energy cutoff energy was 500 eV, and the total energy convergence criterion was set to 10 −5 eV.A 4 × 4 × 1 supercell with over 15 Å of lattice constant was chosen to avoid the interactions of adjacent point defects.The K-mesh was optimized based on the minimum energy principle, and 3 × 3 × 1 Monkhorst-Pack grids were found to be sufficient for sampling the Brillouin zone for the supercell.A vacuum space of over 15 Å was applied to avoid the interactions of neighbor images along the z direction.
To solve the defect equilibrium between TM dopants and the compensating V Cl , the formation energy of charged defects was calculated with the following formula [45,46]: where E tot (q,α) and E tot (host) are the total energy of the defect system and pristine BiOCl, respectively; n is the number of Cr, Bi, and TM dopants; q is the number of electrons transferred from supercell to reservoirs in forming the defect; and µ i is the chemical formula of constituent i with energy E i .E f is the Fermi energy with respect to the valence band maximum (VBM) and ranges from the VBM to the conduction band minimum (CBM).
The chemical potentials of Bi (µ Bi ) and TM (µ TM ) were derived from the energies of their corresponding metals.The chemical potentials of O and Cl were limited to avoid the formation of elementary substances and maintain a stable BiOCl compound.µ Cl and µ O are defined as: where ∆H BiOCl is the formation energy of BiOCl.Under Cl-poor conditions, the chemical potential of the Cl atom is defined as the energy in the species Cl + (−5.36 eV).The total energies of the charged systems should be corrected for the interaction of the charged defect with the compensating background and its periodic images.We used Makov-Payne (M-P) corrections, formulated as q2 α/2EL, where L is the linear dimension of the supercell, Eis the static dielectric constant, and α is the Madelung constant.
To evaluate the reaction coordinate of CO 2 reduction, Nørskov's method [47] was used to calculate the Gibbs free energy difference: ∆E, ∆E ZPE , and ∆S denote the absorbed energy difference, the zero-point energy difference, and the entropy difference between the adsorbed state and the corresponding free state for density functional theory calculations, respectively.T is the temperature of the system, 298.15 K.
The reactions of CO 2 adsorption and the *CO 2 plus H formation *COOH process are defined as follows [48]: * + CO 2 → * CO 2 (5) where * represents adsorption sites, and *CO 2 and *COOH represent adsorption intermediate states.Therefore, the CO 2 RR is calculated as: Charge density redistribution is determined by the charge density difference between the pairs and corresponding isolated systems [49].The charge density difference is calculated as follows: where ρ A and ρ B are the charge densities of isolated systems A and B, respectively, and ρ AB is the charge density of the pairs.

Structures and Stability
In stable BiOCl, the [Bi 2 O 2 ] 2+ layer is sandwiched by two Cl − layers (Figure 1), and the lattice constant is calculated to be 3.854 Å, agreeing with previous works [13,50,51].The TM SACs are doped on the BiOX by substituting a bismuth atom based on experiments [16][17][18][19][20].The TM SAC dopants are denoted as the Bi TM on the intrinsic BiOX surface as shown in Figure 1b.The dopant-defect pair (P TM@nCl ) comprises the Bi TM and its surrounding Cl vacancy as shown in Figure 1.The P TM@3Cl is the combination of a TM SAC and Cl trivacancy.In the Cl trivacancy, the proportion of unlocked Bi atoms is up to 19%, which is large enough to comprise the unlocked Bi region.The structures of the Bi TM and P TM@1Cl (P TM for short) are given in Figure 1b,c, and the P TM@nCl (n > 1) is shown in Figure S1.To assess the stability of the pairs, the formation energy (∆H f ) of charged defects is investigated as a function of E f based on Equation (1).The cathode where the CO 2 RR takes place is an electron reservoir and corresponds to the n-type semiconductor.Therefore, we focus on the electron-rich (e-rich) condition where the E f is close to the CBM.The ∆H f calculated following Equation (1) also depends on µ Cl .We focus on the Cl-rich limit (µ Cl = 0 eV), considering that the electrolyte solution usually accelerates the surface Cl ion exchange in photocatalytic experiments.

G
Charge density redistribution is determined by the charge density difference between the pairs and corresponding isolated systems [49].The charge density difference is calculated as follows: where ρA and ρB are the charge densities of isolated systems A and B, respectively, and ρAB is the charge density of the pairs.

Structures and Stability
In stable BiOCl, the [Bi2O2] 2+ layer is sandwiched by two Cl − layers (Figure 1), and the lattice constant is calculated to be 3.854 Å, agreeing with previous works [13,50,51].The TM SACs are doped on the BiOX by substituting a bismuth atom based on experiments [16][17][18][19][20].The TM SAC dopants are denoted as the BiTM on the intrinsic BiOX surface as shown in Figure 1b.The dopant-defect pair (PTM@nCl) comprises the BiTM and its surrounding Cl vacancy as shown in Figure 1.The PTM@3Cl is the combination of a TM SAC and Cl trivacancy.In the Cl trivacancy, the proportion of unlocked Bi atoms is up to 19%, which is large enough to comprise the unlocked Bi region.The structures of the BiTM and PTM@1Cl (PTM for short) are given in Figure 1b,c, and the PTM@nCl (n > 1) is shown in Figure S1.To assess the stability of the pairs, the formation energy (ΔHf) of charged defects is investigated as a function of Ef based on Equation ( 1).The cathode where the CO2RR takes place is an electron reservoir and corresponds to the n-type semiconductor.Therefore, we focus on the electron-rich (e-rich) condition where the Ef is close to the CBM.The ΔHf calculated following Equation (1) also depends on µCl.We focus on the Cl-rich limit (µCl = 0 eV), considering that the electrolyte solution usually accelerates the surface Cl ion exchange in photocatalytic experiments.Based on previous experiments, surface halogen atoms significantly influence the charge transfer between activation sites and absorbed molecules [52,53].The stability of Cl vacancies is evaluated prior to the dopant-defect pairs.Figure S2a illustrates the formation energy of the Cl monovacancy VCl, Cl di-vacancy (V2Cl), and tri-vacancy VCl (V3Cl) with varied charge states as a function of Ef.The Ef ranges from the VBM and CBM.One Based on previous experiments, surface halogen atoms significantly influence the charge transfer between activation sites and absorbed molecules [52,53].The stability of Cl vacancies is evaluated prior to the dopant-defect pairs.Figure S2a illustrates the formation energy of the Cl monovacancy V Cl , Cl di-vacancy (V 2Cl ), and tri-vacancy V Cl (V 3Cl ) with varied charge states as a function of E f .The E f ranges from the VBM and CBM.One can see that the V Cl is energy favorable in these three kinds of Cl vacancies if E f is above the mid-gap (i.e., for the e-rich condition).In the V 3Cl , 3/16 of Bi atoms are uncovered by the surface Cl atoms so that the V 3Cl is large enough with regard to the unlocked Bi region.Consequently, the V Cl is more stable than the unlocked Bi region in the BiOCl cathode.The V Cl is capable of forming a donor defect with q = +1 charge state (V Cl +1 , with the superscript referring to the charge state) when the E f is in most of the gap and becomes a neutral defect only when the E f is closer to the CBM.The charge transition level E(+1/−1) is only 0.07 eV above the CBM, indicating that it is a donor defect.The Fermi level of the V Cl 0 is located in the conduction band in the spin-polarized density of states, which also indicates the donor defect characteristic.
Considering the stability of the V Cl , we first investigate the P TM , which pairs the Bi TM and surrounding single V Cl , as shown in Figure 1c.The formation energy ∆H f is calculated to investigate the stability and ionization of the P TM and Bi TM .If E f is close to the CBM, as shown in Figure 2, the Bi TM −1 dominates the Bi TM for all the Fe, Co, Ni, and Cu SACs, indicating the p-type impurities characteristic of these TM SACs on pristine BiOCl in the e-rich condition.When the positively charged V Cl +1 and negatively charged Bi TM −1 move close to each other, as shown in Figure 1c, n-p codoping is established due to the coulomb attraction and forms the P TM pairs.Due to the electrostatic attraction of the n-p codoping, the ∆H f of the P TM is lower than the relative Bi TM as shown in Figure 2. Consequently, the V Cl facilitates the stability of the TM SACs on the BiOCl.When E f is close to the CBM, neutral pairs (P TM 0 ) are energetically favorable for the Fe, Co, and Ni SACs, while the negatively charged P Cu −1 with charge state q = −1 is most stable for the Cu SACs.The P Cu −1 possesses the lowest E f (as low as −1.56 eV) in all the defects, suggesting that the P Cu −1 can exist in great quantities on the BiOX cathode.The q = −1 charge state indicates that the P Cu −1 further traps electrons from the surrounding lattice besides the inner electrostatic interaction in the pairs.The Cl-poor condition, where ∆H f is elevated by 5.36 eV depending on the µ Cl , as shown in Figure S3, hinders the TM SACs from anchoring on the BiOCl.
Nanomaterials 2024, 14, x FOR PEER REVIEW 5 of 14 can see that the VCl is energy favorable in these three kinds of Cl vacancies if Ef is above the mid-gap (i.e., for the e-rich condition).In the V3Cl, 3/16 of Bi atoms are uncovered by the surface Cl atoms so that the V3Cl is large enough with regard to the unlocked Bi region.Consequently, the VCl is more stable than the unlocked Bi region in the BiOCl cathode.The VCl is capable of forming a donor defect with q = +1 charge state (VCl +1 , with the superscript referring to the charge state) when the Ef is in most of the gap and becomes a neutral defect only when the Ef is closer to the CBM.The charge transition level Ɛ(+1/−1) is only 0.07 eV above the CBM, indicating that it is a donor defect.The Fermi level of the VCl 0 is located in the conduction band in the spin-polarized density of states, which also indicates the donor defect characteristic.
Considering the stability of the VCl, we first investigate the PTM, which pairs the BiTM and surrounding single VCl, as shown in Figure 1c.The formation energy ΔHf is calculated to investigate the stability and ionization of the PTM and BiTM.If Ef is close to the CBM, as shown in Figure 2, the BiTM −1 dominates the BiTM for all the Fe, Co, Ni, and Cu SACs, indicating the p-type impurities characteristic of these TM SACs on pristine BiOCl in the erich condition.When the positively charged VCl +1 and negatively charged BiTM −1 move close to each other, as shown in Figure 1c, n-p codoping is established due to the coulomb attraction and forms the PTM pairs.Due to the electrostatic attraction of the n-p codoping, the ΔHf of the PTM is lower than the relative BiTM as shown in Figure 2. Consequently, the VCl facilitates the stability of the TM SACs on the BiOCl.When Ef is close to the CBM, neutral pairs (PTM 0 ) are energetically favorable for the Fe, Co, and Ni SACs, while the negatively charged PCu −1 with charge state q = −1 is most stable for the Cu SACs.The PCu −1 possesses the lowest Ef (as low as −1.56 eV) in all the defects, suggesting that the PCu −1 can exist in great quantities on the BiOX cathode.The q = −1 charge state indicates that the PCu −1 further traps electrons from the surrounding lattice besides the inner electrostatic interaction in the pairs.The Cl-poor condition, where ΔHf is elevated by 5.36 eV depending on the µCl, as shown in Figure S3, hinders the TM SACs from anchoring on the BiOCl.The electrostatic interaction is the key factor for establishing n-p codoping.We investigate the electrostatic interaction in the pairs by calculating the charge density redistribution.The charge density redistribution is determined as the charge density difference of the BiOCl system after and before forming P TM pairs.In Figure 3, the yellow and blue isosurfaces refer to the electron accumulation and dissipation regions, respectively.Before combining into P TM pairs, isolated TM SACs are all negatively charged, and V Cl is positively charged based on the ∆H f in Figures 2 and S3a.The negative charge density localizes around the TM SACs in the Bi TM , while the positive charge localizes around the vacancy site in the V Cl as shown in Figure 3e,f.In all the P TM pairs, negative and positive charges also separated at the adjacent TM SACs and V Cl site as shown from Figure 3a-d.An obvious coulomb interaction is elicited between the closer positive and negative charges, forming more stable n-p codoping.
Nanomaterials 2024, 14, x FOR PEER REVIEW 6 of 14 The electrostatic interaction is the key factor for establishing n-p codoping.We investigate the electrostatic interaction in the pairs by calculating the charge density redistribution.The charge density redistribution is determined as the charge density difference of the BiOCl system after and before forming PTM pairs.In Figure 3, the yellow and blue isosurfaces refer to the electron accumulation and dissipation regions, respectively.Before combining into PTM pairs, isolated TM SACs are all negatively charged, and VCl is positively charged based on the ΔHf in Figures 2 and S3a.The negative charge density localizes around the TM SACs in the BiTM, while the positive charge localizes around the vacancy site in the VCl as shown in Figure 3e,f.In all the PTM pairs, negative and positive charges also separated at the adjacent TM SACs and VCl site as shown from Figure 3a-d.An obvious coulomb interaction is elicited between the closer positive and negative charges, forming more stable n-p codoping.

Modulation of CO2 Reduction
CO2 absorption (*CO2) is usually the rate-determining step in the CO2RR [15,18].We investigate CO2 absorption on defective BiOCl systems to evaluate the effects of n-p codoping on the photocatalytic CO2RR.From Figure 4a, one can see that the CO2 is physically absorbed on the pristine BiOCl surface, and hardly any charges transfer between the BiOCl and the CO2.Only the surface charge density distribution of the BiOCl is disturbed by absorbed CO2.The charge transfer between the VCl +1 site and absorbed CO2 is enhanced, but still only physical absorption occurs with over 3.5 Å of distance as shown in Figure 4b.The Gibbs free energy difference (ΔG) is also calculated based on Equations ( 5) and ( 6) to evaluate the CO2 absorption.Negative ΔG refers to the chemisorption of CO2 on the photocatalyst.In Figure S4, one can see that the ΔG is 0.16 eV on the VCl +1 , lower than that on pristine BiOCl (0.32 eV).Consequently, the VCl +1 mildly improves the CO2 absorption.For the CO2 absorption on the dopant-defect pairs, dramatic charge transfers are found in the neutral PFe 0 , PCo 0 , and the negatively charged PCu −1 as shown in Figure 4c, 4d, and 4e, respectively.The lengths of C-TM bonds between the absorbed CO2 and these pairs are 1.90 Å, 1.97 Å, and 1.97 Å as shown in Table 1, while the ΔG sharply declines to −0.52 eV, −0.83 eV, and −0.54 eV for the PFe 0 , PCo 0 , and PCu −1 , respectively, as shown in Figure

Modulation of CO 2 Reduction
CO 2 absorption (*CO 2 ) is usually the rate-determining step in the CO 2 RR [15,18].We investigate CO 2 absorption on defective BiOCl systems to evaluate the effects of n-p codoping on the photocatalytic CO 2 RR.From Figure 4a, one can see that the CO 2 is physically absorbed on the pristine BiOCl surface, and hardly any charges transfer between the BiOCl and the CO 2 .Only the surface charge density distribution of the BiOCl is disturbed by absorbed CO 2 .The charge transfer between the V Cl +1 site and absorbed CO 2 is enhanced, but still only physical absorption occurs with over 3.5 Å of distance as shown in Figure 4b.The Gibbs free energy difference (∆G) is also calculated based on Equations ( 5) and ( 6) to evaluate the CO 2 absorption.Negative ∆G refers to the chemisorption of CO 2 on the photocatalyst.In Figure S4, one can see that the ∆G is 0.16 eV on the V Cl +1 , lower than that on pristine BiOCl (0.32 eV).Consequently, the V Cl +1 mildly improves the CO 2 absorption.For the CO 2 absorption on the dopant-defect pairs, dramatic charge transfers are found in the neutral P Fe 0 , P Co 0 , and the negatively charged P Cu −1 as shown in Figure 4c, Figure 4d, and Figure 4e, respectively.The lengths of C-TM bonds between the absorbed CO 2 and these pairs are 1.90 Å, 1.97 Å, and 1.97 Å as shown in Table 1, while the ∆G sharply declines to −0.52 eV, −0.83 eV, and −0.54 eV for the P Fe 0 , P Co 0 , and P Cu −1 , respectively, as shown in Figure 5.As a result, the n-p codoping facilitates the absorption of the CO 2 on the P Fe 0 , P Co 0 , and P Cu −1 .The charge transfers in the P Ni 0 are slightly weaker, and the distance between the absorbed CO 2 and Ni SACs (3.31 Å) is larger than that in other P TM s, but the ∆G still declines to −0.31 eV.Consequently, the CO 2 absorption is also accelerated by the n-p codoping in the P Ni 0 .
5. As a result, the n-p codoping facilitates the absorption of the CO2 on the PFe 0 , PCo 0 , and PCu −1 .The charge transfers in the PNi 0 are slightly weaker, and the distance between the absorbed CO2 and Ni SACs (3.31 Å) is larger than that in other PTMs, but the ΔG still declines to −0.31 eV.Consequently, the CO2 absorption is also accelerated by the n-p codoping in the PNi 0 .For comparison, we also evaluated the CO 2 absorption on the isolated Bi TM 0 and Bi TM −1 .On the neutral Bi TM 0 , the CO 2 is only absorbed physically on the TM SAC sites with fewer charge transfers as seen in Figure S5.The distances between the CO 2 and Bi TM 0 sites are all beyond 4 Å and are obviously larger than the distances between the CO 2 and the dopant-defect pairs as shown in Table 1.The ∆Gs of CO 2 absorption on the Bi TM 0 are all positive as shown in Figure 5.The CO 2 is therefore still physically absorbed at the neutral Bi TM 0 site.For the CO 2 absorption on the negative Bi TM −1 , we find that the CO 2 is only chemically adsorbed at the Bi Fe −1 but is still physically absorbed at the Bi Co −1 , Bi Ni −1 , and Bi Cu −1 .In detail, an intense charge transfer is formed between the CO 2 and Bi Fe −1 site as shown in Figure S6a, while the length of the C-Fe bond is 1.90 Å as shown in Table 1.The ∆G of CO 2 absorption is −0.18 eV in the Bi Fe −1 .On the other hand, relatively few charge transfers are formed between the CO 2 and Bi Co −1 , Bi Ni −1 , and Bi Cu −1 sites as shown in Figure S6b, 6c, and 6d, respectively.The distances between CO 2 and these activation sites are over 4.0 Å as shown in Table 1.The relative ∆Gs of CO 2 absorption are 0.15 eV, 0.23 eV, and 0.12 eV in the Bi Co −1 , Bi Ni −1 , and Bi Cu −1 systems.Based on previous experiments, the carboxylate pathway is the most common route for the CO 2 RR on Fe-, CO-, and Cu-doped bismuth oxyhalide [16][17][18][19][20]. *CO 2 hydrogenation (*COOH) follows close behind CO 2 absorption (*CO 2 ) in the carboxylate pathway and is reported to be also an important step influencing CO 2 reduction [15,18].We therefore also assess the ∆G of both *CO 2 and *COOH in all defective BiOCl systems.In Figure 5, one can see that n-p codoping facilitates both the CO 2 absorption and *CO 2 hydrogenation with negative ∆G in the BiOCl systems by Fe, Co, and Cu SAC doping.Such an enhancement of CO 2 photoreductions corresponds with previous experiments [18][19][20].For the P Ni 0 -doped BiOCl, the CO 2 absorption is exothermic, but the *CO 2 hydrogenation is endothermic, which is related to the fact that the CO 2 activation prefers the carbide pathway, in which the *COOH is difficult to generate [17].

Orbital Alignments and d-Band Center of the TM_3d
The TM_3d orbital is usually crucial to the reducing capacity of TM SACs [54,55].To account for the improvement in the CO 2 absorption by the n-p codoping, we investigate the d-band center of the TM_3d orbital states on the BiOCl systems.Based on the d-band theory, the d-band center close to the Fermi level usually facilitates the reducing capacity of TM SACs [30,31].For comparison, the d-band centers of the Bi TM 0 and Bi TM −1 are also calculated as shown in Table 1.One can see that the d-band centers of all Bi TM 0 dopants are located far away from the Fermi level.For instance, the d-band centers of the Bi Fe 0 , Bi Ni 0 , and Bi Cu 0 are as low as −2.54 eV, −3.85 eV, and −2.81 eV, respectively.Low d-band centers correspond to the low reducing capacity of the Bi TM 0 and result in physical absorption of the CO 2 .For the Bi TM −1 , the d-band center of Bi Fe −1 is significantly shifted from −2.54 eV to −1.26 eV compared with Bi Fe 0 , resulting in the enhancement from physical absorption to chemical absorption.The d-band center in the Bi Co −1 decreases from −1.81 eV to −2.45 eV, so the CO 2 remains physically absorbed on the Bi Co −1 site.Although the d-band centers of the Bi Ni −1 and Bi Cu −1 are shifted toward the Fermi level compared with the Bi Ni 0 and Bi Cu 0 , they are still lower than the Fermi level by −2.09 eV and −2.31 eV.Consequently, the CO 2 is also physically absorbed on these systems.As shown in Table 1, the d-band centers are upshifted to −1.27, −0.83, −1.46, and −1.45 eV for the P Fe 0 , P Co 0 , P Ni 0 , and P Cu −1 , respectively, and are comparable to those (−1.56eV) of the Pt SACs anchored at the edge of graphene, which are excellent photocatalysts for H 2 reduction [56].Consequently, the modulation of CO 2 absorption on the dopant-defect pairs derives from the enhancement of the reducing capacity of TM SACs with the upshift in the d-band center by the n-p codoping.
The d-band center is determined by the orbital alignment and occupation of the TM-3d.To account for the effects of n-p codoping on the d-band centers, we further calculated the spin-resolved states of density (PDOSs) of the BiOCl systems before and after n-p codoping to illustrate the orbital alignment and occupation.Based on the crystal field theory, the degeneracy of the TM-3d orbital states is broken in a process dependent on the molecular symmetry.As shown in Figure 1b,c In the neutral Bi TM 0 , the oxidation states of the TM dopant are all +3.As shown in Figure 6a, the d electron configuration obeys Hund's rule in Bi Fe 0 .Five spin-up d orbitals are fully occupied, while five spin-down d orbitals are entirely empty.In the Bi Co 0 , Bi Ni 0 , and Bi Cu 0 , the d electron configurations break Hund's rule as shown in Figure 6a.In detail, four spin-up d orbitals and only the spin-down d xz + d yz orbital are occupied in the Co_3d, when the electron configuration of Co is changed from 3d 7 4s 2 to 3d 6 4s 0 after losing three valence electrons in the Bi Co 0 .In the Bi Ni 0 , four spin-up d orbitals and three spin-down d orbitals (d xz + d yz orbitals and d xy orbital) are occupied to form the 3d 7 4s 0 configuration.In the Bi Cu 0 , the Cu_3d orbitals are spin-degenerated, and spin-degenerated d xz + d yz , d xy and d z 2 orbitals are occupied to change the 3d 10 4s 1 to a 3d 8 4s 0 configuration.
In the dopant-defect pairs P TM 0 , negatively charged TM dopants and positively charged V Cl +1 are stabilized by the electrostatic interaction of the n-p codoping.The electron accumulation changes the oxidation state of the TM dopant from +3 to +2.Compared with the electron configuration in the Bi Tm 0 , electron accumulation is located around the TM while charge depletion forms around the V Cl +1 and durable electrostatic attraction is established based on the Coulomb's law as shown in Figure 4. Accumulated electrons occupy the localized empty TM_3d orbital as shown in Figure 6b.In detail, spin-down d z 2 is further occupied in the P Fe 0 to form the 3d 6 4s 0 configuration.In the P Co 0 , the electron from the V Cl +1 further occupies the spin-down d xy to form the 3d 7 4s 0 configuration.Spin-down Ni_d z 2 is further occupied, leading to spin-degenerated Ni_3d states in the P Ni 0 , while the spin degeneration of the Cu_3d states is broken along with the spin-up d x The electron accumulation around the TM SACs plays an important role in changing the d-orbital alignment and shifting the d-band center in the PTM.In contrast to the BiFe −1 , the reducing capacity in the BiTM −1 is weaker than that in the PTM.To account for the reduced capacity in the BiTM −1 , we investigate its d orbital alignment and electron configuration in detail.The electron accumulation is also formed around the BiFe −1 as shown in Figure S7a.The spin-polarized PDOS of TM_3d in the BiTM −1 is shown in Figure S8a.The BiFe −1 keeps in line with the PFe 0 on the d-orbital alignment, the 3d 6 4s 0 configuration, and even Hund's rule, resulting in a similar d-band center and ΔG for the CO2 absorption.Although electron accumulation occurs around the BiCo −1 as shown in Figure S7b, the 3d 7 4s 0 configuration stays the same as that in the PCo 0 , but the spin-polarized d-orbital alignments exhibit a dramatic change compared with those in the PCo 0 .Five spin-up orbitals are occupied in the BiCo −1 as shown in Figure S8a, while only four spin-up orbitals are occupied in the PCo 0 (the same as in the BiCo 0 ) as shown in Figure 6.Such variation changes the high spin state (3µB) to a low spin state (1µB) as well as induces greater broadening of Co_3d orbital states and then a more negative d-band center in the BiCo −1 .The electron localization is weakest in the BiNi −1 as shown in Figure S7c.As a result, although the electron configuration in the Ni_3d is same as the cases in the PCo 0 , the energy level and broadening of Ni_3d, which also influence the d-orbital alignment, are extended to a low energy level in the BiNi −1 , resulting in a much lower d-band center than that in the PCo 0 .The BiCu −1 possesses the same orbital alignment, orbital occupation, and d-band center as the PCu 0 as shown in Figure S8d.Therefore, it possess a similar ΔG for CO2 absorption.The stronger reducing capacity of the PCu −1 derives from the accumulation of additional electrons around the Cu SACs.On the other hand, compared with the BiTM 0 , BiFe −1 , BiNi −1 , and BiCu  The electron accumulation around the TM SACs plays an important role in changing the d-orbital alignment and shifting the d-band center in the P TM .In contrast to the Bi Fe −1 , the reducing capacity in the Bi TM −1 is weaker than that in the P TM .To account for the reduced capacity in the Bi TM −1 , we investigate its d orbital alignment and electron configuration in detail.The electron accumulation is also formed around the Bi Fe −1 as shown in Figure S7a.The spin-polarized PDOS of TM_3d in the Bi TM −1 is shown in Figure S8a.The Bi Fe −1 keeps in line with the P Fe 0 on the d-orbital alignment, the 3d 6 4s 0 configuration, and even Hund's rule, resulting in a similar d-band center and ∆G for the CO 2 absorption.Although electron accumulation occurs around the Bi Co −1 as shown in Figure S7b, the 3d 7 4s 0 configuration stays the same as that in the P Co 0 , but the spinpolarized d-orbital alignments exhibit a dramatic change compared with those in the P Co 0 .Five spin-up orbitals are occupied in the Bi Co −1 as shown in Figure S8a, while only four spin-up orbitals are occupied in the P Co 0 (the same as in the Bi Co 0 ) as shown in Figure 6.Such variation changes the high spin state (3µ B ) to a low spin state (1µ B ) as well as induces greater broadening of Co_3d orbital states and then a more negative d-band center in the Bi Co −1 .The electron localization is weakest in the Bi Ni −1 as shown in Figure S7c.As a result, although the electron configuration in the Ni_3d is same as the cases in the P Co 0 , the energy level and broadening of Ni_3d, which also influence the d-orbital alignment, are extended to a low energy level in the Bi Ni −1 , resulting in a much lower d-band center than that in the P Co 0 .The Bi Cu −1 possesses the same orbital alignment, orbital occupation, and d-band center as the P Cu 0 as shown in Figure S8d.Therefore, it possess a similar ∆G for CO 2 absorption.The stronger reducing capacity of the P Cu −1 derives from the accumulation of additional electrons around the Cu SACs.On the other hand, compared with the Bi TM 0 , Bi Fe −1 , Bi Ni −1 , and Bi Cu −1 upshift the d-band center because additional electrons accumulate around the TM SACs and occupy more 3d orbitals.In contrast to the Bi Co 0 , the Bi Co −1 moves the d-band center downward with negative energy following the transition of spin states.Consequently, the electron accumulation plays a key role in the d-orbital alignment and the d-band center.The transition of spin states also influences the d-band center.

Discussion
Finally, the effect of Cl multivacancies on the n-p codoping is discussed.Taking the most sable P Cu −1 as an example, the P Cu@2Cl and P Cu@3Cl are pairs with double and ternary Cl vacancies, as shown in Figure S1.Both the P Cu@2Cl and P Cu@3Cl exhibit lower formation energy than the P Cu −1 .However, a high ∆H f hinders the stability and solubility of V 2Cl and V 3Cl , as shown in Figure S2.The solubility of both the P Cu@2Cl and P Cu@3Cl is lower than that of the P Cu −1 .In the e-rich condition, both the P Cu@2Cl and P Cu@3Cl exhibit the charge state of q = 0, as shown in Figure S2b (denoted as P Cu@2Cl 0 and P Cu@3Cl 0 ).The ∆G profiles in Figure 7a exhibit obvious high free energy (0.61 eV and 0.43 eV) for CO 2 absorption at both the P Cu@2Cl 0 and P Cu@3Cl 0 sites, which is higher than that at the P Cu −1 and even higher than that at the V Cl +1 .Relatively high d-band centers lead to a weak reducing capacity as shown in Table 1.Obviously positive charge accumulation around the P Cu@2Cl 0 and P Cu@3Cl 0 is also observed as shown in Figure 7c,d, shielding the electron transfer to absorbed CO 2 and decreasing the reducing capacity.Consequently, the P Cu −1 enhances the CO 2 absorption, and doping the TM on the unlocked Bi region is not an efficient strategy to facilitate the CO 2 RR.

Discussion
Finally, the effect of Cl multivacancies on the n-p codoping is discussed.Taking the most sable PCu −1 as an example, the PCu@2Cl and PCu@3Cl are pairs with double and ternary Cl vacancies, as shown in Figure S1.Both the PCu@2Cl and PCu@3Cl exhibit lower formation energy than the PCu −1 .However, a high ΔHf hinders the stability and solubility of V2Cl and V3Cl, as shown in Figure S2.The solubility of both the PCu@2Cl and PCu@3Cl is lower than that of the PCu −1 .In the e-rich condition, both the PCu@2Cl and PCu@3Cl exhibit the charge state of q = 0, as shown in Figure S2b (denoted as PCu@2Cl 0 and PCu@3Cl 0 ).The ΔG profiles in Figure 7a exhibit obvious high free energy (0.61 eV and 0.43 eV) for CO2 absorption at both the PCu@2Cl 0 and PCu@3Cl 0 sites, which is higher than that at the PCu −1 and even higher than that at the VCl +1 .Relatively high d-band centers lead to a weak reducing capacity as shown in Table 1.Obviously positive charge accumulation around the PCu@2Cl 0 and PCu@3Cl 0 is also observed as shown in Figure 7c,d, shielding the electron transfer to absorbed CO2 and decreasing the reducing capacity.Consequently, the PCu −1 enhances the CO2 absorption, and doping the TM on the unlocked Bi region is not an efficient strategy to facilitate the CO2RR., PCu@2Cl 0 , and PCu@3Cl 0 .(b-d) depict the charge density redistribution of the PCu −1 , PCu@2Cl 0 , and PCu@3Cl 0 , respectively.The isovalue is set to 0.03 e/Å.

Conclusions
In summary, n-p codoping engineering is introduced to account for the modulation of photocatalytic CO2 reduction on a BiOCl-based cathode by using first-principles calculation.n-p codoping is established via the Coulomb interactions between positively charged TM SACs and the negatively charged Cl vacancy (VCl) in the dopant-defect pairs.Based on the formation energy of charged defects, we find neutral dopant-defect pairs for the Fe, Co, and Ni SACs (PTM 0 ) and the q = −1 charge state of the Cu SAC-based pair (PCu −1 ).

Figure 1 .
Figure 1.Schematic diagrams of atomic structures for 2D BiOCl with (a) V Cl , (b) Bi TM , and (c) P TM .Blue circle refer to the P TM .

Figure 2 .
Figure 2. Formation energies of BiTM (red lines) and PTM (blue lines) at the Cl-rich limit for Fe, Co, Ni, and Cu.The numerical notations refer to charge states.The dotted lines refer to the formation energies of meta-stable charged defects.

Figure 2 .
Figure 2. Formation energies of Bi TM (red lines) and P TM (blue lines) at the Cl-rich limit for Fe, Co, Ni, and Cu.The numerical notations refer to charge states.The dotted lines refer to the formation energies of meta-stable charged defects.

Figure 3 .
Figure 3. Charge density redistribution of (a) P Fe 0 , (b) P Cu 0 , (c) P Ni 0 , and (d) P Cu 0 .(e,f) The charge density distribution of the isolated Bi Cu −1 and V Cl −1 .The yellow (blue) isosurface refers to the charge accumulation (dissipation) region.The isovalue is set to 0.03 e/Å.
, the TM ion is located at the center of a square plane formed by four nearest neighboring oxygen atoms.In this local symmetry, the TM_3d orbitals split into four groups d xz + d yz , d xy , d z 2 , and d x 2 −y 2 .
−1 upshift the d-band center because additional electrons accumulate around the TM SACs and occupy more 3d orbitals.In contrast to the BiCo 0 , the BiCo −1 moves the d-band center downward with negative energy following the transition of spin states.Consequently, the electron accumulation plays a key role in the d-orbital alignment and the d-band center.The transition of spin states also influences the d-band center.