Ab initio calculation of the binding energy of impurities in semiconductors: Application to Si nanowires
aa r X i v : . [ c ond - m a t . m t r l - s c i ] A p r Ab initio calculation of the binding energy of impurities insemiconductors: Application to Si nanowires
Y. M. Niquet, ∗ L. Genovese, C. Delerue, and T. Deutsch CEA-UJF, INAC, SP2M/L Sim, 38054 Grenoble Cedex 9, France European Synchrotron Radiation Facility,6 rue Horowitz, BP 220, 38043 Grenoble, France Institut d’ ´Electronique, de Micro-´electronique et de Nanotechnologie (UMR CNRS 8520),D´epartement ISEN, 41 boulevard Vauban, F-59046 Lille Cedex, France (Dated: October 16, 2018)
Abstract
We discuss the binding energy E b of impurities in semiconductors within density functionaltheory (DFT) and the GW approximation, focusing on donors in nanowires as an example. Weshow that DFT succeeds in the calculation of E b from the Kohn-Sham (KS) hamiltonian of theionized impurity, but fails in the calculation of E b from the KS hamiltonian of the neutral impurity,as it misses most of the interaction of the bound electron with the surface polarization charges ofthe donor. We trace this deficiency back to the lack of screened exchange in the present functionals. E b of donors and acceptors is a key quantity in semiconductor physicsbecause it determines the doping efficiency. In semiconductor nanostructures for example,confinement and electrostatics tend to shift the impurity levels deeper in the gap, whichdecreases dopant activity [1, 2]. Therefore, ab initio calculations of impurity binding energiesare highly desirable to assess the performances of ultimate nanodevices. Besides, donorsand acceptors are the prototypes of charged defects in semiconductors, and a fundamentalunderstanding of the strengths and weaknesses of present ab initio approaches such as densityfunctional theory (DFT) and the GW approximation [3, 4] would open the way to a moreaccurate modeling of complex defects.So far, the calculation of E b in bulk semiconductors has been possible only with semi-empirical methods [5, 6]. However, calculations based on DFT have become practicable inultimate nanostructures with a smaller number of atoms. Recently, the case of donors in Sinanowires (Si NWs) has been adressed with both semi-empirical methods and DFT, withcontradictory results. Tight-binding [7, 8] and effective mass calculations [9], supportedby experiments [1, 2], indeed suggest that E b increases as 1 /R with decreasing wire radius R , due to the interaction of the bound electron with the surface polarization (or “image”)charges of the impurity, resulting in a significant decrease of the doping efficiency in the R <
10 nm range. In contrast, DFT calculations [10, 11] predict that E b decreases muchfaster than 1 /R , and is about 3 − E b correctly in bulk and nanostructures, because they missmost of the interactions of the carriers with the polarization charges of the impurity. Wepropose an alternative strategy based on the KS hamiltonian of the ionized donor whichcircumvents this deficiency.For a donor, E b is the energy needed to ionize the neutral impurity and bring the electronto the conduction band edge far away. It can be defined as the difference E b = I d ( N + 1) − A p ( N ) between the ionization energy I d ( N +1) of the neutral impurity (with N +1 electrons)and the affinity A p ( N ) of the pristine system (with N electrons and no dopant). Since A ( N ) = I ( N + 1), the binding energy can also be computed as an isoelectronic difference ofionization energies, E b = I d ( N +1) − I p ( N +1), or affinities, E b = A d ( N ) − A p ( N ). In practice,the ionization energies and affinities can be computed either as total energy differences I ( N ) = E ( N − − E ( N ) and A ( N ) = E ( N ) − E ( N + 1) [10], or as “quasiparticle” energies211], i.e., as the highest occupied (HOMO) and lowest unoccupied (LUMO) molecular orbitalenergies. However, the quasiparticle problem should in principle be adressed with many-body perturbation theories (MBPTs) such as the GW approximation, since DFT is knownto miss the HOMO-LUMO gap [3, 4, 12, 13]. We actually show hereafter that the abovedefinitions of E b are consistent in the GW approximation, but not in DFT. Using the insightgained from many-body theory, we conjecture that DFT should succeed in the calculationof E b from the KS LUMO of the ionized impurity, but fails in the calculation of E b fromthe KS HOMO of the neutral impurity, due to the lack of explicit screened exchange in thepresent functionals. We support these conclusions with DFT calculations on Si NWs. The binding energy in many-body theory – In MBPT, the quasiparticle energies E n andwave functions ϕ n of the N -electron system are the solutions of the quasiparticle equation: −
12 ∆ r ϕ n ( r ) + v ion ( r ) ϕ n ( r ) + v h ( r ) ϕ n ( r )+ Z d r ′ Σ xc ( r , r ′ , E n ) ϕ n ( r ′ ) = E n ϕ n ( r ) (1)where v ion ( r ) is the ionic potential, v h ( r ) = R d r ′ ρ ( r ′ ) / | r − r ′ | is the Hartree potentialcreated by the ground-state electronic density ρ ( r ), and Σ xc ( r , r ′ , E n ) is the “self-energy”that describes exchange and correlation effects. The ionization energy is I ( N ) = − E N , whilethe affinity is A ( N ) = − E N +1 .The GW approximation has become the reference for the calculation of the band structureof semiconductors [3, 4]. For illustrative purposes, we shall use hereafter the simpler staticCOHSEX form (COulomb Hole and Screened EXchange) of the GW self-energy [3]. Σ xc canthen be split in two parts Σ COH + Σ
SEX :Σ COH ( r , r ′ ) = 12 δW N ( r , r ) δ ( r − r ′ ) (2a)Σ SEX ( r , r ′ ) = − W N ( r , r ′ ) γ ( r , r ′ ) (2b) γ ( r , r ′ ) = X n ′ ∈{ σ } ϕ ∗ n ′ ( r ) ϕ n ′ ( r ′ ) , (2c)where the sum runs over the occupied states with a given spin σ . W N ( r , r ′ ) is the screened Coulomb interaction, i.e. the total potential created at point r ′ by a test unit charge atpoint r [14]. It can also be split in the bare potential v ( r , r ′ ) = 1 / | r − r ′ | created by this testcharge, plus the response δW N ( r , r ′ ) = W N ( r , r ′ ) − v ( r , r ′ ) of the valence electrons. Σ SEX has the same functional form as the Hartree-Fock exchange, but with a screened instead of a3are Coulomb interaction. Σ
COH ( r , r ′ ) describes the interaction of a carrier at point r withthe valence electrons which dynamically respond to its motion.Before adressing the impurity problem, we shall discuss the form of W N ( r , r ′ ) in bulkmaterials and nanowires. In a solid, a test charge q t = +1 at point r attracts valenceelectrons in a small “cloud” around (over ∼ a bond length). This cloud contains a totalcharge q c = − (1 − /κ ), where κ is the static dielectric constant of the material. Theelectrons are actually dragged from the surface of the system, where they leave an oppositepolarization (or “image”) charge q s = − q c . In bulk, these image charges are infinitely faraway, so that the long-range potential created by the test charge is simply W N ( r , r ′ ) ∼ ( q t + q c ) / | r − r ′ | = 1 / ( κ | r − r ′ | ). In a nanowire, however, the electrons are dragged within afew R ’s only from q t , so that the transfer of charges from the surface to the cloud becomesshorter and shorter-ranged with decreasing R . The screening is therefore reduced by q s andthe potential ultimately tends to W N ( r , r ′ ) ∼ / | r − r ′ | when R → R ’s). This simple picture is consistent with classicalelectrostatics (where the surface polarization charges are given by the discontinuity of theelectric field), and fully supported by quantum calculations [15, 16].As discussed previously, the many-body binding energy of an impurity can be computedas E b = A d ( N ) − A p ( N ), the difference between the affinities of the ionized impurity andpristine systems. They fulfill the equation H p,d ( N ) ϕ p,dN +1 = − A p,d ( N ) ϕ p,dN +1 , where H p ( N )and H d ( N ) are the respective quasiparticle hamiltonians: H p,d ( N ) = −
12 ∆ + v p,d ion + v p,d h + Σ p,d SEX + Σ p,d
COH . (3)The physics of the impurity is most easily brought out from the difference between H p ( N )and H d ( N ). On one hand, the extra proton of the ionized impurity is screened by the valenceelectrons through the Hartree potential v d h . Neglecting short range chemical corrections ina first approximation [17], we can therefore write:[ v d ion + v d h ] − [ v p ion + v p h ] ≃ − W dN ( r i , r ) , (4)where r i is the impurity position. On the other hand, we do not expect significant differencesbetween the screened Coulomb interactions W dN and W pN , nor between the one-particle den-sity matrices γ d and γ p , except possibly right around the donor and surface (image charges),on length scales much shorter than the Bohr radius of the impurity. Hence, Σ d SEX ≃ Σ p SEX ,4 d COH ≃ Σ p COH , and: H d ( N ) ≃ H p ( N ) − W N ( r i , r ) . (5) In a first approximation, the quasiparticle hamiltonian of the ionized impurity is the quasi-particle hamiltonian of the pristine system plus the screened Coulomb potential of a unitcharge at the impurity position . This is the usual “hydrogenic model” [5] used in Refs. [7–9]to calculate E b in Si NWs.The electron is therefore bound to the impurity by the screened Coulomb interaction W N ( r i , r ). In bulk silicon, W N ( r i , r ) ∼ / [ κ | r i − r | ] and E b ≃
50 meV. In a nanowire,however, the electron also interacts with the image charges of the donor. Since the totalsurface polarization charge is q s = (1 − /κ ) ≫ /κ , this leads to a large ∝ /R enhancementof E b with decreasing R [7].Let us now compute the binding energy E b = I d ( N + 1) − I p ( N + 1) from the ionizationenergy of the neutral impurity. I d ( N + 1) and I p ( N + 1) fulfill the equation H p,d ( N +1) ϕ p,dN +1 = − I p,d ( N +1) ϕ p,dN +1 , where, as before, H p ( N +1) and H d ( N +1) are the quasiparticlehamiltonians of the ( N +1)-electron pristine and impurity systems. In the latter, the HOMO ϕ dN +1 is the occupied bound state of the impurity. The neutral impurity as a whole nowintroduces a localized perturbation of the pristine system which is screened by the valenceelectrons. We can therefore write:[ v d ion + v d h ] − [ v p ion + v p h ] ≃ − W dN +1 ( r i , r ) + v b ( r ) , (6)where: v b ( r ) = Z d r ′ W dN +1 ( r , r ′ ) | ϕ dN +1 ( r ′ ) | (7)accounts for the screening of the bound state potential. Assuming again that W pN +1 ≃ W dN +1 ,and that the valence band wave functions ϕ , ..., ϕ N are little affected by the neutral impurity,we further get: Σ d SEX ( r , r ′ ) − Σ p SEX ( r , r ′ ) ≃ − W N +1 ( r , r ′ ) × [ ϕ d ∗ N +1 ( r ) ϕ dN +1 ( r ′ ) − ϕ p ∗ N +1 ( r ) ϕ pN +1 ( r ′ )] . (8)The second term can be neglected in bulk and nanowires where ϕ pN +1 is an extended state.The first term cancels v b ( r ) when applied to the HOMO ϕ dN +1 . The effective hamiltonianfor the bound electron therefore reads: H d ( N + 1) ≃ H p ( N + 1) − W N +1 ( r i , r ) . (9)5n principle, I ( N + 1) = A ( N ) and we should have recovered the same equation as before[Eq. (5)]. Here H p ( N ) is however replaced with H p ( N + 1) and W N with W N +1 . Since ϕ pN +1 is an extended state, H p ( N + 1) and H p ( N ) also primarily differ by the substitution W N → W N +1 . The appearance of W N +1 introduces a residual “self-correlation” error inthe GW ionization energies [18], which is however expected to be limited in solids. Wecan therefore conclude that GW provides a consistent description of the binding energies,whether computed from A d ( N ) or I d ( N + 1).This paragraph clearly demonstrates the importance of screened exchange in the calcu-lation of I d ( N + 1). Screened exchange indeed cancels the unphysical screened interactionof the bound electron with itself which arises from v b ( r ) in Eq. (6). H d ( N + 1) is there-fore the hamiltonian of a charged system as expected (the bound electron interacts with N + 1 ionic charges but N electrons). Such spurious self-interactions are a serious issue inself-consistent descriptions of occupied localized states. In this respect, we would like topoint out that the Hartree-Fock (HF) bare exchange Σ x ( r , r ′ ) = − v ( r , r ′ ) γ ( r , r ′ ) does notproperly correct the screened self-interactions appearing in solids. Following the same linesas before, the HF hamiltonian of the HOMO of the neutral impurity can indeed be written H d HF ( N + 1) ≃ H p HF ( N + 1) − W N +1 ( r i , r ) + v sr b ( r ), where: v sr b ( r ) = Z d r ′ [ W N +1 ( r , r ′ ) − v ( r , r ′ )] | ϕ dN +1 ( r ′ ) | . (10) v sr b ( r ) is the spurious potential created by the valence electrons in response to the boundstate density | ϕ dN +1 ( r ) | , i.e. the potential created by a diffuse charge ρ eff ( r ) ≃ (1 − /κ ) ×| ϕ dN +1 ( r ) | plus the opposite surface polarization charge q s = − (1 − /κ ). These surfacepolarization charges balance those embedded in the impurity potential W N +1 ( r i , r ) (equiv-alent, as discussed before, to the potential of a net charge 1 /κ at r i and q s = (1 − /κ ) atthe surface). H d HF ( N + 1) is therefore approximately equal to the hamiltonian of the pristinesystem plus the bare Coulomb potential of a unit charge spread around the impurity (thecharge 1 /κ at the impurity position plus the diffuse charge ρ eff around). As a consequence, ρ eff plays the role in the HF approximation of an effective polarization charge, mislocalizedwithin the scale of the Bohr radius instead of the surface. The relative error on E b should belimited in thin nanowires (where the Bohr radius is comparable to R ), and maximum in bulk.The implications for hybrid functionals in DFT will be discussed in the next paragraph. The binding energy in DFT – We now discuss the binding energy of the donor within6FT. For the sake of simplicity, we first focus on the local density (LDA) and generalizedgradients (GGA) approximations, then address the case of hybrid functionals. In DFT, theground-state density ρ ( r ) of the N -electron system is computed from the eigenstates of theKohn-Sham hamiltonian [19]: −
12 ∆ r ϕ n ( r ) + [ v ion + v h + v xc ]( r ) ϕ n ( r ) = E n ϕ n ( r ) (11)where v xc ( r ) is the exchange-correlation potential. In LDA and GGA, v xc ( r ) ≡ v xc ( ρ ( r )) isa function of the local density ρ ( r ) and of its derivatives. DFT is known to underestimatethe band gap energy of semiconductors [13]. Still, we show below that DFT should succeedin the calculation of the binding energy from the LUMO of the ionized impurity, but thatpresent functionals fail on the neutral impurity.Let us first compute E b = A d ( N ) − A p ( N ) from the LUMOs of the KS hamiltonians H p KS ( N ) and H d KS ( N ). The previous arguments are also valid in DFT: The extra protonof the donor is screened by the valence electrons, so that Eq. (4) still holds. We do not,moreover, expect much differences between v p xc ( r ) and v d xc ( r ), except possibly right aroundthe impurity. Therefore, in a first approximation: H d KS ( N ) ≃ H p KS ( N ) − W N ( r i , r ) . (12)We hence recover the hydrogenic model as before [Eq. (5)]. The KS hamiltonian of theionized impurity thus embeds the same extra physics (with respect to the hamiltonian ofthe pristine system) as the GW approximation: Although the LUMO energies are typicallyunderestimated by DFT, the binding energies computed as the difference between the KSLUMOs of the ionized impurity and pristine systems should be reasonably accurate . Thisonly holds, of course, as long as the binding energy is not too large with respect to the DFTband gap.Let us now compute E b = I d ( N + 1) − I p ( N + 1) from the HOMOs of the KS hamiltonians H p KS ( N + 1) and H d KS ( N + 1). The KS wave function ϕ dN +1 is the occupied bound state ofthe impurity. The neutral impurity as a whole is again screened by the valence electrons[Eq. (6)]. The exchange-correlation potential v xc ( r ) is also affected by the extra bound statedensity around the impurity. We hence get : H d KS ( N + 1) ≃ H p KS ( N + 1) − W N +1 ( r i , r )+ v b ( r ) + ∆ v xc ( r ) , (13)7here ∆ v xc ( r ) = v d xc ( r ) − v p xc ( r ). At variance with the GW approximation, ∆ v xc ( r ), a localdensity correction within the Bohr radius, can not be expected to cancel v b ( r ) [Eq. (7)], along-range Coulomb term. This results in i ) a self-interaction error, and ii ), an almost com-plete cancellation of the interaction of the electron with the image charges of the impurity.Indeed, W N +1 ( r i , r ) and v b ( r ) are the potentials created by two opposite charges (the ionizedimpurity and bound electron), leaving no net charge in the hamiltonian. Both errors giverise to an increase of the impurity level and to a decrease of the binding energy. Althoughthis is especially sensitive in thin nanowires, where the enhancement of E b is mostly due tothe interaction with the image charges, the LDA and GGA would fail up to the bulk [wherethe impurity potential decreases exponentially instead of 1 / ( κ | r − r i | )]. We stress that thecalculation of the ionization energy or affinity of the impurity as a difference of total energies[10], which involves the neutral impurity as the initial or final state, suffers from the samedeficiencies in the LDA or GGA. Application to Si nanowires – The binding energies of dopant impurities in Si NWs havebeen previously computed from the KS HOMO of the neutral impurity using GGA and ahybrid functional (HGH) [11], i.e. a mixture of Hartree-Fock bare exchange with GGA. Asdiscussed previously, bare echange does not localize the polarization charges properly, theerror being however likely limited in thin nanowires (the total charge being correct). TheGGA results of Ref. [11] are therefore expected to completely miss image charge effects, whilethe HGH results, which include 12% bare exchange, are expected to account for ≃
12% ofthe interactions with the image charges (even though mislocalized). As a consequence, thedifference between the GGA and HGH results of Ref. [11] should be approximately 12% ofthe image charge correction given by Eq. (3) of Ref. [7], that is, 0 .
12 eV for R = 1 nm, 0 . R = 0 .
75 nm, and 0 .
25 eV for R = 0 . d = 1 .
73 nm, either as E Ib = I d ( N + 1) − A p ( N ), or as E Ab = A d ( N ) − A p ( N ), using KSHOMOs and LUMOs as ionization energies and affinities. The LDA was used in a waveletbasis set as implemented in the BigDFT code [20]. The neutral impurity was first relaxed ina 660 atoms supercell. Since the treatment of a charged system is still problematic withinsuch a supercell approach, A d ( N ) (as well as I d ( N + 1) and A p ( N ) for consistency) were8ctually computed in finite rods with lengths l up to 9 . E Ab = 0 .
93 eV and E Ib = 0 .
06 eV for l = 9 . E Ab is much larger than E Ib , and in good agreement with the semi-empirical model of Ref. [7]( E b = 0 .
92 eV when l → ∞ ). This confirms that present functionals are able to predict thebinding energies of impurities or defects from the KS hamiltonian of the charged defect .To conclude, we have shown, by a formal comparison with the GW approximation, thatthe donor binding energies computed from the Kohn-Sham hamiltonians of neutral impu-rities can be strongly underestimated in semiconductor nanostructures (even with hybridfunctionals). This is due to the lack of screened exchange in the present functionals, andexplains the discrepancies between Refs. [10, 11] and previous works [7, 9]. The bindingenergy of a donor should preferably be computed as the difference between the Kohn-ShamLUMOs of the ionized impurity and pristine systems. This provides a reasonable substitutefor much more expensive GW calculations of defect bound states in solids.We thank L. Wirtz, X. Blase and H. Mera for fruitful discussions. This work was sup-ported by the french ANR project “QuantaMonde” (contract ANR-07-NANO-023-02). Thecalculations were run at the CCRT and CINES. ∗ Electronic address: [email protected][1] M. T. Bj¨ork et al ., Nature Nanotechnology , 103 (2009).[2] J. Yoon et al ., Appl. Phys. Lett. , 142102 (2009).[3] L. Hedin and S. Lundqvist, Solid State Physics , ed. by H. Ehrenreich, F. Seitz, andD. Turnbull (Academic Press, New York, London 1969).[4] G. Onida, L. Reining, and A. Rubio, Rev. Mod. Phys. , 601 (2002).[5] W. Kohn and J. M. Luttinger, Phys. Rev. , 915 (1955).[6] A. S. Martins et al ., Phys. Rev. B , 245205 (2002).[7] M. Diarra et al ., Phys. Rev. B , 045301 (2007).[8] M. Diarra et al ., J. Appl. Phys. , 073703 (2008).[9] B. Li et al ., Phys. Rev. B , 115335 (2008).
10] C. R. Leao, A. Fazzio, and A. J. R. da Silva, Nano Lett. , 1866 (2008).[11] R. Rurali et al ., Phys. Rev. B , 115303 (2009).[12] F. Bruneval, Phys. Rev. Lett. , 176403 (2009).[13] M. Lannoo, M. Schl¨uter, and L. J. Sham, Phys. Rev. B , 3890 (1985).[14] More precisely, the screened coulomb interaction in the random phase approximation (RPA).[15] C. Delerue, M. Lannoo, and G. Allan, Phys. Rev. B , 115411 (2003).[16] F. Trani et al ., Phys. Rev. B , 245430 (2006).[17] We also neglect non-linear and exchange-correlation effects beyond the RPA in Eq. (4). Theywould not change the conclusions drawn in this letter.[18] P. Romaniello, S. Guyot and L. Reining, J. Chem. Phys. , 154111 (2009).[19] W. Kohn and L. J. Sham, Phys. Rev. , 1133 (1965).[20] L. Genovese et al ., J. Chem. Phys. , 014109 (2008), 014109 (2008)