Unravelling the Impact of Metal Dopants and Oxygen Vacancies on Syngas Conversion over Oxides: A Machine Learning-Accelerated Study of CO Activation on Cr-Doped ZnO Surfaces

As a critical component of the OX–ZEO composite catalysts toward syngas conversion, the Cr-doped ZnO ternary system can be considered as a model system for understanding oxide catalysts. However, due to the complexity of its structures, traditional approaches, both experimental and theoretical, encounter significant challenges. Herein, we employ machine learning-accelerated methods, including grand canonical Monte Carlo and genetic algorithm, to explore the ZnO(1010) surface with various Cr and oxygen vacancy (OV) concentrations. Stable surfaces with varied Cr and OV concentrations were then systematically investigated to examine their influence on the CO activation via density functional theory calculations. We observe that Cr tends to preferentially appear on the surface of ZnO(1010) rather than in its interior regions and Cr-doped structures incline to form rectangular islands along the [0001] direction at high Cr and OV concentrations. Additionally, detailed calculations of CO reactivity unveil an inverse relationship between the reaction barrier (Ea) for C–O bond dissociation and the Cr and OV concentrations, and a linear relationship is observed between OV formation energy and Ea for CO activation. Further analyses indicate that the C–O bond dissociation is much more favored when the adjacent OVs are geometrically aligned in the [1210] direction, and Cr is doped around the reactive sites. These findings provide a deeper insight into CO activation over the Cr-doped ZnO surface and offer valuable guidance for the rational design of effective catalysts for syngas conversion.


Computational methods
As shown, the value of U was determined by benchmarking with HSE data and was found to be 3.0 eV for Cr and Zn better, which adequately described the  OV for various systems containing Zn or Cr. OV was selected for benchmarking due to its significance in evaluating the equilibrium OV concentrations under the reaction condition, which has a great impact on the C-O bond activation.

ML-accelerated global optimization methods
Table S1.Energy and force root mean squared error (RMSE) on DFT-labelled structures.The Nframes is the number of the training structures, as distinguished by the chemical formula (system) and the number of atoms per cell (Natoms).

System
Natoms Nframes ERMSE (eV/atom) where  DFT can be obtained from the DFT total energy.Vibrational frequency calculations were performed to determine the zero-point energy (∆ ZPE ).() is the enthalpy correction and () is the entropy, taken from the NIST database. 1r CO-assisted OV formation, the pressure of CO2 was evaluated at 673 K, 0.0765 Mpa, according to the experimental results (CO conversation at ~17% and CO2 selectivity at ~45%). 2 The detailed results are listed in Table S2.Under the reaction conditions, these surfaces should undergo reduction by reacting with CO/H2 to form CO2/H2O.The effect of  O on the equilibrium  OV was investigated through the calculation of  OV (Figures S9 and S10).On the unreconstructed surfaces, more OVs can be observed on the surfaces with low  Cr (1.00 M OVs for 0%  Cr and 0.00 M OVs for 12.5%  Cr under ∆ O = -3.5 eV ) and  OV increases gradually with the decline of  O , as shown in Figure S11a.An intriguing observation is that 0.33 M OVs can be found when  Cr is between 4.2% to 8.3% under ∆ O =-3.1 eV, whereas no OV formation is observed in the lower  OV ranging from 0% to 4.2%.This finding can be attributed to the promotion of OV generation around the Zn ions due to the incorporation of Cr, as mentioned previously.After structural reconstruction, the overall trends remain similar, except that the reduction degree may be greater under the same ∆ O (0.00 and 0.33 M OVs for 12.5%  Cr under ∆ O = -3.3eV), as can be seen in Figure S11b.These results suggest, therefore, that one may manipulate  OV of the surface by adjusting  O (temperature and pressure of the experiment).

CO activation
Table S3.Energy differences between the most stable structures produced by the offgrid and on-grid methods.The negative sign means that the structure from off-grid is more stable.
OV  cr 0.00 (ML) 0.33 (ML) 0.67 (ML) 1.00 (ML) 0.0 (at.%) 0.00 0.00 0.00 0.00 2.1 (at.%) 0.00 0.00 0.00 0.00 4.2 (at.%) 0.00 0.00 0.79 1.31 6.3 (at.%) 0.00 0.00 0.55 1.04 8.3 (at.%) 0.00 0.00 0.86 1.61 10.4 (at.%) 0.00 0.92 1.40 1.17 12.5 (at.%) 0.00 1.03 0.80 0.76 We find that CO-induced stabilization effect for on-grid structures outweighs the energetic favourability of the off-grid-derived surface, leading to a reversal in the relative stability of the two systems.Firstly, we calculated the energy difference between most stable structures produced by the off-grid and on-grid methods (Table S3), showing that the structures identified by the off-grid methods exhibit superior stabilities compared to those identified by the on-grid methods.Then we investigated the energy difference between the total energy of CO adsorption on these two types of structures, as summarized in Table S4.Our findings demonstrate that upon CO adsorption, the on-grid structures exhibit greater stabilities than the off-grid structures.It can be seen from Figure 8e in the main text that  a does not fit into a perfect linear relationship with  OV .Especially, for some grey points that covers a range of  Cr from 0 to 8.3% with 0.33 ML OVs, as  Cr increases,  OV decreases (see Figure 4b in the main text).This is because the surfaces with higher  Cr possess the ability to reduce the barrier for C-O bond dissociation: the surfaces with higher  Cr provide more opportunities for the cleavage products to be stabilized by the presence of Cr ions on the surfaces (Figure S14).

𝜇Figure S3 .
Figure S3.Changes of the total energies (black line) and  Cr (red line) of the surface region(Figure 2a) as a function of neural network accelerated GCMC simulation step for ZnO(101 ̅ 0) with 1%  Cr at 798 K.

Figure S4 .
Figure S4.Cr distribution function in z direction (()) for the structures obtained from GCMC simulations.()_0 refers to the () in the original structure, and subsequent () analyses were conducted at intervals of 2000 simulation steps.

Figure S7 .
Figure S7.Calculated density of states of ZnO with 2.1% Cr before(a) and after(b) introducing 0.33 ML OVs with respect to the Fermi level.

Figure S8 .
Figure S8.Most stable structures identified by on-grid and off-grid strategies.

Figure S9 .
Figure S9.Gibbs free energy of OV formation( OV ) with various  Cr and  OV based on the on-grid strategy(a,c) and the off-grid strategy(b,d) under various ∆  (-2.9, -3.1 eV).

Figure S10 .
Figure S10. OV with various  Cr and  OV based on the on-grid strategy(a,c) and the off-grid strategy(b,d) under various ∆ O (-3.3, -3.5 eV).

Figure S11 .
Figure S11.Equilibrium  OV of Cr-doped ZnO identified by on-grid (a) and off-grid (b) methods under different ∆ O .

Figure S13 .
Figure S13.Effective energy barrier( a ) differences of the C-O bond dissociation on the structures identified by the off-grid and on-grid structure search methods.

Figure S15 .
Figure S15.Configuration with a connected 3-coordinated surface OV and a 4coordinated subsurface OV.The blue dashed circle indicates the OV.

Figure S16 .
Figure S16.Optimized transition state of C-O bond dissociation on ZnO[10 1 ̅ 0 surfaces with various arrangements of OVs.

Table S2 .
Correction to gibbs free energies of the molecules.Under the equilibrium condition, the  O = -8.29 eV.In our work, to see the influence of reducing environment more clearly, we defined ∆ O =  O −  ̃O, where  ̃O is the chemical potential of oxygen under the standard condition (298.15K, 0.1 Mpa).