An ab initio based approach to optical properties of semiconductor heterostructures
L. C. Bannow, P. Rosenow, P. Springer, E. W. Fischer, J. Hader, J. V., Moloney, R. Tonner, S. W. Koch

TL;DR
This paper introduces an ab initio method combining density functional theory and the envelope function approach to predict optical properties of semiconductor heterostructures, validated on a quantum well system.
Contribution
It develops a novel procedure integrating DFT and kp-models to accurately predict optical responses of semiconductor heterostructures from first principles.
Findings
Luttinger parameters from TB09 DFT match extrapolated values.
The method accurately predicts absorption spectra for (AlGa)As/Ga(AsP) quantum wells.
The approach enables ab initio design of semiconductor optical devices.
Abstract
A procedure is presented that combines density functional theory computations of bulk semiconductor alloys with the semiconductor Bloch equations, in order to achieve an ab initio based prediction of the optical properties of semiconductor alloy heterostructures. The parameters of an eight-band kp-Hamiltonian are fitted to the effective band structure of an appropriate alloy. The envelope function approach is applied to model the quantum well using the kp-wave functions and eigenvalues as starting point for calculating the optical properties of the heterostructure. It is shown that Luttinger parameters derived from band structures computed with the TB09 density functional reproduce extrapolated values. The procedure is illustrated by computing the absorption spectra for a (AlGa)As/Ga(AsP)/(AlGa)As quantum well system with varying phosphide content in the active layer.
Click any figure to enlarge with its caption.
Figure 1
Figure 2
Figure 3
Figure 4
Figure 5
Figure 6
Figure 7
Figure 8| (eV) | (eV) | (eV) | |||||
|---|---|---|---|---|---|---|---|
| PBE0 | 1.71 | 0.38 | 4.238 | 0.014 | 0.003 | 0.301 | 59.546 |
| HSE06 | 1.11 | 0.37 | 4.456 | 0.014 | 0.003 | 0.293 | 50.541 |
| TB09 | 1.44 | 0.32 | 1.119 | 0.657 | 0.017 | 0.828 | 19.876 |
| TB09a | 1.52 | 0.32 | 1.130 | 0.656 | 0.018 | 0.813 | 20.661 |
| Ref. [46] | 1.52 | 0.34 | 1.13 | 0.759 | 0.0405 | 0.538 | 26.06 |
| (Å) | ||||||||
|---|---|---|---|---|---|---|---|---|
| GaAs | 1.44 | 0.32 | 1.087 | 0.673 | 0.002 | 0.876 | 20.022 | 5.68942 |
| GaP | 2.95 | 0.08 | 1.150 | 0.479 | 0.109 | 0.628 | 19.875 | 5.47747 |
| Ga(As26P1) | 1.48 | 0.31 | 1.167 | 0.633 | 0.001 | 0.877 | 19.884 | 5.68157 |
| Ga(As25P2) | 1.52 | 0.30 | 1.021 | 0.649 | 0.001 | 0.994 | 19.876 | 5.67371 |
| Ga(As24P3) | 1.59 | 0.29 | 0.919 | 0.755 | 0.016 | 1.055 | 20.025 | 5.66586 |
| Ga(As23P4) | 1.62 | 0.28 | 1.124 | 0.639 | 0.055 | 0.750 | 19.817 | 5.65810 |
| Ga(As22P5) | 1.68 | 0.28 | 0.972 | 0.690 | 0.007 | 1.020 | 20.229 | 5.65016 |
| Ga(As21P6) | 1.72 | 0.28 | 1.186 | 0.478 | 0.000 | 0.552 | 19.026 | 5.64230 |
| (Al9Ga18)As | 1.86 | 0.30 | 0.792 | 0.640 | 0.001 | 1.046 | 19.492 | 5.68219 |
| Scissor(b) | ||||||||
|---|---|---|---|---|---|---|---|---|
| GaAs | 1.42 | 0.02 | 0.32 | 1.087 | 0.672 | 0.003 | 0.877 | 19.805 |
| GaP | 2.78 | 0.17 | 0.08 | 1.048 | 0.530 | 0.059 | 0.781 | 19.657 |
| Ga27(As26P1) | 1.47 | 0.01 | 0.31 | 1.168 | 0.632 | 0.001 | 0.877 | 19.780 |
| Ga27(As25P2) | 1.51 | 0.01 | 0.30 | 1.018 | 0.649 | 0.001 | 0.998 | 19.604 |
| Ga27(As24P3) | 1.55 | 0.04 | 0.29 | 0.929 | 0.750 | 0.023 | 1.040 | 19.594 |
| Ga27(As23P4) | 1.60 | 0.02 | 0.28 | 1.127 | 0.637 | 0.058 | 0.745 | 19.579 |
| Ga27(As22P5) | 1.64 | 0.04 | 0.28 | 0.965 | 0.692 | 0.004 | 1.030 | 19.899 |
| Ga27(As21P6) | 1.69 | 0.03 | 0.28 | 1.184 | 0.478 | 0.000 | 0.554 | 18.800 |
| (Al9Ga18)As27 | 1.88 | 0.02 | 0.30 | 0.792 | 0.641 | 0.001 | 1.046 | 19.674 |
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
An ab initio based approach to optical properties of semiconductor heterostructures
L. C. Bannow1, P. Rosenow2, P. Springer1111Current address: Institute of Technical Physics, German Aerospace Center, Paffenwaldring 38-40, 70569 Stuttgart, Germany, E. W. Fischer2, J. Hader3, J. V. Moloney3, R. Tonner2, and S. W. Koch1
1 Department of Physics and Material Sciences Center, Philipps-Universität Marburg, Renthof 5, 35032 Marburg, Germany
2 Department of Chemistry and Material Sciences Center, Philipps-Universität Marburg, Hans-Meerwein-Str. 4, 35032 Marburg, Germany
3 Nonlinear Control Strategies Inc, 7040 N. Montecatina Dr., Tucson, AZ 85704, USA and College of Optical Sciences, University of Arizona, Tucson AZ 85721, USA
Abstract
A procedure is presented that combines density functional theory computations of bulk semiconductor alloys with the semiconductor Bloch equations, in order to achieve an ab initio based prediction of the optical properties of semiconductor alloy heterostructures. The parameters of an eight-band -Hamiltonian are fitted to the effective band structure of an appropriate alloy. The envelope function approach is applied to model the quantum well using the -wave functions and eigenvalues as starting point for calculating the optical properties of the heterostructure. It is shown that Luttinger parameters derived from band structures computed with the TB09 density functional reproduce extrapolated values. The procedure is illustrated by computing the absorption spectra for a (AlGa)As/Ga(AsP)/(AlGa)As quantum well system with varying phosphide content in the active layer.
††: Modelling Simulation Mater. Sci. Eng.
- March 2017
Keywords: density functional theory, semiconductors, -theory, semiconductor Bloch equations, absorption
1 Introduction
The majority of modern electronic devices are based on semiconductors. Therefore, a lot of research aims to improve their optical and electronic properties. Band gap engineering by means of semiconductor alloying has been proven a powerful tool for this purpose. Today, up to five different semiconductors are mixed resulting in quinary alloys, e.g. (GaIn)(NAsSb) [1, 2, 3, 4, 5] or (AlInGa)(AsSb) [6, 7, 8]. In addition, modern growth techniques allow the production of very pure compounds which is essential for optical properties. With this ever increasing number of available materials, it becomes increasingly important to predict the optical properties of new compounds to determine their usability for application in opto-electronic devices such as semiconductor lasers or solar cells.
In this work, we present a method to calculate the optical properties of III-V heterostructures based on a first principles approach. For this purpose, the band structure of the semiconductor heterostructure is obtained within -theory [9] and the envelope function approach [10] as described in section 2.2. However, this approach heavily depends on material parameters such as effective masses which ultimately need to be extracted from experiments. To cope with this issue, we use density functional theory (DFT) [11] to calculate the bulk band structures of each quantum well (QW) within the sample. To overcome the common shortcoming of DFT to underestimate the band gap of semiconductors [12], an advanced functional is used. A bulk -band structure [9] is then fitted to its DFT counterpart. This approach results in a set of effective material parameters which are used to obtain the band structure of the QW sample within the envelope function approach.
The resulting energy bands and wave functions from the -calculation serve as starting point for the calculation of the optical properties. We calculate the absorption using the semiconductor Bloch approach [13] as outlined in section 2.4. The computation of other quantities such as the refractive index or photo luminescence is also possible within our microscopic approach. However, in this work we only aim at demonstrating the possibility of combining DFT and the semiconductor Bloch approach. Therefore, we chose a Ga(AsP) QW with (AlGa)As barriers, a well known III-V material system, for demonstration purposes.
It is worth noting that the -Hamiltonian applied in this work can be extended to include conduction band anti-crossing (CBAC) allowing the description of the band structure of dilute nitrides [14]. Furthermore, it can be extended to include valence band anti-crossing (VBAC) that allows to describe the band structure of alloys such as Ga(AsBi) [15] or even both CBAC and VBAC for Ga(NAsBi) resulting in a Hamiltonian [16] (10 valence and 4 conduction bands).
2 Methods
2.1 DFT calculations
All density functional theory calculations were performed with the Vienna ab initio Simulation Package (VASP) [17, 18, 19, 20]. The lattice constants of ternary alloys were calculated with Vegard’s rule based on computed lattice constants of the binary constituents at the accuracy described below. These alloy lattice constants do not change during relaxation. Atomic positions were relaxed at the generalized gradient approximation level (GGA) with the PBE functional [21, 22], including Grimme’s D3(BJ) dispersion correction [23, 24]. The relaxation was performed with a SCF convergence of , while forces were converged to . Band structures were computed with the TB09 functional (also known as mBJLDA) which combines a meta-GGA exchange part with LDA correlation [25]. Here, eigenvalues were converged to . A plane wave basis set with a kinetic energy cut-off of was used in conjunction with the projector augmented wave (PAW) method [26, 27]. Sampling of -space was performed with a -centered Monkhorst-Pack grid [28] using six points per direction for primitive zinc-blende type unit cells and two points per direction for 54 atom supercells. For band structure calculations spin-orbit coupling (SOC) was enabled and usage of all symmetries, including -space, was turned off.
Supercell (SC) band structures were projected onto the primitive cell (PC) Brillouin zone (“unfolded”) with BandUP, yielding an effective band structure (EBS) [29, 30]. The EBS assigns a spectral weight to PC-wavevectors and energies. This corresponds to the number of states at the given point and is a measure for the Bloch character. Since the translational symmetry of the primitive lattice is broken by substitutions and relaxations in the SC, the spectral weight can take non-integer values. The energy was resolved to in the unfolding procedure.
As SCs we always used special quasirandom structures (SQS) [31] with 54 atoms in this work. The SQS were computed with the Alloy-theoretic automated toolkit (ATAT) [32, 33, 34] where we set the pair length to include third nearest neighbors, the triplet length to include second nearest neighbors, and the quadruplet length to include nearest neighbors.
2.2 calculations
Starting point for our calculations is the fully coupled bulk -Hamiltonian [35]. Included are the twice spin-degenerate -like conduction band, -like heavy-hole, light-hole, and split-off bands as well as the influence of remote bands via Löwdin renormalization of the Luttinger parameters , and [36, 37]. The set of renormalized Luttinger parameters and the Kane energy for bulk materials are obtained from DFT as described in the next section. Since we are dealing with heterostructures, we employ the envelope function approach [10, 38] to account for the quantum confinement of the charge carriers and include epitaxial strain using the Pikus-Bir formalism [39, 40]. Diagonalizing the resulting Hamiltonian yields the energy dispersion and eigenvectors which are used to obtain the dipole and Coulomb matrix elements. All these are required as input for the calculation of the optical properties (see section 2.4).
2.3 Extraction of -parameters from DFT-band structures
In contrast to the spectral weight of the EBS, diagonalizing the -Hamiltonian will yield a set of discrete bands following a dispersion relation . Thus, the EBS must first be reduced to a band structure of the same form. This is carried out in two steps. First, the EBS is averaged over in order to reduce the smearing of the effective bands. Averaging is performed in a small energy window. Since smearing tends to increase with larger distance to for many systems, the averaging can and sometimes needs to be performed over smaller windows close to and over larger windows further away. This leads to an averaged band structure as an intermediate which is less smeared out than the EBS, but still contains spectral weights with generally non-integer values. In the second step, the spectral weights of the averaged band structure are rounded to the next nearest integer and a band structure with an appropriate number of bands at a given wave-vector and energy is obtained. The relevant number of conduction and valence bands constitute the reduced band structure. This band structure has the desired form and can easily be compared to a calculation. A caveat should be mentioned: this procedure requires that the averaged spectral weight is sufficiently close to an integer value, meaning that the Bloch character is not reduced too strong. For example, in the treatment of dilute nitrides and bismides the described procedure needs to be modified because the band anti-crossing leads to sub-band states that can have a small spectral weight.
In order to compare the reduced band structure with its analog, a measure for the deviation between them must be established. This is achieved by summing over the squared difference of energies between both band structures, weighted by the distance to the -point:
[TABLE]
The weighting ensures that momenta close to the center of the Brillouin zone influence the outcome of the fit stronger resulting in a good agreement of the fit in this region. This is desirable since these -points are more important for optical properties and are better described by -theory than -points further away from the center. Only -points within a sphere with radius are taken into account for the fitting procedure. In this work, .
The actual fit is performed by varying the starting parameters in a two-step procedure in order to reduce . In the first step, the parameters are varied by a fixed relative amount (here: of the most recent value) until the difference of to the last step is smaller than (with being for the initial parameters). A variation is accepted only if it leads to a smaller . It turned out that this step is more stable if not all parameters are varied at the same time. Thus, at first only , , and are varied, then , and finally . Especially the last two parameters showed to be unstable if varied simultaneously.
In the second step, all parameters are varied by small random numbers. Again, new values are only accepted if they lead to a smaller . This is performed for 2000 variations of each parameter which showed to lead to a well converged parameter set. Also, no instability of the solution was observed when all parameters were simultaneously changed during this second step.
2.4 Calculation of optical properties
The absorption spectra are calculated using a microscopic theory which is based on the multi-band semiconductor Bloch equations [41, 42, 43, 13]. As input parameters single particle properties are obtained as described in section 2.2. The equation of motion for the microscopic interband polarization between an electron in band and a hole in band is
[TABLE]
where are the renormalized electron (superscript ) and hole (superscript ) energies, is the electron (hole) density, and is the renormalized field. The renormalization is due to the Coulomb interaction of charge carriers. Microscopic electron-electron and electron-phonon scattering that cause a dephasing of the polarization are included via p^{\lambda\nu}_{\boldsymbol{k}}\big{|}_{\textnormal{corr}}. We treat these correlation terms in second order Born-Markov approximation [13]. The equation of motion for the polarization couples to the equations of motion for the electron and hole densities. From the microscopic polarization the absorption and refractive index can be calculated. For a more detailed explanation of the terms and methods, the reader is referred to [42, 43, 13]. In this work, all absorption calculations are performed at K.
3 Results and Discussion
3.1 Method validation
The computational approach was tested on the well studied system GaAs, since it offers reliable reference values and can be computed at low computational cost. Furthermore, it is the host for most alloys of interest. The band structure was computed with the global hybrid functional PBE0 [44], the range-separated hybrid functional HSE06 [45], and the TB09 functional [25]. The resulting data are listed in table 1. The band gap as an experimentally accessible parameter is well reproduced by TB09 with compared to an experimental value of at K [46]. In contrast, PBE0 yields an overestimated band gap of , while the range-separated HSE06 underestimates the gap stronger than TB09 (). A G0W0 calculation of the band gap yielded , showing that TB09 can reproduce features of much more demanding many-body methods, in accordance with other reported findings. Likewise, the -point band gap of the minority component GaP could be reproduced with similar accuracy ( with TB09, exp. [46]).
Applying the fitting routine to the computed band structure, the hybrid functionals deviate significantly – and consistently with each other – from the reference values for the -parameters, indicating that neither reproduces the band curvature sufficiently well. TB09 performs better overall, showing its most significant deviations for the and parameters. The latter directly influences the curvature of the conduction bands which are notoriously difficult to describe with DFT, while the former shows up in matrix elements coupling valence states with each other. The curvature of all bands is influenced rather strongly by the Kane energy , which may compensate the deviation of other parameters.
Shifting the conduction band to the experimental band gap by applying a scissor operator prior to fitting only has a small influence on the -parameters. Overall, the agreement to the reference values is improved slightly. Applying a scissor operator at this point may be desirable, since the band gap directly controls optical properties like the onset of absorption and the photo luminescence wavelength. The agreement between “scissored” and “unscissored” band parameters shows that this step can be applied without changing the result of the fit.
3.2 Determination of -parameters for alloys
For all DFT band structure calculations of ternary compounds we used SQS. The EBS of Ga27(As26P1) () and of Ga27(As21P6) (), which are the extremes of the concentration range studied here, are shown in figure 1. In both cases, a strong Bloch character is retained. Some breaking of the translational symmetry is reflected by the “smearing” of the bands. This increases with the P concentration due to enhanced symmetry breaking. Especially in the vicinity of the -point, which is most important for the optical properties, the smearing is small and a high Bloch character prevails. Visible smearing occurs at the Brillouin zone edge and especially at the folding lines, i. e. at and of the path to the edge.
Averaging over these EBS and reducing them to eight effective bands is straightforward. The reduced band structures are then used for fitting to extract the -parameters. The resulting parameters reproduce the DFT-band structure within the fitting range well (see figure 2) and are listed in table 2. While the band gap increases consistently with the P concentration (and the spin-orbit splitting decreases), no clear trend can be assigned to any of the other parameters.
The DFT calculations are at K whereas we want to calculate the optical properties at K. In order to obtain band structures at K we apply a scissor operator to shift the band gaps to values obtained from a temperature dependent virtual crystal model [46] before fitting. This only weakly affects the band parameters obtained by fitting, keeping the lack of a clear trend and the strong scattering of values (see table 3). Again, the scissor operator can be safely applied without changing the characteristics of the -bands.
An alternative approach to the virtual crystal model is to calculate the band gaps of the ternary alloys using DFT. The only information we require is the experimental band gap of GaAs at the respective temperature. By variation of the parameter in the TB09 functional [25] as done by Kim et al [47] it is possible to accurately fit the GaAs band gap. When using the original TB09 implementation, is calculated from the electron density and its gradient [25]. Here, we calculate the band gap of GaAs with the original and then find a that reproduces the experimental GaAs band gap. We then use the difference as a correction in the calculation of the band gap of the ternary compound Ga(AsP): First, we calculate the band gap with the original TB09 implementation and then we correct the respective value of by . The band gaps calculated with the corrected reproduce the band gaps from the virtual crystal model within a range of meV for all concentrations considered here.
The EBS and fitting results of the barrier material (Al9Ga18)As27 are shown in figure 3. While the valence band states in figure 3(b) mostly retain a strong Bloch character similar to the dilute phosphides, the conduction band states are strongly smeared, an exception being the lowest ones close to the -point. This is exactly the region required for fitting, so averaging the EBS and reducing it to eight bands was again straightforward. The resulting parameters are shown in table 2 and table 3 without and with scissor shift, respectively. As for the other cases described in this work, the -parameters are hardly affected by the scissor shift.
In addition to the -parameters obtained from fitting the ternary Ga(AsP)-SQS band structures, we calculated -parameters for Ga(AsP) by linearly interpolating those of the binary compounds GaAs and GaP from table 3. Interpolation was performed for the Kane energies and un-renormalized Luttinger parameters . The latter are connected to the renormalized ones listed in table 3 via (cf. [36, 37])
[TABLE]
The comparison of the un-renormalized Luttinger parameters and the Kane energies is shown in figure 4. It can be seen that the linearly interpolated -parameters generally overestimate the SQS -parameters. However, a very good agreement is obtained for %.
3.3 Quantum well optical properties
To demonstrate our approach, we calculate the absorption spectra of a (AlGa)As/Ga(AsP)/(AlGa)As QW structure. The barrier material is chosen to be Al0.33Ga0.67As while for the nm QW material we consider six different Ga(AsP) compositions by varying the P content from % to %. From the parameter sets obtained by fitting the bulk DFT-band structures (see table 3) the energy dispersion and wave functions for the eight spin degenerate bands of the QW are calculated (see section 2.2). These serve as input for the computation of the absorption spectrum (see section 2.4) at K.
For each phosphide composition two sets of absorption calculations were performed. We used the -parameters of the Ga(AsP) SQS from table 3 and we used interpolated -parameters that were obtained as described in the last paragraph of section 3.2 to set up the -Hamiltonian.
In figure 5 we present the comparison between the absorption spectra obtained from calculations with linearly interpolated -parameters (lines) and calculations using SQS -parameters (grey areas). Overall, the two sets of calculations deliver comparable results and the shape of the absorption spectra of both sets of calculations coincides. However, in case of %, %, %, and % the spectra from the interpolated calculations are shifted slightly to higher energies by about meV and the peak absorption in the spectra from the SQS calculations is larger. The good agreement of the absorption spectra between both sets of calculations at % and % is consistent with a good agreement of the input Luttinger parameters (see figure 4).
In general, good agreements are expected for very dilute concentrations. Higher concentrations increase the disorder which results in a smearing of states in the EBS. This in turn leads to a higher inaccuracy of the fitted -parameters.
4 Conclusion
We have presented a new method to calculate optical properties of semiconductor heterostructures from first principles using a combination of -theory, density functional theory, and the semiconductor Bloch approach. We demonstrate the applicability of the method by calculating the absorption spectra of (AlGa)As/Ga(AsP)(AlGa)As QWs where the phosphide concentration is varied between .
The only external data we use is the temperature dependent band gap of the alloys. However, we also show that the band gap of Ga(AsP) at the concentrations considered here can be calculated with DFT very accurately using the TB09 functional.
We also show that by interpolation of the Luttinger parameters from the binary compounds GaAs and GaP the absorption spectra obtained are comparable to those calculated for the truly ternary GaAs1-xPx compounds. Both spectra match very well for small concentrations of phosphide. Adding more phosphide increases disorder effects and thus the smearing of the effective band structure.
Our method can, in principle, also be used to describe more complex compounds, e.g. quaternaries, dilute nitrides, and dilute bismides. The difficulty lies in fitting a disordered DFT band structure. For example, the incorporation of nitride strongly perturbs the conduction band while the incorporation of bismide into a host material mainly influences the valence bands. However, at dilute concentrations the distortion is small so a treatment within the here proposed method seems possible.
This work was funded by the DFG via the GRK 1782 “Functionalization of Semiconductors”. LB, PR, EF, RT, and SWK thank the HRZ Marburg, CSC Frankfurt, and HLRS Stuttgart for providing computing time. The work at NLCSTR was supported by the Air Force Office of Scientific Research under the STR Phase II grant FA9550-16-C-0021. LCB and PR contributed equally to this work.
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1[1] Li L H, Sallet V, Patriarche G, Largeau L, Bouchoule S, Merghem K, Travers L and Harmand J C 2003 Electron. Lett. 39 519–520
- 2[2] Yuen H B, Bank S R, Wistey M A, Jr J S H and Moto A 2004 J. Appl. Phys. 96 6375–6381
- 3[3] Harris Jr J S 2005 J. Cryst. Growth 278 3–17
- 4[4] Jackrel D B, Bank S R, Yuen H B, Wistey M A, Jr J S H, Ptak A J, Johnston S W, Friedman D J and Kurtz S R 2007 J. Appl. Phys. 101 114916
- 5[5] Harris J S, Kudrawiec R, Yuen H B, Bank S R, Bae H P, Wistey M A, Jackrel D, Pickett E R, Sarmiento T, Goddard L L, Lordi V and Gugov T 2007 phys. stat. sol. (b) 244 2707–2729
- 6[6] Bugge R and Fimland B O 2006 Phys. Scr. 2006 15
- 7[7] Hosoda T, Belenky G, Shterengas L, Kipshidze G and Kisin M V 2008 Appl. Phys. Lett. 92 091106
- 8[8] Suchalkin S, Jung S, Kipshidze G, Shterengas L, Hosoda T, Westerfeld D, Snyder D and Belenky G 2008 Appl. Phys. Lett. 93 081107
