1.1

Methods

The present work utilizes several electronic structure tools, which are implemented within the projector augmented wave (PAW) framework in the VASP code: (i) Density functional theory (DFT) calculations using the generalized gradient approximation (GGA) [13][14]. (ii) Total energy calculations within the random phase approximation (RPA) based on the adiabatic connection dissipation fluctuation theorem [15][16]. (iii) GW quasiparticle energy calculations [17][18]. The polymorph energies are determined for the rocksalt (RS), zincblende (ZB), and wurtzite (WZ) structures using the RPA, which has recently been shown to provide an improved description of polymorphic energy ordering compared to standard DFT [19]. Since the ZB and WZ polymorphs have rather similar energies and properties, the RS and ZB structures are used as prototypes for the octahedral and tetrahedral coordination symmetries. The respective band-structures are calculated in GW, following the same general recipe that has been published before [20][21], i.e., maintaining the DFT wavefunctions, but iterating the eigenenergies to self-consistency. Brillouin zone sampling was performed on a Γ centered 12×12×12 k-mesh. In order to improve the d-orbital energies, which tend to be overestimated in GW [22], an additional on-site potential for Zn-d (V_d = −1.5 eV) and Cd-d (V_d = −0.5 eV) is applied in the GW calculation, as described previously [23]. For both the RPA and the GW calculations, the GGA+U [24] method is used (U_Zn = 6 eV, U_Cd = 5 eV) to generate the wavefunctions on which the total and quasiparticle energies are evaluated. Following the rationale for combining DFT and GW results as recently established for defect [25] and surface calculations [26], the potential step at the nonpolar (100) and (110) surfaces of the RS and ZB structures is calculated in DFT with the same functional as used for generating the wavefunctions for the GW calculations, so that the GW quasiparticle energies can be directly referred to the vacuum level. The typical overestimation of the lattice constant by about 1% in GGA has been accounted for by a scaling of the lattice volume prior to the calculations of the surface potential step and the quasi-particle energies. Finally, the alloy mixing enthalpies are determined by using the polymorph energies from RPA, and random alloy models, calculated in GGA using the special quasi-random structures (SQS) of Ref. [27] for both the RS and ZB lattices.

2.1

Polymorphism

The polymorphic energy ordering of the II chalcogenides has been calculated on the DFT level in the previous literature [7][9]. In a recent study [19], we have shown that the RPA is able to overcome the difficulty of various DFT functionals to describe the energy ordering in the transition metal oxide MnO [28]. Less dramatic, but nevertheless significant differences between DFT and RPA were also observed for MgO and ZnO. Thus, the energy ordering of the II chalcogenides is reevaluated in the present work using the RPA, and the results are summarized in table 1. Besides an increase of the bulk modulus by about 20% and a reduction of the lattice volume by about 3% (which improves the agreement of the lattice constants with experiment), the most important change due to the RPA is the relative stabilization of the octahedral (RS) over the tetrahedral (ZB and WZ) coordination: The energy difference [E_RPA(RS) – E_RPA(ZB)] – [E_DFT(RS) – E_DFT(ZB)] ranges from −31 meV for ZnS to as much as −139 meV for CdS. The energy difference between the two tetrahedral polymorphs (ZB and WZ) is comparatively small, both in DFT and RPA, as expected. Notably, the RPA predicts that the ZB structure of CdS is slightly lower in energy than the WZ phase, which is generally considered to be the ground state [29], even though the zinc-blende modification of CdS is well known. However, the energy difference of 8 meV/fu is small enough to be considered essentially degenerate within the limits of the calculations, which do not include zero point motion and finite temperature effects. Note that MgO is not stable in the WZ structure [11], and is therefore omitted in table 1. Due to the small energy difference between ZB and WZ, and the generally rather similar properties of the II chalcogenides in these structures, the ZB polymorph is taken as a prototype for the tetrahedral modification in the following. Table 1 reveals a number of polymorph/substrate combinations that are lattice matched and could interesting for the study of polymorphism in epitaxial growth: ZB-MgO/CdO, ZB-ZnO/CdO, ZB-CdO/MgS, RS-ZnS/MgS (note, however, the exceptionally high energy of RS-ZnS), RS-CdS/ZnS.

Table 1.

Polymorphism of group II chalcogenides. The DFT (GGA) and RPA energies are given with respect to the ground state structure (bold). The bulk moduli (B0) and lattice volumes V per formula unit (fu) are given are also given for the RPA.

2.2

Band-structure and band lineup

With exception of MgS, the band-structure properties of the II chalcogenides in their ground state structure are fairly well established [30]. The direct band gap of CdO has been studied in detail recently [31], but the indirect band gap is not well known experimentally. Theoretically, state-of-the-art GW band-structure calculations have been reported mostly for most of the II chalcogenides [12] [22], and some band-structure predictions exist for polymorphism [7][11][32]. The purpose of the present work is to provide a comprehensive comparison between oxides and sulfides, and between the RS and ZB polymorphs, based on a consistent set of data calculated by the same approach. In addition, the combination of GW bulk and DFT surface calculations allows for a quantitative prediction of the band-lineups with respect to the vacuum level at the non-polar (100) and (110) surfaces of the RS and ZB structures, respectively. The results are summarized in table 2, and the band-structures are plotted in figure 1. We observe the following trends: (i) ZB compounds are all direct semiconductors with a band gap at the Brillouin zone center. (ii) In the RS structure, the IIb chalcogenides have the valence band maximum (VBM) at the L point, which can be traced back to the fact that in the octahedral symmetry, the p-d repulsion is forbidden at the Γ point, but allowed at the L point [33]. For example, the calculated direct gap of CdO is 2.0 eV, whereas the indirect gap is only 0.9 eV (see table 2). Note that the Zn and Cd sulfides in the RS structure have valence band energies inside the Brillouin zone that are essentially degenerate with the L point. (iii) The RS sulfides have indirect gaps with the conduction band minimum (CBM) at the X point. (iv) The IPs of the ZB structures are about 1 eV larger that those in the RS structure. (v) The IPs do not show the same chemical trends with the anion type as observed in band-offset calculations [34]: The sulfides have similar or even larger IPs than the oxides, i.e., a lower VBM energy. The present predictions agree rather well with experimental values for ZnO, ZnS, and CdS, where the IPs have been determined as 7.8, 7.5, and 7.3 [35]. However, one should consider that the IP are often sensitive to the surface orientation and termination [26][36]. The difference in the trends between band-offsets and IP is notable, in particular since the band offsets were constructed to eliminate strain effects that inevitably occur at the coherent epitaxial interface in a slab supercell (in contrast, there is no strain at the free surface used to calculate the IP). However, the construction of such a ”natural band offset” [37] relies of the definition of a potential reference, which introduces some ambiguity [38][39]. Further insight, as to whether the band offset in lattice-mismatched non-epitaxial interfaces is controlled by the natural band-offsets or the difference in IPs, is desirable.

Table 2.

The predicted ionization potentials (IP), electron affinities (EA), band-gaps (Eg), and surface energies (Esurf) of II chalcogenides. The IP, EA, and Esurf are for the (100) and (110) surfaces of the RS and ZB structures, respectively, and the direct (d) or indirect (i) character of the band gap is indicated.

Figure 1.

Band-structure energies calculated in the GW approximation for the group II oxides and sulfides. The energy zero is the vacuum level at the (100) surface of the rocksalt structure and at the (110) surface of the zincblende structure.

Shown in Table 2 are also the calculated surface energies (E_surf) for the (100) and (110) surfaces of the Rs and ZB structures, respectively. The chemical trends are that oxides have a larger E_surf than sulfides, and that E_surf decreases along the cation series Mg → Zn → Cd. The structural trend is that RS surfaces have a lower energy than ZB surfaces.

2.3

Alloy mixing enthalpies

In many alloy systems, the mixing enthalpies are well described to quadratic order by a single alloy interaction parameter Ω [40], which implies that the mixing enthalpy of an A_1−xB_x alloy, ΔH_m(x) = E(A_1−xB_x) − (1−x)·E(A) − x·E(B). is symmetric with respect to x = 0.5. Including the possibility of an asymmetric mixing enthalpy, and considering that we are dealing here with alloys between heterostructural end compounds, the mixing enthalpy can be expressed up to 3^rd order as

Here, E₀ and E₁ are the energies of the end compounds relative to the ground state structure, Ω is the regular alloy interaction parameter, and Ω’ is the 3^rd order parameter describing the asymmetry of ΔH_m. Using the DFT (GGA) energies of SQS alloy models for x = 0.25, 0.5, 0.75, and taking E₀ and E₁ as the RPA energy differences between the respective RS and ZB structures (see table 1), the Ω and Ω’ parameters were fitted as shown in table 3. The full ΔH_m(x) diagrams are shown in figure 2.

Figure 2.

Calculated alloy mixing enthalpies ΔH_m for both the octahedral RS and the tetrahedral ZB lattices. The data points are obtained from the RPA polymorph energies (see table 1) and the DFT mixing enthalpies within the RS/ZB lattices for SQS cells at compositions x = 0.25, 0.5, and 0.75. The lines show the results of the fitting according to eq. (1). Note the difference energy scale for oxide-sulfide alloys.

Table 3.

The alloy interaction parameters Ω and Ω’, fitted according to eq. (1), for both the RS and ZB structures.

A number of observations can be made about the alloy interaction parameters: (i) In the RS structure, the mixing enthalpies are larger (Ω) and more asymmetric (Ω’) than in the ZB structure. This behavior is likely related to the much higher density of the RS structure (see table 1). In the more open structure of the ZB lattice, strain due to atomic size mismatch is more easily relieved at a lower energy cost. In cases with strong asymmetries, the signs of Ω’ are positive, indicating that it costs more energy to substitute a large atom for a small atom than vice versa. This behavior can be expected from the asymmetry of atomic pair potentials, where the energy rises more quickly for a distance reduction relative to the energy minimum than for an increase [41]. (ii) The alloy interaction parameters are larger in the mixed anion alloys than in the mixed cation alloys, reflecting the larger size mismatch between O and S, as compared to that between the Mg, Zn, and Cd (cf. the lattice volumes in table 1). (iii) The small values of Ω for the Mg/Zn alloys reflect their small size mismatch. In these cases, the mixing enthalpy is dominated by the polymorphic energy difference. As seen in figure 2, composition dependent structure transitions are expected for Mg_1−xZn_xO, Zn_1−xCd_xO, Mg_1−xZn_xS, Mg_1−xCd_xS, MgO_1−xS_x, and CdO_1−xSx. For Mg_1−xZn_xO and Zn_1−xCd_xO, these transitions have been observed experimentally at similar compositions [4][42]. Note that the DFT values for the polymorph energies of CdO and ZnO (table 1) would imply a minute stability range of the RS structure in Zn_1−xCd_xO alloys in the Cd rich regime, whereas experiments indicate a stability range 0.7 < x ≤ 1.0 that is even wider than predicted from the RPA results. This observation further corroborates the finding that DFT (GGA) overly favors the tetrahedral symmetry. An interesting situation is found for MgO_1−xS_x, where the end compounds have the RS structure, but due to the smaller value of Ω_ZB, the ZB structure is predicted to be stable within the composition interval 0.25 < x < 0.73.