Computational study of positron annihilation parameters for cation mono-vacancies and vacancy complexes in nitride semiconductor alloys

We calculate positron annihilation parameters, namely the S and W parameters from the Doppler broadening spectroscopy and the positron lifetime , for defect-free states as well as cation mono-vacancies and vacancy complexes in nitride semiconductor alloys Al0.5Ga0.5N, In0.5Ga0.5N and Al0.5In0.5N. The obtained distributions of these parameters differ from compound to compound. Especially, the S–W relation for In0.5Ga0.5N is very different from that for Al0.5Ga0.5N. For the cation mono-vacancies, introducing local structural parameters, their correlations with S, W and are investigated. The S and variations are well described with the size distributions of the vacancies while the W variation is related to the presence of localized d electrons. For the vacancy complexes as well as the cation mono-vacancies, multiple-linear-regression models to describe S, W and are successfully constructed using the local structural parameters as descriptors. The S–W and S– relations are also compared with those for AlN, GaN and InN.


Introduction
Nitride semiconductors such as AlN, GaN and InN attract increasing attention not only for applications in light-emitting devices but also for those in power devices. Point defects therein often play a critical role in device performances. Since positrons tend to be trapped at vacancy-type defects, the positron annihilation technique has been utilized to investigate them by measuring positron annihilation Doppler broadening spectra and/or positron lifetime spectra [1,2].
We have been successfully applying the positron annihilation technique in identifying defect species in various nitride semiconductor systems. For example, vacancy type defects in AlN specimens prepared by metal-organic vapor phase epitaxy (MOVPE), molecular beam epitaxy (MBE) and lateral epitaxial overgrowth (LEO) were investigated [3]. For MOVPE and MBE, the major defect species can be identified as vacancy cluster such as V Al -nV N (n = 1-3), while, for LEO, vacancy-oxygen complexes such as V Al -nO N (n = 3, 4) seem to have significant contributions. In Mg-doped GaN, V Ga -nV N (n = 2, 3) were observed [4] and they were identified as the major non-radiative recombination centers [5]. In the above description, V A represents a vacancy at the A site. In practical applications of nitride semiconductors, the alloying technique is used to control the band gap. We have experimentally investigated (Al,Ga)N and (In,Ga)N systems and observed a significant concentrations of vacancy-type defects in high-Al content samples for the former [6], and in samples close to the In 0.5 Ga 0.5 N composition for the latter [7].
Theoretical predictions of positron states and annihilation parameters have often been found to be crucial in enabling specific defect species to be identified from the results of experiments. In consequence, we have been calculating positron lifetime and positron-electron momentum distributions for various types of defects in nitride semiconductors [8]. For alloy systems, in general, each of the nominally similar atomic sites are in fact not equivalent due to the random distribution of the two (or more) cations on the neighbor sites. The calculated positron annihilation parameters show finite distributions as described in our previous studies for cation monovacancies (V c ) and cation-nitrogen di-vacancies (V c -V N ) in Al 0.5 Ga 0.5 N and In 0.5 Ga 0.5 N [8][9][10] and in the previous work by Tuomisto et al for In x Ga 1−x N [11]. As mentioned above, experimental studies provide evidence for the presence of larger vacancy complexes, hence theoretical predictions of the positron annihilation parameters for such defects are required.
In the present study, we calculate positron annihilation Doppler broadening spectra and positron lifetime (τ ) values for defect-free (DF) states as well as for cation mono-vacancies and vacancy complexes V c -nV N (n = 0-3) in nitride semiconductor 'disordered' alloys Al 0.5 Ga 0.5 N, In 0.5 Ga 0.5 N and Al 0.5 In 0.5 N. The Doppler broadening spectra are analyzed by evaluating S and W parameters. For the vacancy-type defects V c -nV N , there are a variety of configurations reflecting the disordered nature, and the obtained S, W and τ values show finite distributions which differ from compound to compound. We use several data-scientific methods to analyze the calculated results. The differences among the three compounds are explained in terms of local structural properties of the vacancies by evaluating linear correlation coefficients between the positron annihilation and local structural parameters. The characteristic behaviors of S, W and τ are successfully modeled by the multiple linear regression using the local structural parameters. Furthermore, to figure out overall relationships among the parameters, we perform principal component analyses. For each compound and n, the averaged values of S, W and τ are evaluated and compared with the values for the DF states as well as those for AlN, GaN and InN.

Computational methods
For each calculation in the present study, an orthorhombic supercell equivalent to 4 × 4 × 2 hexagonal unit cells of the wurtzite structure is used. There are 64 cations and 64 N atoms in the supercell if there is no vacancy. The alloy structures for Al 0.5 Ga 0.5 N, In 0.5 Ga 0.5 N and Al 0.5 In 0.5 N are generated using the special quasirandom structure (SQS) approach [12]. The numbers of the possible vacancy complexes V c -nV N in the above-mentioned supercell are 64, 256, 384 and 256 for n = 0, 1, 2 and 3, respectively. The lattice parameters for the alloy supercells are evaluated by averaging the experimental lattice parameters of the appropriate pair of wurtzite binary structures (AlN, GaN and InN).
Electron and positron states are calculated with our computational code QMAS [13], which adopts the projector augmented-wave (PAW) method [14] and the planewave basis. The generalized gradient approximation for the electronic exchange-correlation energy is used [15]. Internal atomic positions are optimized with evaluating the Hellmann-Feynman force. It has been shown previously that the effects of the positron localization on the atomic configurations and positron annihilation parameters [16][17][18] are insignificant for AlN, GaN and InN [9,19], hence they are not considered in the present study. The vacancy-type defects V c -nV N are thought to be in different charge states depending on the Fermi level [22]. For example, in GaN, the charge state of V Ga varies from 0 to 3− according to density-functional-theory calculations by Neugebauer and Van de Walle [23]. Recent hybrid functional calculations by Diallo and Demchenko [24] suggest +, 0, −, 2− and 3− charge states for V Ga while 3+, +, 0, − and 2− states for V Ga -V N . To investigate the effect of the charge state on the annihilation parameters, we have calculated positron annihilation parameters for 64 types of V c in Al 0.5 Ga 0.5 N using supercells with the charge states of 0 and 3 − . For the 3 − cases, additional three electrons are distributed on N atoms in the vicinity of V c similarly to the previous calculations for GaN [23,24]. Although the obtained S, W and τ for each V c are slightly different between the two charge states, their overall distributions show rather simple parallel shifts. We therefore use neutral supercells throughout the present study. First we obtain the electronic structure with the atomic positions optimized, and then, fixing the electronic structure and atomic positions, the positron state is calculated. The electron-positron correlation energy is included using the local-density-approximation expression obtained for a positron in a homogeneous electron gas [16]. The planewave energy cutoff and the force convergence criterion are set to be 20 Ha and 5 × 10 −5 Ha/bohr, respectively. Only the Γ point is used for the Brillouin zone sampling. We also calculate positron annihilation parameters based on electronic structures obtained using 2 × 2 × 2 k points for 32 cases of V c in Al 0.5 Ga 0.5 N. The results are very similar to those obtained with the Γ-point sampling. Positron states are calculated with the positron effective potential constructed using core electron wavefunctions from isolated atom calculations and valence electron wavefunctions described in the PAW formalism [20]. With the obtained positron and electron wavefunctions, positron annihilation Doppler broadening spectra ρ 1D ( p z ), which corresponds to the one-dimensional projection of the electron-positron momentum density ρ(p), and positron lifetimes τ are calculated as follows: where ψ +, ψ −j , n + and n − are the positron wavefunction, the j th electron wavefunction, the positron density and the electron density, respectively. The direction of γ-ray emission corresponds to the z direction in Cartesian coordinates, r 0 represents the classical electron radius and c is the speed of light, and γ is the enhancement factor. In the present study of semiconductors, we use the expression where ε ∞ is the dielectric constant of the target system and γ BN is the Boroński-Nieminen enhancement factor [16] expressed with the density parameter r s = [3/(4πn − )] 1/3 as The electron-positron correlation energy is corre spondingly modified by scaling [21]. As for ε ∞ , averages of the experimental values for AlN, GaN and InN are used: 5.10 for Al 0.5 Ga 0.5 N, 6.95 for In 0.5 Ga 0.5 N and 6.55 for Al 0.5 In 0.5 N. The Doppler broadening spectra ρ 1D are calculated using spherically averaged momentum densities ρ( p), since directionally averaged analyses are applied in the following. For spherical averaging, first, interpolated values of ρ(p) are obtained with an interval of 0.02 a.u. along 16 directions using the tricubic interpolation method by Lekien and Marsden [25], and then averaged. The 16 directions derive from the 32   icosahedral-group points for numerical integration on a sphere [26]. The obtained spectra are convoluted with a Gaussian whose full width at half maximum (FWHM) is 1.0 keV to mimic the experimental resolution. The resultant spectra are hereinafter referred to as ρ * 1D ( p z ). The positron annihilation S parameter is defined as the ratio of the center area in the momentum region between ±0.419 a.u. to the total spectrum area while the W parameter is that for the higher momentum areas of ±(1.828 ∼ 3.688) a.u.: (7) For pairs of data (x i ; y i ), i = 1, ..., N , the linear correlation coefficient r is given as where x and y are the mean values of x i 's and y i 's, respectively. The value r lies between −1 and 1. If r = 1, the data points (x i ; y i ) line on a straight line with a positive slope while, for r = −1, on that with a negative slope. To evaluate r values, we use the cor() function in the R software package [28], which is an integrated software package to perform data-scientific analyses. For multiple-linear-regression modeling, we use the lm() function in R. The regression results are investigated in terms of the coefficient of determination R 2 adj (adjusted R 2 ), which is given as where y i 's and f i 's are observed (in the present case, theoretically calculated) values and predicted (by multiple linear regression) values. N and p are the numbers of (y i , f i ) pairs and explanatory variables in the regression model, respectively. R 2 adj = 1 indicates that the model is perfect. In a summary of application of the lm() function in R, the R 2 adj value is given. Principal component analysis is a procedure to construct a set of linearly uncorrelated variables called principal components from a possibly correlated observed variables.

Defect-free states
In figure 1, the calculated S, W and τ for the DF states in Al 0.5 Ga 0.5 N, In 0.5 Ga 0.5 N and Al 0.5 In 0.5 N are shown as (a) S-W and (b) S-τ plots. The results for AlN, GaN and InN at their experimental lattice constants as well as those at the lattice constants of the alloys are also plotted for comparison. For AlN, GaN and InN, it is clearly shown that the S, W and τ values change significantly depending on the lattice constants. For example, the W values are in order of GaN InN > AlN at the experimental lattice constants while the W values for GaN and InN are close to each other at the In 0.5 Ga 0.5 N lattice constants. The actual local environment in an alloy system is thought to be in between the two extremes. Anyway the S, W and τ values for an alloy system can be well approximated by the average of the values for their constituent binary systems at the experimental lattice constants or by the average of those at the alloy lattice constants. The present results also indicate that the epitaxial strain can affect the positron annihilation parameters significantly.   In 0.5 Ga 0.5 N (2.90 × 10 −4 ). In parentheses, the standard deviations are shown. We try to understand the reasons of these differences in terms of local structures of the cation mono-vacancies. As shown in figure 2(a), a cation mono-vacancy V c is surrounded by four N atoms and twelve cations. Three structural parameters that characterize the local environment of the vacancy defect in the There are several common characteristics in the r values among the three compounds. The pairs (S, τ ), (S, d N ), and (τ , d N ) show strong correlations (r > 0.9) although the S and τ distributions differ from compound to compound. Here, we define an additional parameter d DF N , which is defined as the average distance between the cation and four N atoms in the DF systems. There is one-to-one correspondence between d DF N and d N since they are evaluated before and after introduction of a vacancy at a certain site and optimization of atomic positions. Relations between d DF N and the corresponding d N are plotted together with their frequency distributions in figure 5. It is shown that the distribution profiles of S, τ and d N for each compound have some similarity to one another, corresponding to the large r values between these parameters listed in table 1. Thus, the distribution characteristics of S and τ reflect the size distribution of the N tetrahedrons surrounding the cation mono-vacancies. The distributions of d DF N (figure 5) show two peaks reflecting the difference in the cation size. According to Shannon [30], the effective ionic radii of fourcoordinated trivalent Al, Ga and In are 0.39, 0.47 and 0.62 Å , respectively. It is clearly shown that the width of the d N distribution is closely related to the separation between the two peaks in the d DF N distribution. It seems that the d N value does not depend on whether A or B is removed to create the vacancy

Cation mono-vacancies
We cannot tell whether V c is V A or V B from its local configuration although it can be practically determined by counting the number of A or B atoms in the supercell. The center of d N distribution is larger than that of d DF N as shown in figure 5. For cation mono-vacancies in the binary compounds AlN, GaN and InN, computational studies report the outward relaxation of the neighboring N atoms [9]. For In 0.5 Ga 0.5 N, the r values for the pairs (S, x B ), (τ , x B ) and (d N , x B ) are opposite in sign to those for the Al 0.5 Ga 0.5 N and Al 0.5 In 0.5 N cases. As for the cation size, Al < Ga, In > Ga and Al < In. Only for the In 0.5 Ga 0.5 N case, B is smaller. Since it is thought that S, τ and d N are closely related to the defect size, the difference of r sign can be attributed to the size relation between the cations.
The correlation between W and x B is also strong although r = 0.865 (<0.9) for In 0.5 Ga 0.5 N. As mentioned above, the width of W distribution is in order of Al 0.5 Ga 0.5 N > Al 0.5 In 0.5 N In 0.5 Ga 0.5 N. We can interpret these results as follows. Al has no d electrons while Ga and In have localized d electrons. The 3d orbitals of Ga are more localized than the 4d orbitals of In. The more localized electronic state, the wider distribution in the momentum space, which contributes to an increase of W. Al 0.5 Ga 0.5 N consists of Al atoms with no d electrons and Ga with the most localized 3d electrons. Depending on the cation site, the x B value is different. These reasons explain why the distribution of W is the widest ( figure 4) and why there is a strong correlation between W and x B for Al 0.5 Ga 0.5 N. As shown in figure 6, in terms of W, the difference between Al and Ga is the most drastic and that between Al and In is the next. The difference between Ga and In is small since both elements have d electrons. Since the relations of the valenceelectron localization are Al < Ga, In < Ga and Al < In, the r values have the same sign among the three compounds.

Vacancy complexes
We expand our calculations to the vacancy complexes V c -nV N (n 3). Figure 6 represents the relations between the calculated S and W values for the vacancy complexes (the results for the cation mono-vacancies are also shown). For each compound, the distributions are similar to that for the  mono-vacancies although their positions shift toward the right (a higher S). This is quite reasonable since the defect size increases with n. Figure 7 represents the relation between S and τ for vacancy complexes. Good correlations are observed even for In 0.5 Ga 0.5 N. With n increasing, the distributions shift toward higher values of S and τ . Al 0.5 In 0.5 N shows the largest distribution for each n. This is reasonable in light of the ionicradius difference.
Similarly to the mono-vacancy case, we try to analyze the obtained positron annihilation parameters S, W and τ in relation to local structural parameters. To deal with vacancy complexes V c -nV N with various n on equal terms, we introduce new local structural parameters d N and d c , which are measured from the redefined defect center, namely the center of twelve cations. For each case of the cation mono-vacancies, the difference between d N and d N and that between d c and d c are marginal. Correlations between S, W, τ , x B , d N and d c are investigated and the obtained r values are listed in appendix.
The r values concerning d N and d c for n = 0 are very similar to the corresponding values for d N and d c in table 1. With n increasing, the correlation between S and τ remains strong while the correlations for the (S, d N ) and (τ , d N ) pairs weaken. Since the number of N atoms decrease with n increasing, the observed variations are thought to be reasonable.
For each compound and n value, we construct multiplelineal-regression models based on the following equations using all the sets of parameters (64, 256, 384 and 256 sets for n = 0, 1, 2 and 3, respectively): Although, as might be expected from r values, it is revealed that not all the three parameters x B , d N and d c are significant in several regression models by checking the t-statistic and corresponding p -value, all of them are included in the regression models to keep the same footing. In table 2, the resultant adjusted R 2 values are listed. The closer value to unity, the better regression model. For the cation mono-vacancy (n = 0) case, the obtained multiple-linear-regression models are excellent. For the other cases also, the regression models are rather good although, as mentioned above, the correlations between S and d N and between τ and d N weaken with n increasing. It is thought that the three descriptors complement one another. In order to strengthen the above discussion, we perform principal component analyses of the parameters S, W, τ , x B , d N and d c with respect to each compound and n. As mentioned in the previous section, the first principal component reflects the variability of these parameters as mush as possible. The second principal component reflects the remaining variability as much as possible under the orthogonality constraint. Factor loadings, which correspond to correlation values between the original parameters and the principal components, are crucial values to understand the nature of a specific component. In the present study, the cumulative proportion values up to the second component, which describe how much information of the parameters is included in the first and second principal components, run from 83% through 99%. We therefore evaluated factor loadings concerning the first and second principal components and plot in figure 8. The linear correlation coefficients in tables 1, A1-A3 represent the quantitative relationship between each pair of the parameters. With the plots of the factor loadings, we can comprehend the relations among all the parameters visually at a time. It should be noted that the sign of each factor loading is meaningful only in a relative sense. Even if the signs of all the values are inverted concerning a certain principal component axis, the interpretation must be the same. If directions from the origin concerning two parameters are parallel or antiparallel to each other, the correlation between them is strong. Figure 8 enables satistically significant correlation between the objective variables S, W and τ and the descriptors x B , d N and d c to be identified. For example, in the plot for Al 0.5 Ga 0.5 N and n = 0, it is shown that x B is an excellent descriptor for W while d N is an excellent descriptor for τ and a good descriptor for S also. In the plot for In 0.5 Ga 0.5 N and n = 2, it looks that each of x B , d N and d c alone is not a good descriptor for S, W, or τ . By (linearly) combining them, good descriptors are expected to be obtained. The plots for Al 0.5 In 0.5 N are roughly similar to one another. The directions of S, W and τ for Al 0.5 Ga 0.5 N and Al 0.5 In 0.5 N are nearly parallel or antiparallel to the first principal component axis while those for In 0.5 Ga 0.5 N are not. This can be also related to the distinct feature in the S-W plot for In 0.5 Ga 0.5 N.
In actual experiments, except for the case that spectral decomposition is possible in positron lifetime measurements, only average values of the positron annihilation parameters for various annihilation sites are obtained for each sample. Nonetheless, the results shown in figures 6 and 7 are meaningful since there is a possibility that vacancy complexes show preferences for surrounding atoms. Here, we assume that such a preference does not exist and take the average for each compound and n. In figure 9, the S-W relation is plotted together with the results for the DF states. For comparison, the results for AlN, GaN and InN also are plotted. Since there two types of vacancy configurations for n = 1-3 in the wurtzite structure,  Table A1. Linear correlation coefficients (r xy ) between positron-annihilation and structural parameters for vacancy complexes V c -nV N (n = 0-3) in Al 0.5 Ga 0.5 N.   two values are plotted for each of them. Depending on the system, different behaviors are observed with n varying. As already reported in our previous studies [3,4,8,27], with n increasing, the (S, W) values shift toward the right in the S-W plot with a negative slope for AlN or with a positive slope for GaN and InN. The shift directions for the alloy systems are in between the corresponding two of them. In figure 10, the (S, τ ) values are plotted with the results for AlN, GaN and InN.
In contrast with the S-W relation, there are small variations among the systems. A strong correlation between S and τ is clearly shown.

Conclusions
For the DF states as well as the cation mono-vacancy and vacancy complexes V c -nV N (n = 0-3) in Al 0.5 Ga 0.5 N, In 0.5 Ga 0.5 N and Al 0.5 In 0.5 N, the positron annihilation parameters S, W and τ have been calculated for all the possible configurations in the supercells, which contain 128 atoms if there is no vacancy. For the mono-vacancy, the relations between S and W as well as S and τ are presented together with the distributions of S, W and τ . They differ from compound to compound. For example, the S-W plots for Al 0.5 Ga 0.5 N and In 0.5 Ga 0.5 N are very different from each other. We introduce the local structural parameters x B , d N and d c and investigate correlations between the positron annihilation parameters and these local structural parameters. S and τ show strong correlations with d N while W has a strong correlation with x B . The former is understood in terms of the vacancy size while the latter is thought to be related to the localization property of valence electrons. The S-W and S-τ relations are also presented for the vacancy complexes. For each compound, the distributions are similar to that for the mono-vacancies although their positions shift toward a higher S (and also toward a higher τ for S-τ ) with n increasing. Using the local structural parameters x B , d N and d c as the descriptors, multiple-linear-regression models are successfully constructed to describe S, W and τ . By performing principal component analyses, the factor loadings for S, W, τ , x B , d N and d c are obtained. The results are consistent with the other analyses in the present study. The average values of S, W and τ are evaluated and compared with the values for the DF states and those for AlN. GaN and InN in the S-W and S-τ plots, which are expected to be useful in interpreting experimental results.