The ordered structure of the MgO–ZnO alloy system is a versatile tunable optical material, which is promising for diverse optoelectronic applications. However, isovalent and isostructural alloy compositions of MgO–ZnO are generally unstable at ambient conditions. Using state-of-the-art ab initio evolutionary simulations, we predict and study the properties of stable phases of MgO–ZnO. We establish the dynamical stability of the predicted crystal structures through phonon and Raman spectroscopy. Detailed analyses of two of the most stable structures reveal highly tunable properties that could be explored for photonic and optical applications.

Zinc oxide (ZnO) is a wide, direct bandgap Eg3.4 eV II–VI semiconductor that has attracted increasing interest for both fundamental research and diverse optoelectronic applications.1–3 For example, ZnO is a promising material for light-emitting and laser diodes operating in blue and ultraviolet (UV) regimes.2–4 However, low electroluminescence often diminishes the efficiency of ZnO-based devices.5,6 To realize and enhance the efficiency of ZnO-based devices, their electronic structure, through band-engineering, can be systematically tuned by doping. ZnO can be alloyed with MgO to form a ternary MgxZn1xO compound, thus, enabling bandgap engineering and luminescence in the UV regime.7 Often, the combination of these group-II oxides in alloys leads to crystal structure mismatch: undoped ZnO prefers the hexagonal wurtzite (B4) structure or the fourfold coordinated zinc-blende structure (B3), while MgO favors the cubic rock salt (B1) structure at ambient conditions.8,9 Experiments have shown that MgxZn1xO exhibits the B4 structure for high ZnO concentration,10 while preferring the B1 structure at high MgO concentration.11 At intermediate concentrations, it exhibits phase separation12 due to its compositional gradient, which often leads to a thermodynamically unstable crystal structure. In general, the isovalent and isostructural II–VI alloys are thermodynamically unstable because the mixing enthalpy in either the B1, B3, or B4 structure is always positive.12,14

The difficulty in synthesizing stable phases of MgO–ZnO alloys has hindered the adequate exploitation of their promising properties for technological applications. This prompted earlier calculations on the stability of the MgxZn1xO alloy system.12,13 For example, the work of Sanati et al.12 on various alloys of MgO–ZnO found energetic stability, under certain conditions, in the sixfold-coordinated structure for Zn concentrations below 67%. Instead of the often-used traditional doping approach of alloying ZnO with MgO to obtain the phase immiscible MgxZn1−xO alloy, we can apply high-throughput crystal structure prediction techniques to search for energetically and dynamically stable stoichiometric MgO–ZnO alloy systems. A high-throughput experiment to achieve this is, however, daunting because of the ample parameter space that needs to be explored. Computationally, such an approach is now routinely used to predict new crystal structures,15,16 including, recently, the prediction of energetically stable high pressure phase of MgxZn1xO.17 Herein, we report the electronic structure and thermodynamic stability of the stable phases of MgZnO2 crystal structure at ambient conditions obtained from the recently developed evolutionary algorithm for predicting crystal structure.16 We provide detailed physical properties and characterize electronic, mechanical, phononic, and optical spectroscopy. We extensively discuss the results in relation to experimental data of the parent systems and provide details of the crystal structure information (see supplementary material). We hope that they will motivate future experimental investigations of the electronic structure, particularly using optical measurements and photoemission, and vibrational spectroscopy by Raman measurements.

The ground state crystal structures of MgZnO2 are predicted using the ab initio evolutionary algorithm CALYPSO.16 We performed variable compositional simulations for MgnZnmOm+n systems (n, m 30) while maintaining the stoichiometric composition of MgZnO2. The structural optimizations and electronic structure calculations for the selected structures are performed using density functional theory (DFT) as implemented in the Vienna Ab Initio Simulation Package (VASP) code.18 The electronic wavefunctions of the systems are represented by a plane wave basis with an energy cutoff of 550 eV. The exchange-correlation interactions were described with the HSE06 hybrid functional19 with a 10×10×10Γ—centered k—point grid, which ensures that the energy (charge) is converged to within 104 (108) eV. We also obtained the electronic structure using the strongly constrained and appropriately normed (SCAN) semilocal density functional20 (see the supplementary material). The phononic and vibrational spectroscopy is obtained using density functional perturbation theory (DFPT)21 as implemented in VASP in conjunction with the Phonopy package.22 

In this work, we present two stable ground state crystal structures of MgZnO2. The two of the predicted structures belong to the tetragonal crystal symmetry I41/amd space group No. 141 and P4¯m2 space group No. 115, respectively (see, Fig. 1). We designate the former and latter structures as A and B, respectively. The predicted lattice constants, along with other salient parameters, are presented in Table I, and the Wyckoff positions in the crystallographic information file format are included in the supplementary material. The phase stability of the MgO–ZnO systems is studied by calculating the formation energy Eform per atom; Eform is a measure of the energetic stability of a material and for our system, it is obtained as Eform=Etoticiμi, where i denotes different types of atoms present in the unit cell of the material with concentration ci and μi is the chemical potential approximated with the standard state (bulk) energy Ei of the corresponding atom, i. The two structures exhibit negative formation energies with an average Eform2.11eV/atom; this implies an exothermic process confirming the phase and energetic stability of the crystal. To further confirm the energetic stability for the predicted structures, we calculate the enthalpy profile as ΔH=ΔG+TΔS, where H, G, S, and T is enthalpy, free energy, entropy, and temperature, respectively, for DFT, T=0K as a function of generations (see Fig. S1 in the supplementary material). All the predicted structures show a negative formation enthalpy ΔH>5.0 eV, which supports the energetic stability of the structures. The predicted ΔH is consistent with previous theoretical results for the Mg1xZnxO system;12,13 though these previously reported ΔH are negative, they are smaller than our results due to different stoichiometric compositions.

FIG. 1.

Crystal structure of MgZnO2. The crystal lattice of structure A (a) and structure B (b).

FIG. 1.

Crystal structure of MgZnO2. The crystal lattice of structure A (a) and structure B (b).

Close modal
TABLE I.

Some important predicted parameters of MgZnO2.a

Namea/cSGEgme(/)mh(/)BEGMμζΘ Dϵ(/)ϵ0(/)
6.008/8.491 141 4.62(I) 0.296/0.293 2.973/2.361 180.82 223.63 86.45 296.09 0.29 2.09 377.22 2.87/2.87 7.91/7.78 
3.245/4.509 115 3.70(D) 0.304/0.312 1.772/2.728 140.13 116.26 42.80 197.20 0.36 3.27 502.46 2.83/2.83 6.80/5.40 
Namea/cSGEgme(/)mh(/)BEGMμζΘ Dϵ(/)ϵ0(/)
6.008/8.491 141 4.62(I) 0.296/0.293 2.973/2.361 180.82 223.63 86.45 296.09 0.29 2.09 377.22 2.87/2.87 7.91/7.78 
3.245/4.509 115 3.70(D) 0.304/0.312 1.772/2.728 140.13 116.26 42.80 197.20 0.36 3.27 502.46 2.83/2.83 6.80/5.40 
a

The lattice constant includes a (Å), space group SG with crystal symmetry I41/amd (P4¯m2) for structure A (B); the electronic properties include the energy bandgap Eg (eV), electron me and hole Mh effective masses in units of free electron mass m0, where (⊥) denotes parallel (perpendicular) direction, and I (D) is indirect (direct) energy bandgap; the mechanical properties include the bulk modulus B (GPa), Young’s modulus E (GPa), shear modulus G (GPa), P-wave modulus M (GPa), Poisson’s ratio μ, and bulk/shear ratio ζ; the dielectric properties include the high photon-energy ε and static ε0 dielectric constants; and the Debye temperature Θ D (K) representing the thermal property.

While the energetic stability is often a confirmation that a crystal can be experimentally synthesized, in most cases, basing the overall stability on the energy alone may be misleading. In the absence of any phase transition, a crystal may have negative formation energy yet thermodynamically unstable. To ascertain the dynamical stability of the predicted structures, we performed phonon and vibrational spectroscopy calculations (Fig. 2). We also show in Fig. 2 the Raman and Infrared spectroscopy. Both structures show positive definite phonon frequency along the various high symmetry point in the first Brillouin zone, confirming the dynamical stability. From the phonon atomic projected density of states, the acoustic branch is dominated by Zn atoms, while the optical branch is composed of a strong hybridization between the O and Mg atoms. This trend is consistent with the domination of the frequency scale of the acoustic (optical) phonon modes by atoms with larger (smaller) masses similar to a diatomic linear chain model.23 

FIG. 2.

Phonon and vibrational spectroscopy of MgZnO2. Phonon spectra and the corresponding vibrational density of states, the Raman, and Infrared spectra for structure A (a) and structure B (b).

FIG. 2.

Phonon and vibrational spectroscopy of MgZnO2. Phonon spectra and the corresponding vibrational density of states, the Raman, and Infrared spectra for structure A (a) and structure B (b).

Close modal

The predicted crystal structures belong to the tetragonal crystal symmetry, albeit with a different point group symmetry leading to a different irreducible representation of the vibrational modes. Crystal structure A belongs to the point group D4h (4/mmm) with 16 atoms in the primitive cell. The irreducible representation of the corresponding 48 normal vibrational modes at the Brillouin zone center (Γ-point) is Γ=2A1g(R)+3A1u+A2g+6A2u(IR)+2B1g(R)+3B1u+B2g(R)+6B2u+9Eu(IR)+3Eg(R), where IR and R are the infrared and Raman active modes. A1u, A2g, B1u, and B2u are the silent modes. The R-active Eg, B1g, and the A1g modes are mainly due to the symmetric stretching vibration, the symmetric bending vibration, and the antisymmetric bending vibration of the O–Mg/Zn–O, respectively. The phonon spectra, the phonon total, and the atomic projected density of states along with the Raman and Infrared spectra are presented in Fig. 2(a). (Further analysis of the irreducible representation of the phonon modes are shown in Table S1 of the supplementary material.) We hope the presented phonon, vibrational, Raman, and Infrared spectroscopy will guide experimental and other analyses. For the crystal structure B with point group symmetry D2d(42m), the primitive cell contains four atoms; as such, there will be a total of twelve normal vibrational modes at the Γ point [Fig. 2(b)]. The irreducible representation of the vibrational modes at the Γ point could be expressed as Γ=A1(R)+3B2(IR+R)+4E(IR+R). One low-energy B2 and one double E modes are the usual acoustic modes characterized by the transverse acoustic (TA), longitudinal acoustic (LA), and the out-of-plane transverse acoustic or flexural mode (ZA); the remaining belong to the optical mode. Using group theory analysis24 in conjunction with the Raman and Infrared spectra [Fig. 2(b)], we assign the feature at 178.70 cm1 to the degenerate IR- and R-active E mode and the structure at 326.28 cm1 to the IR- and R-active B2 mode. We predict a degenerate IR- and R-active E mode at 413.94 cm1, followed by IR- and R-active E mode at 480.60 cm1. We predict an R-active A1 mode though with small activity around 531.00 cm1 before a degenerate IR- and R-active E mode at 548.81 cm1. (Further details of the irreducible representation are presented in Table S2 in the supplementary material.) We note that the predicted features in both the vibrational and Raman spectra depict more of ZnO character than MgO. For example, room temperature Raman spectroscopy experiments reported the R-active modes E583 cm1 and A1574 cm1.25 

The Debye temperature ΘD is an important thermal property of a material. It quantifies the temperature below which modes start to freeze out and above which, the modes gain enough energy to be excited.26 The Debye temperature can be obtained as ΘD3=13i1ΘΓ,i3, where ΘΓ,i=ΩΓ,i/κB and ΘΓ,i corresponds to the frequency at the zone boundary of the ith acoustic mode.27 The predicted ΘD is presented in Table I, which seem to be closer to that of ZnO with ΘD440 K than to MgO with ΘD948 K as reported by experiments.9,25

The elastic properties not only provide the information on the mechanical stability of material but also serve as a guide on the potential device applications of a given material. In this regard, we have calculated the elastic properties of the studied materials using the finite difference method. We used a plane wave energy cutoff of 700 eV to ensure the adequate convergence of the stress tensor. The elastic tensor matrix is determined by performing six finite distortions of the lattice and obtaining the elastic constants from the strain–stress relationship.28,29 The predicted elastic tensor matrices are presented in Table II. Using the calculated elastic tensor matrix, we predict the elastic constants presented in Table I. Aside from the isotropic elastic behavior [the calculated elastic anisotropy ratio A=2C44/(C11C12)1.0], the two studied structures show markedly different elastic behavior. Structure A (B) could withstand linear compressibility of 1.81 (1.65 TPa1). The predicted Bulk modulus of structure A is 41 GPa larger than that of structure B. The same trend is observed for Young’s modulus E, shear modulus G, P-wave modulus M, and the Poisson ratio μ, which is larger by E107.37 GPa, G43.65 GPa, M98.89 GPa, and μ0.07. Though structure A exhibits large elastic constants, structure B is, however, more ductile from the bulk/shear ratio value (Table I).

TABLE II.

Calculated elastic constant tensor Cij for MgZnO2. Only half side is shown by symmetry.

Structure AStructure B
Cij123456123456
312.05 88.39 142.04 0.36 0.0 0.0 225.77 85.46 123.68 0.0 2.32 0.0 
 313.29 143.02 0.38 0.0 0.0  227.90 124.45 0.0 2.41 0.0 
  255.17 0.45 0.0 0.0   151.77 0.0 2.40 0.0 
   114.33 0.0 0.0    70.48 0.0 0.06 
    113.96 0.0     70.61 0.0 
     58.38      20.55 
Structure AStructure B
Cij123456123456
312.05 88.39 142.04 0.36 0.0 0.0 225.77 85.46 123.68 0.0 2.32 0.0 
 313.29 143.02 0.38 0.0 0.0  227.90 124.45 0.0 2.41 0.0 
  255.17 0.45 0.0 0.0   151.77 0.0 2.40 0.0 
   114.33 0.0 0.0    70.48 0.0 0.06 
    113.96 0.0     70.61 0.0 
     58.38      20.55 

The description of the electronic properties of the two forms of MgZnO2 studied herein as obtained using the HSE06 hybrid functional is provided by the band structure, total, and atomic projected density of states in Fig. 3. We predict an indirect energy bandgap Eg4.62 eV along the ΓN point of the high symmetry point of the Brillouin zone for structure A. We note that the direct energy bandgap is 4.52 eV, which is only 0.10 eV smaller than the largest bandgap; as such, the indirect and direct energy bandgaps could be said to be almost degenerate. Structure B is predicted to be a direct bandgap semiconductor with Eg3.70 eV at the Γ-point. The predicated bandgaps are listed in Table I. We note that the energy bandgap predicted with the semilocal SCAN (see Fig. S2 in the supplementary material) and the Perdew–Burke–Ernzerhof generalized gradient approximation (PBE-GGA)30 is significantly smaller in both materials. For example, for structure A, SCAN and PBE-GGA functional led to Eg3.32 and 2.73 eV, respectively. We obtained Eg2.43 and 1.98 eV using SCAN and PBE-GGA functional, respectively, for structure B.

FIG. 3.

Electronic properties of MgZnO2. Electronic band structure and the corresponding electron density of states for structure A (a) and structure B (b). The horizontal dashed line is the Fermi energy, which is set at zero of the energy scale.

FIG. 3.

Electronic properties of MgZnO2. Electronic band structure and the corresponding electron density of states for structure A (a) and structure B (b). The horizontal dashed line is the Fermi energy, which is set at zero of the energy scale.

Close modal

The two studied structures share in common some of the features of ZnO. For example, the predicted energy bandgap is in the same range as the 3.4 eV reported for bulk ZnO by several experiments and computations.9,25,31–33 Aside from the difference in the energy bandgap of the two structures, the character of the states in the proximity of the Fermi energy are similar for the two structures. From the projected density of states, the lower-lying group in the valence states consists of highly hybridized O (from s- and p-states), while the upper group is mostly dominated by Zn (d-states), O (p-states), and tiny contribution from Mg (s- and p-states). The uppermost conduction states are dominated by a strong hybridization between Zn (s-states) and O (s- and p-states). This is evident from the free electron-like (s-like) states at the conduction band minimum (CBM) located at the Γ-point. The total width of the upper valence bands within 10.0 eV from the valence band maximum (VBM) is 6.30 (5.28 eV) for structure A (B). This value is in good agreement with 6.50 eV for MgO and approximately 1–2 eV smaller than 7.2 eV for ZnO reported by experiments.31,34 The total width of the entire valence states is 19.52 (18.25 eV) for structure A (B). These values are within the experimental uncertainty of the total valence width of ZnO.25 We note that the valence bandwidth obtained with the scan functional decreased by 5.2%. With the PBE-GGA functional, the decrease is more significant 14.07%.

The effective mass is a measure of the curvature of the calculated band structure, and it provides insight into the transport properties of materials. To this end, we calculate the carrier band effective mass of the two studied structures. The band effective mass mb is obtained from the band structure by fitting a parabola Ek=22m0kTAk to the states around the band extremum (CBM and VBM) of Fig. 3, where k=(,) is the parallel and perpendicular direction of the k-point measured from the band extremum and m0 is the free electron mass. The predicted electron and hole effective masses along the () direction are listed in Table I. The effective mass decreased by 6.50 and 15.20%, respectively, with the SCAN and PBE-GGA functionals. The large hole effective mass, especially for structure A, is due to the rather flatbands around the top of the valence states. We note that the predicted me/ straddle that of ZnO 0.275 and MgO 0.370 m0.9,25,35,36 While the predicted mh is in good agreement with 0.60 m0 reported for ZnO,25 the calculated mh is in better agreement with the 1–2 m0 reported for MgO.36 

To gain further insights into the physical properties of the predicted structures, we calculate the electrodynamical response as a function of photon energy (ω). The optical spectroscopy provides a detailed and direct information on charge dynamics and “bulk” electronic structure of materials since it probes both free carriers and interband excitation. The optical response is determined using the spectra from our electronic structure calculations obtained with HSE06 hybrid functional without any shift or scissor operator. We used a small broadening parameter 104 eV as such, our optical spectra show a lot of features. Using random phase approximation, we compute the real part of the complex dynamical dielectric function

ϵαβ(2)=4π2e2Ωlimq1q2c,v,kΞ2wkδ(εckεvkω),

where Ξ=uck+eαq|uvkuck+eβq|uvk, ω is the photon frequency, α (β) is band indices, v and c denote valence and conduction band states, and uck is the cell periodicity at the k-point over volume Ω. The dispersive (real) part of the dynamical dielectric function ϵ(1)(ω) is then obtained using the appropriated Kramers–Kronig transformation

ϵαβ(1)=1+2πP0ϵαβ(2)(ω)ωdωω2ω2+iη,

where P is the principal value.21 The calculated optical spectra are presented in Fig. 4. The onset of absorption that is determined by the absorption edge depicts the characteristic optical bandgap. This corresponds to the lowest sharp structure around 4.53 (3.12) eV for structure A (B) due to the direct optical excitation.

FIG. 4.

Optical properties of MgZnO2. Calculated imaginary part ϵ2(ω) (a) and the real part ϵ1(ω) (b) of the dielectric function of MgZnO2. Blue (red) color denotes structure A (B).

FIG. 4.

Optical properties of MgZnO2. Calculated imaginary part ϵ2(ω) (a) and the real part ϵ1(ω) (b) of the dielectric function of MgZnO2. Blue (red) color denotes structure A (B).

Close modal

The predicted optical bandgap is within the numerical uncertainty of the ones obtained from our band structure calculations (Table I). The optical spectra for both materials share almost the same characteristic features. In the ϵ2(ω) spectra [Fig. 4(a)], we found a sharp structure around 6.50 (5.06), followed by a set of close sharp features and shoulders that transcend in both structures up to 17.5 eV, before the spectra systematically decay toward zero photon energy.

In Fig. 4(b), we present the calculated dispersive parts of the dynamical dielectric function ϵ(1)(ω). The main features in both structures are sharp structures around 4.95 (3.50) eV for structure A (B), series of other features between 5.20 eV and 15.10 eV, characterized by a steep decrease that started around 7.80 eV. Then, ϵ(1)(ω) for both structures becomes negative with a minimum at 15.10 eV. The optical spectra obtained using the electronic structure spectra from the SCAN semilocal functional though exhibits the same trend, the critical energies are significantly different (see Fig. S3 of the supplementary material). We attribute this to the nonlocal effects that are inherent in the HSE06 hybrid functional as evident from the sharper structures in the optical spectra. We further determined the high photon energy ϵ and the static photon energy ϵ0 limits of the dielectric constants. ϵ is obtained as the zero-photon energy limit of the ϵ(1)(ω) spectra and the obtained values are listed in Table I. To obtain ϵ0, we note that ϵ0=ϵ+ϵion, where ϵion is the static photon energy (ω=0) contribution of the ionic lattice dynamics to the dielectric constant. We obtained from our phonon calculations ϵion(/)5.04/4.91 (3.97/2.57) for structure A (B). Using the obtained ϵion, we predict the static photon energy of the dielectric constant summarized in Table I. The predicted optical behavior of structure A is generally isotropic while that of structure B could be said to be slightly anisotropic, mainly from the dynamics of ϵion(/). We note that the ϵ(/) obtained using the eigenstates from the SCAN functional as the input in our optical calculations is 3.43/3.43 (3.34/3.31) for structure A (B), which seem to be closer to the values reported for ZnO.25 Overall, compared to the parent materials, the predicted limits of the dielectric constants straddle those of ZnO and MgO. For example, experiments report ZnO has ϵ(/)3.75/3.70 and ϵ0(/)8.75/7.80, while MgO exhibits ϵ(/)2.94 and ϵ0(/)9.83.25 

We report detailed ab initio results of the structural, phononic, electronic, transport, and optical properties of newly predicted crystal structures of the MgZnO2 alloy system. Through the calculations of formation energies and the phonon and vibrational spectra, including the Raman and Infrared spectroscopy, we establish the dynamical stability of the predicted structures. We demonstrate highly tunable electronic and optical properties of two of the predicted structures with optoelectronic features that make them promising for application as a photonic material. We provide detailed crystal structure information to guide further experimental studies.

See the supplementary material for the predicted crystal lattice parameters, the irreducible representation of the phonon active modes at the Γ-point, and the predicted electronic and optical properties obtained with the SCAN functional.

Supercomputer and computational resources were provided by the Lehigh University (LU) HPC Center. K.A. acknowledges support of the LU Physics REU program (No. PHY-1852010) and C.E.E. acknowledges the LU start-up and summer research fellowship.

1.
X.
Yu
,
T. J.
Marks
, and
A.
Facchetti
,
Nat. Mater.
15
(
4
),
383
396
(
2016
);
[PubMed]
M.
Vijaya
,
Piezoelectric Materials and Devices: Applications in Engineering and Medical Sciences
(
CRC Press
,
2016
);
S.
Pearton
,
D.
Norton
,
Y.
Heo
,
L.
Tien
,
M.
Ivill
,
Y.
Li
,
B.
Kang
,
F.
Ren
,
J.
Kelly
, and
A.
Hebard
,
J. Electron. Mater.
35
(5),
862
868
(
2006
).
2.
M. H.
Huang
,
S.
Mao
,
H.
Feick
,
H.
Yan
,
Y.
Wu
,
H.
Kind
,
E.
Weber
,
R.
Russo
, and
P.
Yang
,
Science
292
(
5523
),
1897
1899
(
2001
).
3.
S. C.
Su
,
Y. M.
Lu
,
Z. Z.
Zhang
,
C. X.
Shan
,
B. H.
Li
,
D. Z.
Shen
,
B.
Yao
,
J. Y.
Zhang
,
D. X.
Zhao
, and
X. W.
Fan
,
Appl. Phys. Lett.
93
(
8
),
082108
(
2008
).
4.
D.
Thomas
,
J. Phys. Chem. Solids
15
(
1-2
),
86
96
(
1960
).
5.
A.
Tsukazaki
,
A.
Ohtomo
,
T.
Onuma
,
M.
Ohtani
,
T.
Makino
,
M.
Sumiya
,
K.
Ohtani
,
S. F.
Chichibu
,
S.
Fuke
,
Y.
Segawa
,
H.
Ohno
,
H.
Koinuma
, and
M.
Kawasaki
,
Nat. Mater.
4
,
42
(
2004
).
6.
J.-H.
Lim
,
C.-K.
Kang
,
K.-K.
Kim
,
I.-K.
Park
,
D.-K.
Hwang
, and
S.-J.
Park
,
Adv. Mater.
18
(
20
),
2720
2724
(
2006
).
7.
T. A.
Wassner
,
B.
Laumer
,
S.
Maier
,
A.
Laufer
,
B. K.
Meyer
,
M.
Stutzmann
, and
M.
Eickhoff
,
J. Appl. Phys.
105
(
2
),
023505
(
2009
);
F. K.
Shan
,
B. I.
Kim
,
G. X.
Liu
,
Z. F.
Liu
,
J. Y.
Sohn
,
W. J.
Lee
,
B. C.
Shin
, and
Y. S.
Yu
,
J. Appl. Phys.
95
(
9
),
4772
4776
(
2004
);
D.
Jiang
,
C.
Shan
,
J.
Zhang
,
Y.
Lu
,
B.
Yao
,
D.
Zhao
,
Z.
Zhang
,
X.
Fan
, and
D.
Shen
,
Cryst. Growth Des.
9
(
1
),
454
456
(
2009
).
8.
D. C.
Reynolds
,
D. C.
Look
, and
B.
Jogai
,
Solid State Commun.
99
(
12
),
873
875
(
1996
).
9.
O.
Madelung
,
U.
Rössler
, and
M.
Schulz
,
Landolt-Brnstein Database, Vol. 10
(
1999
), see http://www.springermaterials.com.
10.
Y.
Ohno
,
D. K.
Young
,
B.
Beschoten
,
F.
Matsukura
,
H.
Ohno
, and
D. D.
Awschalom
,
Nature
402
(
6763
),
790
792
(
1999
).
11.
J.
Narayan
,
A.
Sharma
,
A.
Kvit
,
C.
Jin
,
J.
Muth
, and
O.
Holland
,
Solid State Commun.
121
(
1
),
9
13
(
2001
).
12.
M.
Sanati
,
G. L. W.
Hart
, and
A.
Zunger
,
Phys. Rev. B
68
,
155210
(
2003
).
13.
X. F.
Fan
,
H. D.
Sun
,
Z. X.
Shen
,
J.-L.
Kuo
, and
Y. M.
Lu
,
J. Phys. Condens. Matter
20
(
23
),
235221
(
2008
);
[PubMed]
Y.-S.
Kim
,
E. C.
Lee
, and
K. J.
Chang
,
J. Kor. Phys. Soc.
39
,
92
(
2001
).
14.
S.-H.
Wei
and
A.
Zunger
,
Phys. Rev. B
37
,
8958
8981
(
1988
);
A.
Seko
,
F.
Oba
,
A.
Kuwabara
, and
I.
Tanaka
,
Phys. Rev. B
72
,
024107
(
2005
);
A.
Malashevich
and
D.
Vanderbilt
,
Phys. Rev. B
75
,
045106
(
2007
).
15.
R.
Gautier
,
X.
Zhang
,
L.
Hu
,
L.
Yu
,
Y.
Lin
,
T. O. L.
Sunde
,
D.
Chon
,
K. R.
Poeppelmeier
, and
A.
Zunger
,
Nat. Chem.
7
,
308
(
2015
);
[PubMed]
A. R.
Oganov
and
C. W.
Glass
,
J. Chem. Phys.
124
(
24
),
244704
(
2006
).
[PubMed]
16.
Y.
Wang
,
J.
Lv
,
L.
Zhu
, and
Y.
Ma
,
Phys. Rev. B
82
,
094116
(
2010
);
Y.
Wang
,
J.
Lv
,
L.
Zhu
, and
Y.
Ma
,
Comput. Phys. Commun.
183
(
10
), pp.
2063
2070
(
2012
).
17.
F.
Tian
,
D.
Duan
,
D.
Li
,
C.
Chen
,
X.
Sha
,
Z.
Zhao
,
B.
Liu
, and
T.
Cui
,
Sci. Rep.
4
(
1
),
5759
(
2014
).
18.
G.
Kresse
and
J.
Furthmüller
,
Comput. Mater. Sci.
6
(
1
),
15
50
(
1996
).
19.
A. V.
Krukau
,
O. A.
Vydrov
,
A. F.
Izmaylov
, and
G. E.
Scuseria
,
J. Chem. Phys.
125
(
22
),
224106
(
2006
);
[PubMed]
J.
Heyd
,
G. E.
Scuseria
, and
M.
Ernzerhof
,
J. Chem. Phys.
124
(
21
),
219906
(
2006
);
J.
Heyd
,
G. E.
Scuseria
, and
M.
Ernzerhof
,
J. Chem. Phys.
118
(
18
),
8207
8215
(
2003
).
20.
J.
Sun
,
A.
Ruzsinszky
, and
J. P.
Perdew
,
Phys. Rev. Lett.
115
,
036402
(
2015
).
21.
M.
Gajdoš
,
K.
Hummer
,
G.
Kresse
,
J.
Furthmüller
, and
F.
Bechstedt
,
Phys. Rev. B
73
,
045112
(
2006
).
22.
A.
Togo
,
F.
Oba
, and
I.
Tanaka
,
Phys. Rev. B
78
,
134106
(
2008
).
23.
C.
Ekuma
,
S.
Najmaei
, and
M.
Dubey
,
Mater. Today Commun.
19
,
383
392
(
2019
).
24.
E.
Kroumova
,
M.
Aroyo
,
J.
Perez-Mato
,
A.
Kirov
,
C.
Capillas
,
S.
Ivantchev
, and
H.
Wondratschek
,
Phase Transit.
76
(
1-2
),
155
170
(
2003
).
25.
O.
Madelung
,
U.
Rössler
, and
M. E.
Schulz
, Numerical Data and Functional Relationships in Science and Technology, Landolt-Börnstein, New Series, Group III (Springer, 2006), Vols. 17a, and 22a.
26.
N.
Ashcroft
and
N.
Mermin
,
Solid State Physics
, HRW International ed. (
Holt
,
Rinehart and Winston
,
1976
). ISBN 9780030839931; T. Nakashima and Y. Umakoshi, Philos. Mag. Lett. 66(6), 317–321 (1992).
27.
D. T.
Morelli
and
J. P.
Heremans
,
Appl. Phys. Lett.
81
(
27
),
5126
5128
(
2002
).
28.
Y.
Le Page
and
P.
Saxe
,
Phys. Rev. B
65
,
104104
(
2002
).
29.
X.
Wu
,
D.
Vanderbilt
, and
D. R.
Hamann
,
Phys. Rev. B
72
,
035105
(
2005
).
30.
J. P.
Perdew
,
K.
Burke
, and
M.
Ernzerhof
,
Phys. Rev. Lett.
77
,
3865
3868
(
1996
).
31.
U.
Özgür
,
Y. I.
Alivov
,
C.
Liu
,
A.
Teke
,
M. A.
Reshchikov
,
S.
Doǧan
,
V.
Avrutin
,
S.-J.
Cho
, and
H.
Morkoç
,
J. Appl. Phys.
98
(
4
),
041301
(
2005
).
32.
L.
Franklin
,
C. E.
Ekuma
,
G. L.
Zhao
, and
D.
Bagayoko
,
J. Phys. Chem. Solids
74
(
5
),
729
(
2013
).
33.
J. F.
Muth
,
R. M.
Kolbas
,
A. K.
Sharma
,
S.
Oktyabrsky
, and
J.
Narayan
,
J. Appl. Phys.
85
(
11
),
7884
7887
(
1999
).
34.
S.
Kowalczyk
,
F.
McFeely
,
L.
Ley
,
V.
Gritsyna
, and
D.
Shirley
,
Solid State Commun.
23
(
3
),
161
169
(
1977
).
35.
R. G.
Evans
,
E. L.
Clark
,
R. T.
Eagleton
,
A. M.
Dunne
,
R. D.
Edwards
,
W. J.
Garbett
,
T. J.
Goldsack
,
S.
James
,
C. C.
Smith
,
B. R.
Thomas
,
R.
Clarke
,
D. J.
Neely
, and
S. J.
Rose
,
Appl. Phys. Lett.
86
(
19
),
191505
(
2005
).
36.
Q.
Yan
,
P.
Rinke
,
M.
Winkelnkemper
,
A.
Qteish
,
D.
Bimberg
,
M.
Scheffler
, and
C. G.
Van de Walle
,
Appl. Phys. Lett.
101
(
15
),
152105
(
2012
).

Supplementary Material