Band alignment of β-(Al_xGa_1−x)₂O₃ alloys via atomic solid-state energy scale approach

The band lineup or alignment of semiconductors is very important for new device design. When semiconductors get in contact, band-offset values are also important physical and chemical properties. Therefore, to determine band alignments for various semiconductors, several theoretical approaches have been presented.^1,2 For example, the band lineup relative to the charge neutrality level (also known as the branch point energy)³ or the band alignment relative to the vacuum level calculated by the GW approach was demonstrated.⁴ However, the “natural band alignments” that have been calculated are most often typical compound semiconductors without strain or dipole effects.^5,6 The calculation methods of their natural band alignment were useful guidelines for screening of the heterojunction type or doping limit due to low computational cost.^7–10 Recently, Wetson et al. established a method to obtain the band-offset values by calculating the electrostatic potential alignment of the heterostructure superlattice and the constituent semiconductors.¹¹ This method has high computational efficiency, but it did not provide a case for estimating the band offset of compound semiconductor alloys. As far as we know, there have been few reports about the determination of band alignment of alloy materials by a simple theoretical model. Therefore, a simple method to estimate band alignment in compound alloy semiconductors for high throughput screening is needed.

Toroker et al. proposed a first principles approach with a simple equation that can be used to estimate the band edge positions.¹² The valence band maximum (VBM) and conduction band minimum (CBM) relative to the vacuum are determined by the energy of the bandgap center (E_BGC) and energy gap (E_g) of the material. This method is effective to determine the band edge positions of nonpolar crystals. A clear advantage of the method is that it explains the band edge positions just in terms of slab model calculations. However, this approach has limitations; the surface of the slab model cannot have a net electrical dipole moment. In contrast, the classical approach provided a similar empirical equation based on the electron affinity (EA) rule.¹³ Fortunately, the empirical method was found to be good in estimating the band edges of oxides and metal-sulfide materials without considering the surface state effect.¹⁴ However, this method is based on the Mulliken electronegativity (χ) of a material. χ is the absolute energy scale, which is defined by the arithmetic mean of the atomic first ionization energy (IE) and electron affinity (EA).¹⁵ Therefore, when one considers the crystalline state, the χ concept cannot suitably represent condensed atoms.

More recently, Pelatt et al. presented an atomic solid-state energy (SSE) scale, which is an alternative approach to electronegativity.^16,17 SSEs are calculated from a given atom by assessing an average EA (for a cation) or an average IE (for an anion) for solid-state binary inorganic compounds. We focused on the SSE concept, which is more suitable for the atomic energy state in a condensed matter than χ. Clarifying whether SSE can be applied to estimate the band edge of a semiconductor alloy is particularly important for developing a new heterodevice in the future. Moreover, it can also provide useful guidelines for the development of photocatalysts and photoelectrochemical cells.^18,19

In this study, we present a theoretical approach for estimating the band alignment of β-$(AlxGa1−x)2O3$ alloys. We observed that the use of the SSE scale in combination with density functional theory (DFT) calculations reproduced natural band alignment in accordance with the conventional electrostatic potential approach. Additionally, we show that the calculated band offsets for selected heterostructures were in agreement with the experimental values than the previous χ scheme. Our results indicated that simple and practical prediction of the band offset and alignment of ionic bonded semiconductors can be attained.

We used the monoclinic β-gallia structure, where the alloy models contained 20 atoms, 8 Ga atoms and 12 O atoms. We employed the “Supercell program” for the atomic substitutions in crystals to reduce the number of potential configurations by the algorithm for detecting symmetry-equivalent structures.²⁰ The geometry optimization was performed within the DFT of the plane wave based on the Quantum ESPRESSO package.^21,22 We used the Perdew–Burke–Ernzerhof (PBE) form within the generalized gradient approximation (GGA) with projector-augmented-wave potentials.^23,24 The electronic wave-functions and charge densities were expanded in plane waves with cutoff energies of 50 Rydberg and 350 Rydberg (Ry), respectively. The k-points were generated using the Monkhorst–Pack (MP) scheme with a mesh size of 4 × 8 × 6 for the structure optimized calculation.²⁵ For geometry optimization, the criteria for stress and force on atoms were set to be 0.05 GPa and 0.001 Ry/a.u., respectively. After the geometry optimization of the models, the most favorable alloy structures were determined using formation enthalpy. The formation enthalpy (ΔH) was calculated by

where E$Ga2O3$ and E$Al2O3$ are the energies per atom of the lowest structures, monoclinic for Ga₂O₃ and corundum for Al₂O₃. The results shown in Fig. 1 include all nonequivalent atomic substitution alloy models. The monoclinic structure has two cation sites: an octahedrally coordinated and a tetrahedrally coordinated site. However, previous studies have shown that the substitute Al prefers an octahedral site over a tetrahedral site.^26–28 Our results show similar trends only for the formation enthalpy minima structures. In the case of Al with 50% content, the energetically stable AlGaO₃ ordered structure appeared. Such an ordered alloy structure may have been realized in the synthesis of oxide powders.²⁶ For all alloys, the energetically preferable configurations are illustrated in Fig. 2. The crystal structures were drawn with VESTA.²⁹ These results clearly show that Al prefers to occupy the octahedral sites for x ≤ 0.5. When the octahedral site is filled, Al will be substituted on the tetrahedral site (x > 0.5). To attain high accuracy DFT calculations, the preferable structures were used to calculate the electronic structures by using the full potential (L)APW + lo method implemented in WIEN2k code.^30,31 For self-consistent field (SCF) calculations, we selected the muffin-tin radii of O, Al, and Ga to be 1.64–1.66 a.u., 1.66–1.76 a.u., and 1.81–1.83 a.u., respectively. The plane-wave cutoff was determined using R_MTK_max = 7.0, where R_MT is the muffin-tin radius and K_max gives the magnitude of the largest K vector in the plane-wave expansion. The Brillouin zone (BZ) was sampled on the MP mesh with 2000 k points in each BZ for SCF calculations within the GGA-PBE exchange–correlation energy functional. After SCF calculations, the E_g values were corrected by the Tran–Blaha modified Becke Johnson (TB-mBJ) potential approach with Koller et al. presented parameter.³² We chose the TB-mBJ approach because it shows better performance for predicting E_g values than the hybrid functional.³³ 

To estimate the band alignment in β-$(AlxGa1−x)2O3$ alloys, we followed Butler and Ginley who presented a simple equation for band edge determination using χ and E_g.¹³ For our studies here, we replaced the χ with the SSE of the alloys. The positions of the energy of the valence band maximum (E_VBM) and conduction band minimum (E_CBM) are obtained using the following relation:

where E_g is the bandgap and ξ_M is the SSE of the material, which is defined by the geometric mean of the SSE of constituent atoms. For a β-$(AlxGa1−x)2O3$ alloy, the ξ_M can be calculated as

Here, x is the composition of aluminum and ξ_Al, ξ_Ga, and ξ_O are the SSE of absolute values for atoms. We used the following SSE of absolute values for O, Al, and Ga: 7.98 eV, 3.14 eV, and 3.82 eV, respectively.¹⁷ 

Figure 3 shows the calculated bandgap of β-$(AlxGa1−x)2O3$ alloys as a function of Al content. The experimental bandgap values of β-Ga₂O₃ varied between 4.5 eV and 4.9 eV.^34–38 Our calculated value was 4.93 eV, which was slightly larger than the reported results. For monoclinic Al₂O₃ (x = 1), the crystal structure denoted θ-Al₂O₃.³⁹ The calculated bandgap for θ-Al₂O₃ was 7.0 eV, which falls within the range of 6.9–7.1 eV (amorphous to crystal) reported in the results of x-ray photoelectron spectroscopy.^40,41 The bowing parameter b is estimated from

We obtained a bowing parameter of about b ≃ 1.4 eV for the β-$(AlxGa1−x)2O3$⁠. The parameter values were not consistent with previous theoretical results.^27,28 The different E_g values of the start (4.69 eV, 4.895 eV) and endpoints (7.03 eV, 7.24 eV) caused a variety of b values. However, our calculated E_g results were in good agreement with the experimental data.^26,41–44 In addition, the bandgap values as a function of Al were also similar to the values of other hybrid functional calculation results except for the bowing parameter. Therefore, our calculated E_g values were available for band alignment estimation. The position of the band edge will vary depending on whether the bandgap is underestimated or overestimated. However, from Eqs. (2) and (3), the effect of underestimating (or overestimating) the bandgap is halved with respect to the positions of the VBM and CBM. For example, in the case of β-Ga₂O₃, the bandgap may be overestimated by ∼0.4 eV, but it only shifts by ∼0.2 eV each when determining the position of the VBM and CBM. This is considered accurate enough to account for the experimental errors that will be discussed later. The result is shown in Fig. 4, and the band edge positions were determined using Eqs. (2) and (3). The CBM positions monotonically increased with increasing Al content. The VBM positions only reduced slightly compared to the CBM positions. This was mainly because the O 2p states dominated the VBM. The only exception was the x = 0.125 case, for which the VBM position was slightly larger than that of Ga₂O₃. The overestimation for the Ga₂O₃ bandgap likely caused the out of linear trend. The differences of each band edge position due to the Al content were consistent with previous experimental results.⁴² However, other DFT calculation results showed a slightly upward VBM trend.^27,28 Distinct differences exist between our method and the previous approach based on the electrostatic potential.⁴⁵ For the bandgap correction method, a previous study used the computationally costly hybrid functional of Heyd, Scuseria, and Ernzerhof (HSE).^46,47 Notably, the estimation method for band edge positions using the averaged electrostatic potential procedure with 120 atoms was completely different from our presented scheme. The results of band alignment, despite different techniques being used for both band edge positions related to the vacuum, corresponded with our study and previous DFT calculations. Our presented method has the advantage that band alignment can be determined using only the E_g and SSE of materials. To assess our band alignment prediction scheme, we compared the band-offset value obtained by the SSE approach using the experimental E_g.

We inspected the available literature for experimental data on band offsets for Ga₂O₃ and related oxide.^{42–44,48–82} The calculated valence band offsets (VBOs) and conduction band offsets (CBOs) obtained by combining the ξ_M from Eq. (4) with E_g values from the experiment are given in the supplementary material (see Table S1). Moreover, another χ scheme dataset is available as Table S2. The comparison between calculated and experimental VBOs and CBOs is presented graphically in Fig. 5 for selected heterostructures. The conventional χ scheme tends to overestimate the VBOs, whereas the CBOs are underestimated in many cases. Conversely, the SSE approach does not show a pronounced trend toward under- or overestimation for VBOs and CBOs. To quantitatively evaluate the performance of different approaches, we calculated the mean error (ME), mean absolute error (MAE), mean squared error (MSE), and root mean square error (RMSE) with respect to the experiment. These statistical quantities are shown in Table I. In comparison with both statistical evaluations, the SSE approach predicts values closer to the experimental offsets than the χ scheme. These results indicated that our scheme provided a good estimation for the VBO and CBO of the selected heterostructures than conventional techniques. However, our presented method has some limitations; for example, we omitted all interface physical properties. Therefore, in some cases, a relatively large difference between the experimental and SSE offsets occurred. To resolve the errors, we required that the correction term (e.g., ±ΔV) include all inter-facial effects, such as dipole, potential, and polarization.^83–90 Unfortunately, the major problem with the correction term was that it could not use the exact numerical value based on physical phenomena in the interface, i.e., the actual interface states were too complex. On the other hand, in several cases, our calculated offsets were in good agreement with the experimental data, which may be explained by interfacial dipole vanish models (dipolar response or dipole screening mechanism).^3,91,92 Thus, our method is suitable for prediction of the offset where a pseudo-natural band alignment is realized, except for the heterostructure that has an interface dipole. In addition, the experimental results may contain errors due to Kraut’s method.^93,94 For example, VBOs for β-Ga₂O₃/SiO₂ heterostructures were quite widespread within a range of ∼1.8 eV (see the supplementary material). These accuracy problems persisted until comparing to internal photoemission measurement results considering the deposition method.^51,58,95,96 Moreover, the differences between physical and chemical states at the interface must also be considered.⁹⁷ Therefore, our method is limited in its applicability to general heterojunctions, but it is sufficient for finding band edges in isolated materials, e.g., III–V and II–VI compounds.

To understand these results obtained using our approach, we considered ξ_M. Originally, the SSE was defined by averaged IEs and averaged EAs from inorganic binary semiconductors. We simplified and expanded this idea, and by considering the case of ionic bonded materials, we identified that anion and cation atoms dominated the VBM and CBM, respectively. Furthermore, we defined the SSE of the material (ξ_M) by the constituent atoms. In the case of the A_aB_b compound, the ξ_M representation can then be written as

where $ξAa$ and $ξBb$ are the absolute SSE values of atoms A and B, respectively. The position of ξ_M is presented as a schematic in Fig. S1 (see the supplementary material). We pointed out that the energetic position of −ξ_M is located on the absolute energy scale relative to the vacuum level. In standard DFT calculations, the VBM and CBM positions are not placed at the vacuum level. Therefore, only by using ξ_M as an index can we put the VBM and CBM on an energy scale with respect to the vacuum level. The location of ξ_M halfway between the VBM and the CBM corresponded to the E_BGC that has been introduced in the literature.³² Therefore, for band edge estimation, Eqs. (2) and (3) were similar to those of the previous studies.⁹⁸ Note that the definition of ξ_M was not applicable to covalent bond materials such as Si and Ge. Furthermore, since the SSE value depends on the oxidation state of the atom, selecting an appropriate value according to the system is necessary.⁹⁹ However, our method has the advantage that natural band alignment can be determined using the ξ_M and E_g of materials.

We have examined a scheme for calculating band alignments and offsets that produces good accuracy by simple equations. The band-offset values were in good agreement with the experiment, and the band alignment reproduced the previous calculation results, which is remarkable considering its little computational efforts compared to conventional approaches. Combining our present methodology with informatics will achieve faster estimation for the natural band lineup aligned to the vacuum level of various ionic bonded materials.¹⁰⁰ Therefore, the approach is expected to be very useful for screening of p-type oxides with Hosono’s empirical rules.¹⁰¹ Furthermore, our proposed method could help predict new catalysts for water splitting and carbon dioxide reduction.

Details about the band edge position and offset values are provided in the supplementary material.

The author would like to thank Enago (www.enago.jp) for the English language review.

The data that support the findings of this study are available within the article and its supplementary material.