Hybrid density functional theory has become a standard method for calculations of defects in semiconductors. The majority of work in this field is done using hybrid functionals tuned to reproduce the experimental bandgap of the host material. This approach usually yields results in reasonable agreement with the experiment. Alternatively, hybrid functional can be tuned to fulfill the generalized Koopmans' condition for defect orbitals, which cancels self-interaction energy and restores the linear behavior of energy with respect to electron occupation. Here, we investigate the methods of hybrid functional tuning, which both satisfy the generalized Koopmans' condition and reproduce the experimental bandgap, using one of the most well-studied defects in GaN, carbon acceptor. We test different charged defect correction schemes, the influence of Ga3d-electrons, and compare the results with accurate photoluminescence measurements. We find that using different charged defect correction methods can lead to substantially different hybrid functional parametrizations. However, the calculated optical properties of the carbon acceptor are found to be weakly dependent on specific parameters.

Gallium nitride is a wide bandgap semiconductor with numerous applications in light-emitting devices,1 lasers,2 solar cells,3 and high-power electronics.4,5 A wide variety of devices depend on the ability to control and take advantage of the properties of defects in GaN. Point defects in GaN can be detrimental for device operation because they lead to lowering of the breakdown voltage in high-power devices. On the other hand, intentional doping by impurities, such as magnesium, leads to the creation of p-type materials, indispensable in solid-state lighting. One of the most common impurities in GaN is carbon. It is often found as a contaminant in GaN samples grown by a variety of methods. Carbon is present in undoped, Fe-doped, and Si-doped GaN grown by metalorganic vapor phase epitaxy (MOVPE), undoped GaN grown by hydride vapor phase epitaxy (HVPE), and undoped GaN grown by molecular beam epitaxy (MBE). Carbon occupying a nitrogen site in the GaN lattice (CN) acts as an acceptor and is energetically favorable in n-type and high resistivity GaN samples with Fermi energy above the middle of the bandgap.6,7 Early density functional theory (DFT) calculations, which used local density approximation (LDA), predicted a relatively shallow acceptor level of 0.2 eV above the valence band maximum (VBM) for CN.6 However, more recent theoretical calculations using hybrid density functional theory consistently yield a deep acceptor 0/– level at 0.9–1.0 eV above the VBM for CN. They also predict CN to be a source of the yellow luminescence (YL or YL1) band with a maximum at about 2.2 eV in GaN.7–13 

The attribution of the YL1 band to the isolated carbon acceptor is firmly supported by recent experiments. In particular, Zvanut et al.14 have obtained p-type GaN by doping with carbon. The temperature dependence of the hole concentration revealed the ionization energy of the CN acceptor of about 1 eV. Electron paramagnetic resonance (EPR) experiments also showed the presence of CN acceptors that can be ionized with 1 eV light.14 In photoluminescence (PL) studies, the concentration of defects responsible for the YL1 band matches the concentration of carbon estimated from secondary ion mass spectrometry analysis.15 Finally, a second transition level (0/+) at about 0.3 eV above the valence band maximum (VBM) for the CN defect, predicted theoretically,8,9 has been detected experimentally.15,16 The YL1 band is shown in Fig. 1.

FIG. 1.

The YL1 band in GaN co-doped with C and Si. At low temperatures, the ZPL is observed at 2.59 eV and the band maximum is at 2.17 eV. The as-measured spectrum, I(λ), after corrections for the spectral response of the measurement system, is multiplied by λ3 to present it in units proportional to the number of emitted photons as a function of photon energy ħω.15 The solid and dashed lines show steady-state PL spectra at 18 and 150 K, respectively. The filled circles show a time-resolved PL spectrum at 18 K taken at ∼10−5 s after a laser pulse.

FIG. 1.

The YL1 band in GaN co-doped with C and Si. At low temperatures, the ZPL is observed at 2.59 eV and the band maximum is at 2.17 eV. The as-measured spectrum, I(λ), after corrections for the spectral response of the measurement system, is multiplied by λ3 to present it in units proportional to the number of emitted photons as a function of photon energy ħω.15 The solid and dashed lines show steady-state PL spectra at 18 and 150 K, respectively. The filled circles show a time-resolved PL spectrum at 18 K taken at ∼10−5 s after a laser pulse.

Close modal

The following properties of the YL1 band and CN acceptor have been established.15,17 At T < 20 K, the YL1 band is caused by transitions from shallow donors to the 0/− level of CN, the so-called donor–acceptor pair (DAP) recombination. At T > 50 K, electron transitions from the conduction band to the same acceptor replace the DAP transitions. The zero-phonon line (ZPL) of the YL1 band is found at 2.57 eV at T < 20 K and 2.59 eV at T ≈ 50 K in unstrained GaN. Analysis of ZPL in a large number of GaN samples reveals that the 0/− level of the CN is located at 916 ± 3 meV above the VBM. The maximum of this broad band is at 2.17 ± 0.05 eV, and the full width at half maximum is about 0.44 eV at low temperature. The hole-capture and electron-capture coefficients for the CN acceptor are (3.7 ± 1.6) × 10−7 and (1.1 ± 0.1) × 10−13 cm3/s, respectively, and their temperature dependence is insignificant. While YL in GaN is sometimes also attributed to other defects, such as VGa, VGa-3H, and VGaON-2H,18 we do not observe any other YL band except for the YL1 band with the above-listed parameters and properties. The reliable identification and accurately determined parameters of the YL1 band make it an excellent case-study object for comparison with theoretical calculations.

The most widespread theoretical approach for calculating the optical properties of defects in semiconductors is hybrid density functional theory. Among different hybrid functional methods, for calculations of defects in semiconductors, Heyd–Scuseria–Ernzerhof (HSE)19,20 hybrid functional is often used, with functional parameters fit to reproduce the host semiconductor bandgap. Recently, it has also been suggested that correct electronic interactions can be reproduced in HSE if the generalized Koopmans' condition for the defect orbitals is fulfilled instead of (or in addition to) tuning HSE to the experimental gap.13,21–23 In this approach, the HSE parameters are adjusted to cancel the self-interaction energy by wavefunction relaxation with the addition (or removal) of an electron to (from) the defect state orbital. This procedure restores the correct linear behavior of the total energy E(N) with electron occupation in HSE. Usually, this is achieved by finding the HSE parameters, fraction of exact exchange α and range separation parameter μ, for which the electron addition energy is equal to the defect state Kohn–Sham eigenvalue E(N+1)E(N)=ei(N).22 However, within HSE, the set of parameters that fulfills this condition is not unique, and a range of α and μ can satisfy it for a given defect orbital. Therefore, recently, it was suggested that this allows one to tune HSE to both obtain the correct bandgap and enforce the generalized Koopmans' condition for the defect state.13 The tuning procedure involves computing total energies of the defect upon addition/removal of an electron, and in periodic supercells, it depends on the correction scheme for the spurious interaction between periodic images of charged defects. The two most commonly used correction schemes are due to Lany–Zunger (LZ)24 and Freysoldt–Neugebauer–Van de Walle (FNV).25,26 Both approaches yield similar results, within less than 0.1 eV for ∼100 atom supercells,27 but it is not yet clear which method is more appropriate for tuning of the HSE since small total energy differences (less than 0.1 eV) can lead to substantially different values of α and μ that appear to satisfy the generalized Koopmans' condition. Here, we test HSE tuning to enforce the generalized Koopmans' condition using LZ and FNV charged defect correction schemes. We also test a commonly used parametrization of HSE which only reproduces the experimental bandgap, as well as the inclusion of the semicore Ga3d-electrons in the valence band.

All calculations were performed using Vienna Ab initio Simulation Package (VASP)28 using the projector augmented-wave method. HSE hybrid functional calculations were carried out with a varying fraction of exact exchange α and range separation parameter μ in 128 atom hexagonal (4 × 4 × 2) supercells at the Γ point, with plane-wave energy cutoffs of 500 eV. During the initial tuning procedure, lattice parameters were kept at values (a = 3.21 Å, c = 5.20 Å) relaxed for a common gap-tuned HSE (α = 0.31, μ = 0.2). During initial tuning, unrelaxed defect geometries were used, which allow more accurate calculations of charged defect corrections in the FNV scheme, circumventing complications due to ionic screening effects. Electron addition energies were computed as total energy differences of the negatively charged and neutral CN acceptors. Total energies of the charged defects were corrected using both LZ and FNV correction schemes. In the LZ model, the artificial interactions between charged defects are calculated as the electrostatic energy due to periodic images. It contains the well-known Madelung energy of the periodic array of point charges and the energy of interaction of the localized charges with the compensating background, which also includes screening effects. Here, the LZ scheme was also applied in combination with the potential alignment procedure.27 In the FNV model, the defect induced potential is obtained directly from ab initio calculations and is used to match the potential of the model charge density. The correction for the charged defect total energy is then calculated from the Coulomb interaction of the periodic array with model charge density. A detailed review of both correction schemes is given in Ref. 27. In both correction schemes, a high-frequency (ion-clamped) dielectric constant ε=5.0 was used for unrelaxed lattice structures. This value was calculated for GaN using the self-consistent polarization response to the applied electric field. The high-frequency dielectric constant was found to be nearly isotropic, 4.95 in-plane and 5.08 along the wurtzite c-axis, and practically independent of HSE parametrization, within the parameter range used in this work.

Figure 2 shows the planar averaged defect induced electrostatic potentials of negatively charged and neutral CN acceptors (solid lines). The defect is located at the origin, and the first periodic image is 4a away at 12.8 Å. The electrostatic potential of the negative acceptor behaves as expected, with the model potential for FNV correction method25,26 calculated using the above dielectric constant (dashed line) reproducing ab initio calculation very well. However, for the neutral CN, the electrostatic potential is not constant between defects, indicating interactions between periodic images. In previous works, interactions between formally neutral defects were analyzed using an effective defect charge.29,30 Similarly, here these interactions are estimated by adopting an effective defect charge of q = 1 in the FNV model and finding the model potential (dashed line) that approximates the ab initio potential (solid line). This estimation suggests that the energy of spurious interactions of neutral defects is roughly 3.6 times lower than that of negatively charged defects. Therefore, during HSE tuning, FNV and LZ corrections were applied to both negatively charged CN and neutral CN0 defects (for the latter, using an effective dielectric constant of 18).

FIG. 2.

Planar averages of the defect induced potential for neutral and negative charge states of the CN acceptor (solid lines), i.e., the difference of the potential of the supercell containing the defect and that of the bulk. Dashed lines show the FNV correction model potentials.25,26 For the negatively charged acceptor CN, the model potential computed with high-frequency dielectric constant ε=5 reproduces the HSE computed potential very well. In the neutral state of CN0, electrostatic interactions between periodic images of defects are also present and are reproduced by the model potential with effective dielectric constant ε=18.

FIG. 2.

Planar averages of the defect induced potential for neutral and negative charge states of the CN acceptor (solid lines), i.e., the difference of the potential of the supercell containing the defect and that of the bulk. Dashed lines show the FNV correction model potentials.25,26 For the negatively charged acceptor CN, the model potential computed with high-frequency dielectric constant ε=5 reproduces the HSE computed potential very well. In the neutral state of CN0, electrostatic interactions between periodic images of defects are also present and are reproduced by the model potential with effective dielectric constant ε=18.

Close modal

Optical transitions via CN were computed from defect energies following the established procedure.31,32 For each parametrization of HSE, lattice constants of bulk GaN were relaxed minimizing forces to 0.01 eV/Å or less. Atoms in supercells containing defects were relaxed to 0.05 eV/Å. The PL maximum is calculated as

PLMAX=Eg+E(CN0)[E(CN)EVBM],
(1)

where E(CN0) is the total energy of the neutral acceptor, E(CN) is the total energy of the negatively charged acceptor in neutral acceptor lattice geometry, and Eg = 3.5 eV is the experimental bandgap of GaN. ZPL is calculated as

ZPL=Eg+E(CN0)[E(CN)EVBM],
(2)

where E(CN0) and E(CN) are total energies of neutral and negatively charged defects in their respective relaxed lattices. All total energies were corrected applying LZ and FNV electrostatic corrections as outlined above. One important caveat is that high-frequency dielectric constant ε is used for corrections of E(CN) (negatively charged defect in neutral lattice geometry) when calculating the PL maxima, since in the neutral defect lattice, ionic screening effects are negligible. On the other hand, the low frequency (relaxed-ion) dielectric constant ε0 was used for the relaxed negatively charged defect, when calculating ZPLs to account for ionic screening. This is in line with the analysis of defect corrections given in Ref. 27. Our tests show that ε0=8 produces the most reasonable FNV model potentials, albeit it is slightly below the experimental value of 9.5.33 

Figure 3 shows the electron addition energies calculated using the LZ correction scheme (solid black lines), FNV correction scheme (solid red lines), and unoccupied single-particle eigenvalues (dashed lines) of the neutral carbon acceptor for a range of HSE parameters α and μ. Higher values of α indicate larger fractions of exact exchange, which increases the value of the bandgap and moves the electron addition energy and eigenvalues with respect to each other. When the generalized Koopmans' condition is fulfilled, the electron addition energy is

Eadd=E(N+1)E(N)=ei(N),
(3)

where E(N) and E(N + 1) are the total energies of the system with N and N + 1 electrons, i.e., the neutral and negatively charged CN defect, respectively; and ei(N) is the (unoccupied) defect state eigenvalue in the neutral state. In Fig. 3, the intersection points (filled symbols) show the parameter sets for which the generalized Koopmans' condition is fulfilled. Using local approximations to the DFT, such as LDA, leads to the violation of this condition, which can be expressed as22 

EaddLDA=ELDA(N+1)ELDA(N)=eiLDA(N)+Πi+Σi,
(4)

where Πi is the positive self-interaction energy of the electron in the orbital i fixed at the initial unfilled state, and Σi is the negative wavefunction relaxation energy upon addition of an electron. In Hartree–Fock theory (exact exchange is unscreened), the self-interaction energy is zero, and from expression (4) it follows that due to the negative wavefunction relaxation energy, the eigenvalue ei(N) is higher than Eadd. On the other hand, in LDA, where self-interaction energy Πi is significant, the situation is reversed, and the eigenvalue ei(N) is lower than Eadd. Therefore, increasing the fraction of exact exchange from zero to one will move the defect state eigenvalue and Eadd in opposite directions. Similarly, Eadd increases and ei(N) decreases with μ because the effective screening length ls of the exact exchange changes as ls=2/μ. At the points where Eadd is equal to ei(N), we obtain parameters where self-interaction is canceled by wavefunction relaxation (Fig. 3), and Eq. (3) is restored.

FIG. 3.

Electron addition energies (solid lines) and defect state eigenvalues (dashed lines), with their intersections shown as filled circles and diamonds computed using (LZ) and (FNV) correction methods, respectively. These intersections trace HSE parameters μ and α along the dotted lines for which the generalized Koopmans' condition is fulfilled for the 0/– transition level of CN.

FIG. 3.

Electron addition energies (solid lines) and defect state eigenvalues (dashed lines), with their intersections shown as filled circles and diamonds computed using (LZ) and (FNV) correction methods, respectively. These intersections trace HSE parameters μ and α along the dotted lines for which the generalized Koopmans' condition is fulfilled for the 0/– transition level of CN.

Close modal

In the case of the FNV correction scheme, the electron addition energies are shifted upward by 0.1 eV, resulting in intersections with the single-particle eigenvalues that are about 0.06 eV higher. This shifts the HSE parameters that satisfy the generalized Koopmans' condition by a significant amount. The two sets of parameters obtained using the two correction schemes are compared in Fig. 4(a). With one parameter fixed (α or μ), another parameter varies by about 0.05 for the two correction schemes. For instance, depending on the preferred electrostatic correction method, one may conclude that the set either α = 0.3, μ = 0.17 (using FNV corrections) or α = 0.3, μ = 0.22 (using LZ corrections) satisfies the generalized Koopmans' condition. The same could be claimed about the uncertainty in the fraction of exact exchange α. This uncertainty also leads to significant variations in the computed bandgap.

FIG. 4.

(a) HSE parameters μ and α that satisfy the generalized Koopmans' condition using LZ and FNV correction methods. (b) Bulk GaN bandgaps for the same sets of parameters μ and α. The dashed horizontal line shows the experimental bandgap of 3.67 eV (with added renormalization energy due to the electron–phonon interactions). The HSE parameters that satisfy the generalized Koopmans' condition and yield the correct bandgap are α = 0.5, μ = 0.428 when using LZ corrections, and α = 0.36, μ = 0.236 when using FNV corrections.

FIG. 4.

(a) HSE parameters μ and α that satisfy the generalized Koopmans' condition using LZ and FNV correction methods. (b) Bulk GaN bandgaps for the same sets of parameters μ and α. The dashed horizontal line shows the experimental bandgap of 3.67 eV (with added renormalization energy due to the electron–phonon interactions). The HSE parameters that satisfy the generalized Koopmans' condition and yield the correct bandgap are α = 0.5, μ = 0.428 when using LZ corrections, and α = 0.36, μ = 0.236 when using FNV corrections.

Close modal

Figure 4(b) shows the dependence of the bulk bandgap of GaN for the set of parameters that satisfy the generalized Koopmans' condition (using LZ and FNV corrections) for CN. HSE hybrid functionals are commonly tuned to the experimental low-temperature bandgap of 3.50 eV for GaN. However, recently it was suggested that since in experiment the bandgap is renormalized by electron–phonon interactions, and calculations are performed in the static lattice, the correct electron interactions are reproduced for the HSE sets of parameters that yield the non-renormalized bandgap.13 The zero-temperature gap renormalization energy was estimated to be 0.17 eV for bulk GaN.34 Therefore, the target bandgap value, to tune the HSE functional, is 3.67 eV, shown as a dashed line in Fig. 4(b). Assuming that the renormalization affects the valence and conduction band edges similarly for calculations of optical transitions, Eqs. (1) and (2), HSE computed EVBM is adjusted upward by half of the renormalization energy. Figure 4(b) illustrates how sensitive the Koopmans' and bandgap tuned HSE parameters are to the adopted electrostatic correction scheme. Using the LZ correction scheme, the HSE parameters that both satisfy the Koopmans' condition and yield the correct bandgap are α = 0.5, μ = 0.428. Similarly, using the FNV correction scheme, we obtain a significantly different HSE with parameters α = 0.36, μ = 0.236. Interestingly, if neutral defect energies are not corrected for the (small) electrostatic interaction, we also obtain quite different HSE parameters (α = 0.4, μ = 0.297) and (α = 0.3, μ = 0.126) using LZ and FNV corrections, respectively.

It should be noted that tuning HSE to the generalized Koopmans' condition and the bandgap may still not be optimal. For example, it has been suggested that tuning HSE to the Koopmans' condition and bulk ionization potential would yield better optical transition energies relative to the valence band.21 This is because in this case, the VBM should be correct with respect to a physical reference, such as the vacuum level. The problem, however, is that the values of the ionization potential in GaN obtained from experiments are scattered between 6.2 and 7.2 eV (see Ref. 21 and references therein), making reliable tuning of HSE to experimental ionization potential difficult. Fitting HSE to fulfill the generalized Koopmans' condition for a well-localized defect state, such as the deep 0/– level of the gallium vacancy, and the value of GaN ionization potential of 6.7 eV (an average value of available experimental results) yields HSE with parameters α = 0.25, μ = 0.161.21 There is significant uncertainty in whether this parametrization accurately reproduces the bulk ionization potential. Nevertheless, it is illustrative to compare this approach to the parametrizations outlined above.

The results of the calculations of optical properties of CN are shown in Table I. Interestingly, even though different approaches in electrostatic corrections lead to significantly different sets of HSE parameters, the resulting overall calculated energies of optical transitions are relatively close to each other. Tuning HSE to satisfy the generalized Koopmans' condition using the LZ correction scheme produces results slightly closer to experiment, by roughly 0.1 eV, compared to those obtained by tuning HSE using the FNV corrections. Including a small correction for the interactions of neutral defects does not affect the computed ZPL's and PL maxima, while values of μ and α are significantly different. On the other hand, values of PL maxima are all systematically underestimated by 0.1–0.35 eV. This indicates overestimated computed lattice relaxation energies. Tests for HSE (α = 0.25, μ = 0.161) using 300 atom supercells yield very similar results (ZPL = 2.58 eV, PLMAX = 2.02 eV), indicating that this overestimation is unlikely due to the supercell size. We also note that a commonly used tuning of α = 0.31, μ = 0.2 produces similar reasonable results. Tests show that for this parametrization, the electron addition energy (using LZ corrections) is 0.07 eV away from the neutral defect eigenvalue, i.e., it is close to being Koopmans' compliant.

TABLE I.

Optical properties of the CN acceptor computed with HSE tuned to satisfy the generalized Koopmans' condition using LZ and FNV electrostatic correction methods. The same correction method was used in each case for optical properties’ calculations. The 0/– transition level is given relative to the VBM. Theoretical bandgaps are calculated in the primitive cells relaxed for each HSE. The experimental bandgap of 3.5 eV is used in optical transition calculations. For HSE parametrizations labeled with an asterisk, only the negative charge state of CN was corrected for both HSE tuning and optical calculations. The values that are closest to the experiment are given in boldface.

HSE α, μ−1)ZPL (eV)PLMAX (eV)0/– level (eV)Eg (eV)a (Å)c (Å)
0.5, 0.428 (LZ) 2.52 1.89 0.98 3.66 3.193 5.181 
0.36, 0.236 (FNV) 2.46 1.82 1.03 3.67 3.200 5.193 
0.4, 0.297 (LZ)* 2.53 1.91 0.97 3.68 3.197 5.188 
0.3, 0.126 (FNV)* 2.40 1.81 1.095 3.68 3.205 5.201 
0.25, 0.161 (LZ)21  2.59 2.06 0.91 3.22 3.212 5.215 
0.26, 0.1 (Ga-d) (LZ)13  2.56 2.05 0.94 3.63 3.181 5.171 
0.31, 0.2 (LZ) 2.50 1.93 1.00 3.49 3.205 5.202 
Experiment 2.587 2.17 0.916 3.67 3.189 5.185 
HSE α, μ−1)ZPL (eV)PLMAX (eV)0/– level (eV)Eg (eV)a (Å)c (Å)
0.5, 0.428 (LZ) 2.52 1.89 0.98 3.66 3.193 5.181 
0.36, 0.236 (FNV) 2.46 1.82 1.03 3.67 3.200 5.193 
0.4, 0.297 (LZ)* 2.53 1.91 0.97 3.68 3.197 5.188 
0.3, 0.126 (FNV)* 2.40 1.81 1.095 3.68 3.205 5.201 
0.25, 0.161 (LZ)21  2.59 2.06 0.91 3.22 3.212 5.215 
0.26, 0.1 (Ga-d) (LZ)13  2.56 2.05 0.94 3.63 3.181 5.171 
0.31, 0.2 (LZ) 2.50 1.93 1.00 3.49 3.205 5.202 
Experiment 2.587 2.17 0.916 3.67 3.189 5.185 

Here, we also compare HSE parametrization used in Ref. 13 (α = 0.26, μ = 0.1), obtained by tuning HSE to Koopmans' condition and the bandgap, including semicore Ga3d-electrons in the valence band. In their respective relaxed primitive cells, the bandgap computed with this parametrization of HSE including Ga3d-electrons is 3.63 eV, while removing them reduces the bandgap to 3.49 eV. The Ga3d-electrons push the valence band up, which changes the position of the VBM; however, the changes in the exact exchange energy also increase the bandgap. Including Ga3d-electrons does not change the position of the defect state eigenvalue above the VBM. However, it moves Eadd down somewhat, which means that without Ga3d-electrons Eadd is higher than the defect state eigenvalue, indicating an excess of self-interaction energy. In other words, such HSE parametrization without Ga3d-electrons is slightly closer to the local approximations of DFT.

Better agreement with the experimental results is obtained with HSE parameters from Refs. 13, 15 and 21. This is especially surprising for HSE (α = 0.25, μ = 0.161), given the fact that this parametrization was obtained without tuning to the bandgap, i.e., it yields Eg = 3.22 eV, significantly lower than the experimental (non-renormalized) value of 3.67 eV. Nevertheless, optical transitions obtained by this HSE parametrization are very similar to those obtained using HSE tuned with Ga3d-electrons, which is a significantly more computationally demanding calculation. This suggests that reproducing the experimental bandgap is less important than satisfying the generalized Koopmans' condition.

Using LZ and FNV electrostatic correction schemes, we obtained several parametrizations of the HSE hybrid functional, which satisfy the generalized Koopmans' condition for the CN acceptor and reproduce the experimental bandgap of GaN. Optical properties of the CN acceptor calculated using these parametrizations of HSE are compared with experimental data: (1) The HSE parameters are sensitive to the electrostatic correction scheme used during tuning. Nevertheless, all Koopmans' tuned HSE parametrizations yield reasonable results. (2) The results that agree better with experiment are obtained using Koopmans' gap-tuned HSE with Ga3d-electrons. (3) Good results can also be obtained without Ga3d-electrons, by tuning HSE to the Koopmans' condition and ionization potential. (4) The results of the latter HSE parametrization (α = 0.25, μ = 0.161) for CN are in close agreement with the experiment. This parametrization was also shown to reproduce the experiment for the magnesium acceptor MgGa and fulfill the Koopmans' condition for the gallium vacancy VGa in GaN.21 This suggests that Koopmans' tuned HSE parametrization is transferrable and can be used for a variety of defects in GaN. This is also in agreement with earlier studies that suggested that the same set of hybrid functional parameters satisfies the Koopmans' condition for a series of defects.23 

This work was supported by the National Science Foundation (NSF) (No. DMR-1904861). The calculations were performed at the VCU Center for High Performance Computing.

1.
S.
Nakamura
,
T.
Mukai
, and
M.
Senoh
,
Appl. Phys. Lett.
64
,
1687
(
1994
).
2.
F. A.
Ponce
and
D. P.
Bour
,
Nature
386
,
351
(
1997
).
3.
R.
Dahal
,
J.
Li
,
K.
Aryal
,
J. Y.
Lin
, and
H. X.
Jiang
,
Appl. Phys. Lett.
97
,
073115
(
2010
).
4.
Y.
Saitoh
,
K.
Sumiyoshi
,
M.
Okada
,
T.
Horii
,
T.
Miyazaki
,
H.
Shiomi
,
M.
Ueno
,
K.
Katayama
,
M.
Kiyama
, and
T.
Nakamura
,
Appl. Phys. Express
3
,
081001
(
2010
).
5.
M.-W.
Ha
,
C. H.
Roh
,
D. W.
Hwang
,
H. G.
Choi
,
H. J.
Song
,
J. H.
Lee
,
J. H.
Park
,
O.
Seok
,
J.
Lim
,
M.-K.
Han
, and
C.-K.
Hahn
,
Jpn. J. Appl. Phys.
50
,
06GF17
(
2011
).
6.
A. F.
Wright
,
J. Appl. Phys.
92
,
2575
(
2002
).
7.
J. L.
Lyons
,
A.
Janotti
, and
C. G.
Van de Walle
,
Appl. Phys. Lett.
97
,
152108
(
2010
).
8.
J. L.
Lyons
,
A.
Janotti
, and
C. G.
Van de Walle
,
Phys. Rev. B
89
,
035204
(
2014
).
9.
D. O.
Demchenko
,
I. C.
Diallo
, and
M. A.
Reshchikov
,
Phys. Rev. Lett.
110
,
087404
(
2013
).
10.
M. A.
Reschchikov
,
D. O.
Demchenko
,
A.
Usikov
,
H.
Helava
, and
Y.
Makarov
,
Phys. Rev. B
90
,
235203
(
2014
).
11.
M.
Matsubara
and
E.
Bellotti
,
J. Appl. Phys.
121
,
195701
(
2017
).
12.
S. G.
Christenson
,
W.
Xie
,
Y. Y.
Sun
, and
S. B.
Zhang
,
J. Appl. Phys.
118
,
135708
(
2015
).
13.
P.
Deák
,
M.
Lorke
,
B.
Aradi
, and
T.
Frauenheim
,
Phys. Rev. B
99
,
085206
(
2019
).
14.
M. E.
Zvanut
,
S.
Paudel
,
E. R.
Glaser
,
M.
Iwinska
,
T.
Sochacki
, and
M.
Bockowski
,
J. Electron. Mater.
48
,
2226
(
2019
).
15.
M. A.
Reshchikov
,
M.
Vorobiov
,
D. O.
Demchenko
,
Ü
Özgür
,
H.
Morkoç
,
A.
Lesnik
,
M. P.
Hoffmann
,
F.
Hörich
,
A.
Dadgar
, and
A.
Strittmatter
,
Phys. Rev. B
98
,
125207
(
2018
).
16.
T.
Narita
,
K.
Tomita
,
Y.
Tokuda
,
T.
Kogiso
,
M.
Horita
, and
T.
Kachi
,
J. Appl. Phys.
124
,
215701
(
2018
).
17.
M. A.
Reshchikov
,
J. D.
McNamara
,
F.
Zhang
,
M.
Monavarian
,
A.
Usikov
,
H.
Helava
,
Y.
Makarov
, and
H.
Morkoç
,
Phys. Rev. B
94
,
035201
(
2016
).
18.
J. L.
Lyons
,
A.
Alkauskas
,
A.
Janotti
, and
C. G.
Van de Walle
,
Phys. Status Solidi B
252
,
900
(
2015
).
19.
J.
Heyd
,
G. E.
Scuseria
, and
M.
Ernzerhof
,
J. Chem. Phys.
118
,
8207
(
2003
).
20.
J.
Heyd
,
G. E.
Scuseria
, and
M.
Ernzerhof
,
J. Chem. Phys.
124
,
219906
(
2006
).
21.
D. O.
Demchenko
,
I. C.
Diallo
, and
M. A.
Reshchikov
,
Phys. Rev. B
97
,
205205
(
2018
).
22.
S.
Lany
,
Phys. Status Solidi B
248
,
1052
(
2011
).
23.
G.
Miceli
,
W.
Chen
,
I.
Reshetnyak
, and
A.
Pasquarello
,
Phys. Rev. B
97
,
121112
(
2018
).
24.
S.
Lany
and
A.
Zunger
,
Phys. Rev. B
78
,
235104
(
2008
).
25.
C.
Freysoldt
,
J.
Neugebauer
, and
C. G.
Van de Walle
,
Phys. Rev. Lett.
102
,
016402
(
2009
).
26.
C.
Freysoldt
,
J.
Neugebauer
, and
C. G.
Van de Walle
,
Phys. Status Solidi B
248
,
1067
(
2011
).
27.
H.-P.
Komsa
,
T. T.
Rantala
, and
A.
Pasquarello
,
Phys. Rev. B
86
,
045112
(
2012
).
28.
G.
Kresse
and
J.
Furthmuller
,
Phys. Rev. B
54
,
11169
(
1996
).
29.
F.
Oba
,
A.
Togo
,
I.
Tanaka
,
J.
Paier
, and
G.
Kresse
,
Phys. Rev. B
77
,
245202
(
2008
).
30.
D.
Skachkov
,
A.
Punya Jaroenjittichai
,
L.
Huang
, and
W. R. L.
Lambrecht
,
Phys. Rev. B
93
,
155202
(
2016
).
31.
C. G.
Van de Walle
and
J.
Neugebauer
,
J. Appl. Phys.
95
,
3851
(
2004
).
32.
C.
Freysoldt
,
B.
Grabowski
,
T.
Hickel
,
J.
Neugebauer
,
G.
Kresse
,
A.
Janotti
, and
C. G.
van de Walle
,
Rev. Mod. Phys.
86
,
253
(
2014
).
33.
A. S.
Barker
and
M.
Ilegems
,
Phys. Rev. B
7
,
743
(
1973
).
34.
M.
Cardona
and
M. L. W.
Thewalt
,
Rev. Mod. Phys.
77
,
1173
(
2005
).