Complexes and compensation in degenerately donor doped GaN

Gallium nitride is a direct, wide bandgap semiconductor material widely used in solid state lighting, making significant inroads into the power electronics market.^1,2 Ongoing research is exploring its potential for plasmonics, quantum information applications, and radiation- and high-temperature-tolerant integrated circuitry.^1,3–6 Silicon and germanium are the two primary n-type dopants in this material due to their shallow activation.^3,7–9 However, at very high doping levels, compensation effects arise for both dopants that limit their efficacy and restrict the maximum achievable electron concentration to different amounts.

It has recently been demonstrated in AlN that the silicon self-compensation effect is caused by the geometrically increasing favorability of multi-donor $vAl$-$n$$SiAl$ complexes vs the isolated dopant with increasing Si availability.¹⁰ It is also well-known that defects in the III-nitrides tend to exhibit similar behavior across cation chemistries, with the different dopant behaviors arising primarily from differences in the band edge locations, and that Si and Ge exhibit fundamentally different responses in AlN due to slight differences in defect behavior.¹¹ Based on this, it was hypothesized that (a) such complexes could be responsible for the compensation effects encountered in degenerate doping conditions in GaN and (b) differences between Si and Ge manifest with the above conduction band Fermi level in degenerate conditions, leading to differences in compensation behavior. In this article, this hypothesis is explored using first principles calculations based on density functional theory (DFT) and semiconductor compensation calculations based on the grand canonical thermodynamics ensemble.

All DFT calculations were performed with the HSE06^12,13 screened hybrid exchange correlation functional in VASP 5.3.3 with an exact exchange amount of 0.289, collinear spin polarization, and a plane wave kinetic energy cutoff of 520 eV.^14–17 The amount of exact exchange was chosen to correct the underestimation of the bandgap common to traditional functionals by tuning the bandgap based on the results of Monemar,¹⁸ who carefully examined the cryogenic electronic behavior of GaN to separate the fundamental band edge from excitonic effects. All defect calculations were performed in a 96 atom supercell with a 2 × 2 × 2 k-point grid. The electronic density of states was calculated using a 13 × 13 × 13 k-point mesh. PAW (Projector Augmented Wave) pseudopotentials contained 13, 5, 14, and four electrons for Ga, N, Ge, and Si, respectively.

The formation energy of a defect is given by Eq. (1). This expression balances the energy change from introducing a defect (⁠$EDqtot−Ebulktot$⁠) with exchange of atomic (⁠$−∑iniμi$⁠) and electronic (⁠$q(μe+Ev)$⁠) species with the respective reservoirs. Finally, $EDqcorr$ is a finite size correction, based on the method of Kumagai and Oba and using a relative dielectric constant of 8.9, taken from experiment;^19,20 this yields defect formation energies typically converged to within 200 meV of much larger 360 atom supercells, even for neutral defects,¹⁰ 

With an appropriate thermodynamics model, this expression can be used to examine defect formation energies as a function of processing environment. By imposing charge neutrality, the defect concentrations can also be self-consistently determined for a given processing environment and impurity content, yielding the equilibrium concentrations of all point defects and charge carriers in the system. This process self-consistently determines the Fermi level and impurity chemical potentials and, thus, the distribution of a dopant across all of its defect forms for a given processing/doping combination. The prefactor to the concentration expressions includes terms related to the defects' configurational entropy. These concentrations can be further used to obtain the room temperature point defect and charge carrier concentrations by self-consistently solving for reionization. Further details on these methods may be found in previous work.^10,21–25

The temperature dependence of the bandgap observed by Monemar¹⁸ and Su et al.²⁶ was used to calculate valence edge and conduction edge movement using the method described in Ref. 27. Native chemical potentials (that is, $μGa$ and $μN$⁠) were selected here to correspond to extremely nitrogen-rich growth conditions of the kind expected in an metalorganic chemical vapor deposition (MOCVD) reactor flowing ammonia and hydrogen, swept between low and high $μN$ conditions. This allows for a qualitative evaluation of how growth conditions might impact compensation. Because of the constraints imposed by equilibrium and the stoichiometry of the compound, the increase in $μN$ leads to a 1:1 decrease in $μGa$⁠. It should be stressed that the compensation and defect ensemble behavior discussed in this article are general to most growth methods, but minor variations in predicted properties are expected as each technique has access to unique chemical potential regimes.^28–31 This is because, unlike stoichiometric complexes (e.g., v_Ga–v_N), the donor complex's formation energies depend on the native chemical potentials in the growth environment. For this work, native vacancies, divacancy complexes, donor substitutionals, and complexes between the gallium vacancy and 1, 2, and 3 donor substitutionals were simulated and were managed and analyzed with the asphalt point defect simulation informatics suite.²⁵ 

As discussed in Ref. 10, there are 58 geometrically inequivalent configurations for these second nearest neighbor metal vacancy n-donor complexes up to n = 3. The primary factor in formation energy ordering of vacancy multi-silicon complexes in that work was found to be electrostatic interactions,²⁷ and many of the configurations were consequently found to not contribute strongly to the charge balance and, thus, position of the Fermi level. Subsequent analysis identified a much smaller set of 16 configurations that very closely reproduced the electron vs donor concentration curves calculated using the full dataset. Given these factors, it was assumed that this reduced set plus the isolated donor (which will be referred to as the onsite donor throughout the rest of this manuscript) could describe compensation in GaN to a similarly high degree of numeric precision for both Si and Ge. Formation energy diagrams for the simulated defects and binding energy diagrams for the simulated complexes are shown in Fig. 1, while thermodynamic transition level (TTL) ranges for donor-related defects are listed in Table I. This work adopts the convention of Ref. 32 for binding energies, where a positive binding energy indicates energetic favorability of defect association and the binding energy is taken as the negative of the difference between the minimum formation energies of the complex and its constituents as a function of Fermi level. As such, the complexes with the lowest formation energy will have the highest binding energy, i.e., the uppermost lines of each color in Figs. 1(c) and 1(f) correspond to the lowest lines of each color in Figs. 1(a) and 1(b) and 1(d) and 1(e), respectively.

Each dopant is shown in low $μN$ conditions at low (10¹⁷ cm⁻³) [Figs. 1(a) and 1(d)] and high (⁠$5×1020$ cm⁻³) [Figs. 1(b) and 1(e)] dopant concentrations. The onsite defects are seen to always dominate for Fermi levels within the bandgap for both dopants. Both isolated and complexed donor defects exhibit transitions well above the conduction band edge, which become relevant in certain conditions, and so they have been tabulated along with the in-gap TTLs in Table I and emphasized to indicate their above band edge nature. The onsite donors display what appears in Fig. 1 to be DX transitions (a transition involving both electron capture and offsite displacement of the capturing atom) above the band edge. However, these are not true DX transitions, as they first do not involve the displacement of the donor atom, as the a- and c-DX configurations were found to be less stable than the onsite, and second the 0 charge state has very nearly the same formation energy as the +1 and −1 charge states at the (1 | −1) crossing.

The 16 complex configurations considered consisted of the three configurations of the vacancy-1 donor complexes, seven of the vacancy-2 donor complexes, and six of the vacancy-3 donor complexes. Based on similar locations of the thermodynamic transition levels in the complexes and the Ga vacancy, and the simple donor character of the onsite donors, the states of the multi-donor complexes seem to be derived from the Ga vacancy transitions, perturbed by the presence of the donors. This is most obvious in Fig. 1. Once the Ga vacancy associated levels are filled, the net charge of the complexes over the upper half of the bandgap comes out to the sum of the formal oxidation states of the vacancy and the donors. The single donor complexes have minimal scatter as compared to the higher order complexes. This is primarily geometric: all three single donor complexes are essentially identical aside from the anisotropy of the wurtzite lattice. Because of this, they all have very similar electrostatic interactions and, hence, minimal formation energy scatter. In contrast, the two and three donor complexes each have multiple donors occupying different combinations of Ga sites around the vacancy, leading to more dramatic variation in the electrostatic energy and, thus, to more dramatic scatter in the formation energies as compared to the single donor complexes. These variations are further perturbed by ligand field differences around the vacancy site, leading to the emergence of different TTL groupings for some charge states. To borrow the language of coordination chemistry, the two and three donor complexes are made up of a vacancy and different combinations of 1N-1 donor and 1N-2 donor functional groups. Among those complexes in the reduced dataset, the two donor complexes with different TTLs from the rest have a 1N-2 donor functional group, while the three donor complexes with different TTLs consist solely of 1N-1 donor functional groups. As can be seen in Fig. 1, all complexes have a binding energy above 1 eV for Fermi levels in the upper half of the bandgap and for a considerable range above the conduction band. The Si complex binding energies are generally higher than the corresponding Ge complex binding energies, reflecting their increased stability at a given availability.

However, the complexes do not become energetically relevant until two conditions are met: the Fermi level must reach or exceed the conduction band (i.e., degenerate doping conditions), and the availability of the dopant must be sufficiently high. This regime is shown in Figs. 1(b) and 1(e). The first condition is a consequence of the smaller bandgap in GaN as compared to that in AlN. In AlN, the much higher conduction band means that complexes become favorable at Fermi levels still inside the bandgap, which, in turn, leads to their onset at lower dopant concentrations in that material system.¹⁰ 

The second condition is related to how the impurity concentration is tied to the chemical potentials. Because the different modes of donor incorporation involve between one and three donor dopants counterbalanced by between one and four displaced cations, increasing the concentration geometrically increases the relative favorability of the higher donor-count complexes over the lower donor-count complexes, as a simple consequence of Eq. (1). However, the isolated donor substitutionals have a large head start, energetically speaking, over any of the complexes, and so the complexes will not become more energetically favorable than the isolated donors until quite high availabilities. By the time both conditions are met, the three donor complex is more favorable than either the one donor or the two donor complex, and so when the Fermi level is high enough to make complexes more favorable than the onsite, most of the donor switches from incorporating as onsite directly to incorporating as v_Ga-3Ge_Ga. The three donor complexes themselves are neutral and consume dopants rather than directly compensating them, with the majority of compensation arising instead from the two donor complexes, which act as acceptors.

The collective behavior of these defect ensembles and their interactions with the growth environment govern the donor efficacy of these dopants, which are shown in Figs. 2 and 3 for Ge and Si, respectively. Subfigures 2(a) and 3(a) show predictions of the room temperature electron concentration vs donor doping level, showing how the curves are expected to change under variations in growth conditions when changing $μN$⁠. The insets illustrate the linear, processing-invariant donor activation in most regimes, which is what initially led to the success of these donors in this material system. Note that the electron-vs-donor curves will die off quickly as the donor dopant falls below the level of compensating background impurities (such as in the case of MOCVD carbon), which tend to “clamp” the onset of the electron concentration by counteracting the donor's effect on the Fermi level. Note that at simultaneously high C and donor levels, a large amount of C may also incorporate as donor-C complexes instead of C_N, as discussed for AlN.³³ 

Germanium will be discussed first. The room temperature effective density of states calculated from the band curvatures is marked with a horizontal line in the inset of Fig. 2(a). This point is passed at relatively modest (when compared to what is routinely obtained in experiment) Ge doping levels in the mid 10¹⁸ cm⁻³, with the first compensation effects not showing up for another order of magnitude. Past this point, different growth conditions will lead to different maximum achievable electron concentrations. The best results are obtained for low $μN$ growth conditions (corresponding, for example, to lower ammonia partial pressures in the case of MOCVD).

The self-compensation effect itself is a result of multiple mechanisms activating at similar doping levels, shown in Fig. 2(b). First, the deviation from the linear activation: this is caused by the neutral three donor complex consuming enough of the dopant to cause a nonlinearity in the dopant-vs-onsite ratio, while the charged two donor complex reaches a high enough concentration to contribute to the charge balance. The flattening out of the response is from the charged onsite and the charged two donor complexes having similar shapes but different onsets, while the neutral three donor complex consumes the excess dopant. Finally, the dip in the electron concentration seen on the high side is caused by the single donor complex becoming high enough in concentration to contribute to the charge balance, tipping it in favor of a decreasing electron concentration with the increasing dopant content. While these mechanisms activate, a small fraction of the above band edge defect states fill, neutralizing some of the onsite donor and making the compensating effects of the complexes more severe, although these are very minor effects compared to the other changes, particularly in the case of the complexes. Because of the Fermi level dependence of occupation probability, most of the onsite donor remains ionized as +1, a small amount of the donor becomes neutral, and an even smaller amount becomes an acceptor. Although it may appear from Table I that above band edge states associated with the three donor complex would contribute strongly to the charge balance by filling before the onsite donors, the configurations with TTLs lower than the onsite are of such low concentration in these conditions that they never contribute more than a hundredth of a percent of the charge balance, at least in the conditions simulated. As the nitrogen chemical potential in the growth environment is increased (and $μGa$ is decreased), the onsite donor becomes less favorable with respect to the complex, shifting the charge balance and resulting in a decrease in the maximum achievable electron concentration. This is shown in Fig. 2. Changes with nitrogen chemical potential are expected to be even more dramatic when moving between different growth techniques and could enable swings in the maximum electron concentration between low 10¹⁹ and mid 10²⁰ per cubic centimeter.

The same set of results for Si is shown in Fig. 3, with the highest obtained electron concentration for Ge shown as a gray line in Fig. 3(a) for reference. The same major mechanisms control the behavior, with minor energetic differences in the defects leading to a lower maximum achievable electron concentration than that found for Ge. Specifically, the multi-donor Si complexes are more favorable relative to their onsite than the multi-donor Ge complexes and are more strongly bound (see Fig. 1).

Fundamentally, this difference has to do with silicon's lower electronegativity than germanium and the effect this has on the strain fields in and around the multi-donor complexes (i.e., change in bond lengths in the first and second nearest neighbor shells around the complex defect centers). Because of the electronegativity difference, Si forms stronger and shorter bonds with surrounding N atoms (1.78 $Å$ for $SiGa+1$⁠) than Ge (1.87 $Å$ for $GeGa+1$⁠). On its own, this means that Si should have a higher strain associated with its incorporation than Ge and, thus, a higher formation energy at a given availability. However, while the Ge–N bonds are much closer to GaN's normal Ga–N bond length (1.95 $Å$⁠) than the Si–N bonds, once the donors begin complexing with Ga vacancies, these bonds interact with the latter's positive strain field, which normally warps the surrounding nitrogen atoms outward from the vacancy core. Because silicon's normal Si–N bond is shorter than germanium's, its complexes provide greater strain relief by more fully accommodating the strain fields of both the donor atoms and the Ga vacancies, increasing the favorability of the Si complexes vs the isolated Si as compared to Ge, which is what ultimately leads to germanium's better performance as a donor dopant in this material system. Although these calculations do not include vibrational free energy, preliminary work on including vibrational effects shows that the strain field differences around the dopants cancel most of the mass difference effects, leading to only minor changes in the electron concentrations, while the overall behaviors and findings discussed in this manuscript remain unchanged. Work in Ref. 34 suggests that these preliminary vibrational energies are likely converged to within 10 meV of results obtained with larger 216 atom supercells.

Because of these factors, Ge is intrinsically a more effective donor dopant than Si in GaN, even though the differences do not manifest until the degenerate doping regime. Even then, differences in growth conditions can obscure the differences between the dopants, and so the growth environment must be carefully managed to maximize donor efficacy by minimizing the compensating effects of the multi-donor-vacancy complexes.

These findings also have interesting implications for defect engineering design rules with respect to dopant selection. It is fairly well known that dopants will tend to complex with native defects if they are of opposite charge, and this is typically detrimental to the intended function of the dopant. While this can sometimes be avoided by changing the processing or growth regime to modulate complex favorability, it is clear from this work that choosing a chemically similar dopant (i.e., similar valence shell) with different electronegativity or atomic radius may also be effective, by changing the strain fields associated with the defect complexes.

The authors received financial support from AFOSR Grant No. FA9550-17-1-0225 and computer time from the DoD HPCMP. The authors thank S. Washiyama for fruitful discussions.