Decoupling Many-Body Interactions in CeO2 (111) Oxygen Vacancy Structure: Insights from Machine-Learning and Cluster Expansion

Technical Summary: Decoupling Many-Body Interactions in CeO2 (111) Oxygen Vacancy Structure

Problem Statement
Oxygen vacancies ($V_O$) are critical to the redox and catalytic properties of ceria ($CeO_2$). However, accurately predicting the distribution and stability of these vacancies on the $CeO_2(111)$ surface is computationally prohibitive due to the vast number of possible electronic configurations. As the concentration of vacancies increases, the number of associated $Ce^{3+}$ ions (polarons) rises, leading to an exponential increase in configurations. The stability of these structures depends on complex many-body interactions, including $Ce^{3+}$-$Ce^{3+}$, $V_O$-$Ce^{3+}$, and $V_O$-$V_O$ repulsions and attractions. Traditional methods struggle to sample these configurations effectively, and previous studies have largely focused on isolated vacancies or low concentrations, often overlooking the potential stability of vacancies in deeper oxygen layers (beyond the subsurface).

Methodology
The authors developed a hybrid theoretical framework combining first-principles calculations, machine learning (ML), and cluster expansion (CE) to decouple these many-body interactions.

1. First-Principles Calculations: Spin-polarized DFT+U calculations (PBE functional, $U=5.0$ eV) were performed on $CeO_2(111)$ slabs. The study focused on the top three oxygen layers and the top two cationic layers. $Ce^{3+}$ identification was based on local magnetic moments ($\sim 1 \mu_B$) and $Ce$-$O$ bond lengths.
2. Cluster Expansion Model: A CE model was constructed to represent the formation energy of any configuration as a sum of effective cluster interactions ($J_\alpha$). The initial pool consisted of 3,144 clusters involving up to four points (representing $V_O$ or $Ce^{3+}$) within an 8 Å cutoff.
3. Machine Learning (LASSO): The Least Absolute Shrinkage and Selection Operator (LASSO) regression was employed to select the optimal subset of clusters, preventing overfitting and identifying the most significant interaction terms. The training set was iteratively refined using an active learning loop involving Metropolis Monte Carlo (MC) sampling at 1000 K.
4. Monte Carlo Sampling: Once the CE model was converged (validated by a prediction error threshold of <0.05 eV/$V_O$), it was used to perform extensive MC simulations (over $10^8$ steps) in an $8 \times 8$ supercell. These simulations explored vacancy distributions across various concentrations ($\Theta = 1/16, 1/8, 3/16, 1/4$) and temperatures (300 K to 1100 K).

Key Contributions and Results

• Decoupling Interactions: The study successfully isolated 15 key clusters from the initial 3,144, effectively decoupling the many-body interactions. The model achieved a maximum absolute error (MAE) of 0.07 eV/$V_O$ on the training set and 0.08 eV/$V_O$ on cross-validation.
• Identification of Motifs: The analysis revealed specific local motifs governing stability:
  • Isolated Vacancies: Subsurface vacancies are favored over surface vacancies by $\sim 0.2$ eV.
  • Interactions: Strong repulsion exists between nearest-neighbor vacancies in different layers. Conversely, specific patterns, such as a subsurface vacancy paired with a third-layer vacancy (motif $VO_2VO_3\_5N$), exhibit strong attraction ($-0.16$ eV).
  • Polaron Localization: $Ce^{3+}$ ions generally prefer surface sites but show strong attraction to next-nearest-neighbor (NNN) vacancies.
• Concentration-Dependent Stability:
  • At low concentrations ($\Theta \le 1/8$), vacancies predominantly occupy the subsurface layer, consistent with prior literature.
  • At higher concentrations ($\Theta \ge 3/16$), vacancies begin to populate the third oxygen layer. The most stable configuration at $\Theta = 1/4$ (termed A2233) features two vacancies in the subsurface and two in the third layer. This structure is 0.09 eV more stable than the previously proposed all-subsurface $(2 \times 2)$ ordered structure.
• Mechanism of Stability: The stability of the third-layer vacancies at high concentrations is driven by extensive geometric relaxation. While bond-breaking energy ($E_b$) penalizes third-layer vacancies, the relaxation energy ($E_r$) gained in the A2233 configuration (where surface $Ce^{4+}$ ions can relax upward more freely due to the absence of opposing bonds) outweighs this penalty. In contrast, configurations with all vacancies in the third layer (A3333) are unstable due to prohibitively high bond-breaking energies.
• Temperature Effects: MC simulations show that as temperature and concentration increase, the probability of finding vacancies in the third layer increases relative to the surface layer, although surface vacancies do appear occasionally at high temperatures.

Significance
The paper claims that this work introduces a novel theoretical framework for unraveling complex vacancy structures in metal oxides by decoupling many-body interactions. The primary significance lies in the discovery that oxygen vacancies in reduced $CeO_2(111)$ do not remain confined to the surface or subsurface at high concentrations but migrate to the third oxygen layer to maximize geometric relaxation. This finding suggests that previous interpretations of ceria surface chemistry, which neglected the third layer, may require re-evaluation. The coupled machine-learning and cluster expansion approach is presented as a robust method for studying other metal oxides with well-defined lattice structures, offering a pathway to understand redox and catalytic mechanisms more accurately.