Stochastic approach to phonon-assisted optical absorption
Marios Zacharias, Christopher E. Patrick, Feliciano Giustino
SStochastic approach to phonon-assisted optical absorption
Marios Zacharias, Christopher E. Patrick, ∗ and Feliciano Giustino Department of Materials, University of Oxford, Parks Road, Oxford OX1 3PH, United Kingdom (Dated: February 21, 2018)We develop a first-principles theory of phonon-assisted optical absorption in semiconductors andinsulators which incorporates the temperature dependence of the electronic structure. We showthat the Hall-Bardeen-Blatt theory of indirect optical absorption and the Allen-Heine theory oftemperature-dependent band structures can be derived from the present formalism by retaining onlyone-phonon processes. We demonstrate this method by calculating the optical absorption coefficientof silicon using an importance sampling Monte Carlo scheme, and we obtain temperature-dependentlineshapes and band gaps in good agreement with experiment. The present approach opens the wayto predictive calculations of the optical properties of solids at finite temperature.
PACS numbers: 78.40.-q, 71.15.-m, 71.38.-k
In semiconductors and insulators exhibiting indirectband gaps the optical transitions near the fundamentaledge require the absorption or emission of phonons in or-der to fulfill the crystal momentum selection rule. Thismechanism is discussed in every introduction to solidstate physics [1, 2]. The theory of phonon-assisted indi-rect optical transitions was developed by Hall, Bardeen,and Blatt (HBB) [3, 4], and forms the basis for our cur-rent understanding of phonon-assisted optical processes.Despite the popularity of the HBB theory, only veryrecently was this formalism combined successfully withfirst-principles density-functional theory calculations [5]powered by Wannier interpolation [6, 7]. The work ofRef. 5 stands as the most sophisticated calculation ofindirect optical absorption available today, yet it is notentirely parameter-free since an empirical shift of the ab-sorption onset at each temperature was needed in orderachieve agreement with experiment. This correction wasunavoidable because the HBB theory does not take intoaccount the temperature dependence of band structures.A consistent theory of temperature-dependent bandstructures was developed by Allen and Heine (AH) [8, 9].In recent years this approach was successfully demon-strated and improved within the framework of first-principles density-functional theory calculations [10–13].Given these recent advances it is natural to ask whetherthe HBB theory of indirect absorption and the AH theoryof temperature-dependent band structures could be com-bined in a more general formalism, in view of fully pre-dictive calculations of phonon-assisted optical processesat finite temperature.In this manuscript we show that the quasiclassicalmethod introduced by Williams [14] and Lax [15] (WL)provides a unified framework for calculating optical ab-sorption spectra of solids, including phonon-assisted ab-sorption and electron-phonon renormalization on thesame footing. Indeed we show that the HBB and AHtheories can be derived from the WL formalism by ne-glecting electron-phonon scattering beyond one-phononprocesses. In order to demonstrate the power of the WL approach we calculate from first principles the phonon-assisted optical absorption spectrum of silicon at differ-ent temperatures using a stochastic importance samplingMonte Carlo method [16] and no adjustable parameters.Our calculations are in very good agreement with exper-imental spectra measured at several temperatures. Wealso calculate temperature-dependent band gaps and findgood agreement with experiment.The premise of the conventional HBB theory is thatelectrons in solids experience a time-dependent poten-tial which arises from the oscillatory motion of theatoms around their equilibrium positions. Followingthis premise, indirect electronic transitions are obtainedwithin time-dependent perturbation theory to first orderin the atomic displacements [4]. This amounts to con-sidering optical transitions whereby the absorption of aphoton is accompanied by the emission or absorption ofone phonon.At variance with the HBB point of view, in the WL ap-proach electrons and phonons are described on the samefooting, and optical excitations correspond to transitionsbetween Born-Oppenheimer product states of electronsand quantum nuclei [14, 15]. The quantized final vi-brational states are then replaced by a classical contin-uum, leading to an expression for the optical absorptionwhich only involves the nuclear wavefunction of the ini-tial state [15]. This replacement can be justified usingthe adiabatic approximation. The temperature depen-dence is then obtained as a canonical average over theinitial states of the system. This theory was successfullyemployed to explain the optical properties of cold lithiumclusters [17] and diamondoids [16].The imaginary part of the WL temperature-dependentdielectric function is given by: (cid:15) ( ω ; T ) = Z − (cid:88) n exp( − E n /k B T ) (cid:104) (cid:15) ( ω ; x ) (cid:105) n , (1)where ω is the photon frequency, k B the Boltzmann con-stant, and T the temperature. Z is the canonical par-tition function among the quantum nuclear states with a r X i v : . [ c ond - m a t . m t r l - s c i ] O c t FIG. 1. (a) The absorption coefficient of bulk silicon at 300 K: calculation with the atoms clamped at their equilibriumpositions (blue dashed line), calculation using the WL method [Eq. (2), blue solid line], and experimental data from Ref. 18(grey filled discs). The thin vertical lines indicate the direct and indirect band gaps with nuclei in their equilibrium positions.(b) Temperature dependence of the absorption coefficient of silicon: WL theory (solid lines) and experimental data for 78 K[19], 300 K [18], and 415 K [20] (grey discs). The calculated spectra were broadened using Gaussians of width 30 meV, andtruncated at the smallest excitation energy in order to avoid artifacts. energies E n , and (cid:104) (cid:105) n stands for the expectation valuetaken over the n -th many-body nuclear wavefunction. (cid:15) ( ω ; x ) denotes the imaginary part of the dielectric func-tion evaluated with the nuclei clamped at the positionsspecified by the set of normal coordinates { x ν } , whichwe indicate collectively as x . In order to keep the nota-tion light we label the normal modes of vibration and theelectronic states by integer indices; accordingly the fol-lowing equations will refer to a Born-von K´arm´an (BvK)supercell of the crystal unit cell. An intuitive interpre-tation of Eq. (1) is that in the adiabatic approximationthe electronic and nuclear timescales are decoupled, andthe measured absorption spectrum is described as an en-semble average over instantaneous absorption spectra atfixed nuclear coordinates. In the harmonic approxima-tion Eq. (1) simplifies via Mehler’s formula [21]: (cid:15) ( ω ; T ) = (cid:90) Π ν dx ν G (cid:2) x ν ; (cid:104) x ν (cid:105) T (cid:3) (cid:15) ( ω ; x ) , (2)where G [ u ; σ ] is a normalized Gaussian of width σ inthe variable u . (cid:104) x ν (cid:105) T = (2 n ν + 1) l ν represents the meansquare nuclear displacement at the temperature T , with n ν the Bose-Einstein occupation factor of the mode withenergy ¯ hω ν at the temperature T , and l ν the correspond-ing zero-point amplitude [22].For simplicity we calculate the dielectric functionwithin the independent-particle approximation, althoughthe present formalism is general and can be used with anydescription of optical transitions at fixed nuclei. In theelectric dipole approximation we have [23]: (cid:15) ( ω ; x ) = 2 πmN ω ω (cid:88) cv | p xcv | δ ( ε xc − ε xv − ¯ hω ) , (3) where m is the electron mass, ω p the plasma frequency, N the number of electrons in the unit cell, and the fac-tor of two is for the spin degeneracy. p xcv is the matrixelement of the momentum operator along the polariza-tion direction of the photon, taken between the valenceand conduction Kohn-Sham states | v x (cid:105) and | c x (cid:105) with en-ergies ε xv and ε xc , respectively. The superscripts indicatethat these states are calculated with the nuclei fixed inthe configuration specified by the normal coordinates x ;the same quantities evaluated at the equilibrium atomicpositions will be denoted without superscripts. Equa-tion (2) was evaluated within density functional theoryusing importance sampling Monte Carlo integration in aBvK supercell, as described below in the Methods.In Fig. 1(a) we compare the optical absorption co-efficient of silicon calculated from first principles usingEqs. (2) and (3) with the experimental spectrum, bothat 300 K. The absorption coefficient was obtained as κ ( ω ; T ) = ω (cid:15) ( ω ; T ) /c n ( ω ) where c is the speed of lightand n ( ω ) the refractive index. The spectrum calculatedwith the nuclei clamped in their equilibrium positions[dashed blue line in Fig. 1(a)] exhibits an onset around3.3 eV, corresponding to the direct Γ (cid:48) v → Γ c tran-sition in silicon. The sub-gap absorption between 1.1–3.3 eV observed in experiments [18] is completely miss-ing in this calculation. At variance with this result,our WL spectrum correctly captures indirect absorption[solid blue line in Fig. 1(a)], and exhibits very good agree-ment with experiment without any adjustable parame-ters. Since we are not including excitonic effects, thestrength of the E transition is underestimated in ourcalculations, as can be seen at energies around 3.3 eV inFig. 1(a) [24]. The agreement between theory and ex-periment in Fig. 1 remarkably extends over five orders ofmagnitude.In order to shed light on the ability of the WL theory tocapture indirect optical absorption we express the depen-dence of the optical matrix elements on the atomic po-sitions using time-independent perturbation theory. Tofirst order in the atomic displacements we have: p xcv = p cv + (cid:88) (cid:48) nν (cid:20) p cn g nvν ε v − ε n + g cnν p nv ε c − ε n (cid:21) x ν l ν , (4)where g mnν = (cid:104) m | ∂V /∂x ν | n (cid:105) l ν is the electron-phononmatrix element associated with the Kohn-Sham potential V , and in the primed summation the terms n = c , v areskipped. The spectral range of indirect absorption cor-responds to photon energies ¯ hω < E dg , with E dg = 3.3 eVbeing the direct band gap of silicon. In this range directoptical transitions are forbidden, therefore from Eq. (4)we have p cv = 0. If we retain only one-phonon processesand neglect the dependence of the electron energies onthe nuclear coordinates, Eqs. (2)-(4) yield: (cid:15) ( ω ; T ) = 2 πmN ω ω (cid:88) cvν (cid:12)(cid:12)(cid:12)(cid:12)(cid:88) (cid:48) n p cn g nvν ε n − ε v + g cnν p nv ε n − ε v − ¯ hω (cid:12)(cid:12)(cid:12)(cid:12) × δ ( ε c − ε v − ¯ hω ) (2 n ν + 1) . (5)This expression is essentially the same as given by theconventional HBB theory of indirect optical absorp-tion [4], and employed in the first-principles calculationsof Ref. 5. The only difference is that the HBB theorycontains phonon energies ± ¯ hω ν in the denominators andthe Dirac delta functions, corresponding to phonon emis-sion and absorption processes, respectively. In the WLapproach these terms are neglected since in the adiabaticapproximation ¯ hω ν (cid:28) ε c − ε v . In Figure S1 [25] we showthat the present result agrees well with the indirect ab-sorption spectrum of silicon calculated using the conven-tional HBB theory in Ref. 5.In Fig. 1(b) we compare our calculated temperaturedependence of the indirect optical absorption lineshapeof silicon with experiment. We focus on the energy range1.1–2.3 eV where the effect of excitonic spectral weighttransfer on the dielectric function is negligible. Our cal-culations are in good agreement with experiment. In par-ticular the theoretical spectra capture both the smoothincrease of the absorption coefficient with temperature,and the concurrent redshift of the absorption onset. Westress that the observed redshift arises naturally in ourcalculations, in contrast with the HBB theory where thiseffect needs to be included empirically [5]. The slight lossof intensity near the indirect edge at the highest tem-perature [spectrum at 415 K in Fig. 1(b)] results fromthe incomplete sampling of multi-phonon processes in ourstochastic approach.In order to understand the effect of temperature inthe WL approach we note that temperature enters theformalism in two ways: firstly in the Bose-Einstein fac- tors (2 n ν + 1) in Eq. (5), as in the conventional HBBtheory. This term mainly modifies the absorption inten-sity. Secondly temperature enters in the electron-phononrenormalization of the electronic band structure, leadingto a temperature-dependent shift of the absorption onset.The latter contribution can be analyzed by rewriting theenergies inside the Dirac delta functions in Eq. (3) usingtime-independent perturbation theory. The result accu-rate to second order in the atomic displacements reads: ε xc = ε c + (cid:88) ν g ccν x ν l ν + (cid:88) µν (cid:20)(cid:88) (cid:48) n g cnµ g ncν ε c − ε n + h cµν (cid:21) x µ x ν l µ l ν , (6)where h nµν = (cid:104) n | ∂ V /∂x µ ∂x ν | n (cid:105) l µ l ν / ε xc in Eq. (6) following the same pre-scription as for the dielectric function in Eq. (2) we obtain(up to third order in the displacements): ε c ( T ) = ε c + (cid:88) ν (cid:20)(cid:88) (cid:48) n | g cnν | ε c − ε n + h cνν (cid:21) (2 n ν + 1) . (7)In the first term inside the square brackets we recog-nize the Fan (or self-energy) electron-phonon renormal-ization; the second term is the Debye-Waller renormal-ization [6, 8–10, 12, 13]. Both terms can be derived froma diagrammatic analysis by considering only one-phononprocesses [27]. Equation (7) represents precisely the AHtheory of temperature-dependent band structures, andexplains the temperature shift of the indirect absorptionlineshapes in Fig. 1(b).From the calculated optical absorption spectra we canextract the temperature dependence of the indirect anddirect band gaps of silicon, following the standard ex-perimental procedure. In fact within the HBB theorythe absorption coefficient near the indirect edge goeslike ω − (¯ hω − E g ± ¯ hω ν ) [4, 23], therefore the indirectgap E g is straightforwardly extracted from a linear fitto ω / κ ( ω ) / . As expected Fig. 2(a) shows that ourcalculated spectra follow a straight line when plotted as ω / κ ( ω ) / . The intercept of this line with the horizon-tal axis yields the indirect band gaps for each tempera-ture, and the results are shown in Fig. 2(b) for two fittingranges, 0–1.5 eV and 0–2 eV. Single-oscillator fits to ourdata using E g ( T ) = E g (0) − a B { / [exp(Θ /T ) − } following Ref. 28 gave a zero-point renormalization of a B = 60–72 meV and an effective temperature Θ = 368–494 K for the two ranges considered. These values are ingood agreement with the experimental data 62 meV and395 K, respectively [29].In Figure S3 [25] we show that the WL spectrum canalso be used to extract the temperature dependence ofthe direct band gap of silicon using standard lineshapeanalysis of second-derivative spectra. Also in this casewe obtain good agreement with experiment. Overallthe agreement between theory and experiment in Fig. 1, FIG. 2. (a) Extraction of the temperature-dependent indirect band gap of silicon using lineshape analysis. The calculated[ ωκ ( ω )] / at each temperature are shown as blue lines, and the corresponding linear fits as thin black lines. The interceptsof the straight lines with the horizontal axis give the band gaps. The linear fits were determined in the energy range 0–2 eV.(b) Temperature-dependent indirect band gap of silicon: the band gaps extracted from the lineshape analysis in (a) usinglinear fits in the ranges 0–2 eV and 0–1.5 eV are shown as blue filled discs and open circles, respectively. Grey filled discs areexperimental data from Ref. 26. The solid lines are single-oscillator fits to the calculated data, and the dashed lines are thecorresponding high-temperature asymptotes. Figure S2 [25] shows the sensitivity of the band gaps to the fitting range. Theshading is a guide to the eye and can be taken as the uncertainty of the theoretical band gaps. Fig. 2, and Figure S3 provides strong support to the va-lidity of the WL theory for first-principles calculations ofphonon-assisted optical absorption spectra.In future work it will be important to test the roleof additional correction terms, such as nonadiabaticity[12], quasiparticle corrections [13], and anharmonicity[30]. While these further refinements will modify the pre-cise values of the zero-point renormalization of the bandgap, it is expected that they will not change any of thefeatures of the lineshapes in Fig. 1.The stochastic approach employed here is remarkablyefficient in sampling the vibrational phase space, to thepoint that the optical spectrum can be calculated using a single configuration of the nuclei (Fig. S4 [25]). This is anunexpected finding and warrants separate investigation.While the present method lacks the elegance of standarddensity-functional perturbation theory approaches [31],it comes with distinctive advantages: (i) the electron-phonon coupling is included to all orders, (ii) the methodcan be used in conjunction with higher-level theories,such as hybrid functionals [32, 33] and the GW/Bethe-Salpeter method [34], and (iii) the anharmonicity of thepotential energy surface can be incorporated by using theappropriate nuclear wavefunctions [35].In conclusion we have demonstrated a new theory ofphonon-assisted optical absorption in solids, based on theWilliams-Lax quasiclassical approximation. This theoryincorporates for the first time the temperature-dependentelectron-phonon renormalization of the electronic struc-ture, and enables calculations of optical spectra at finitetemperature over a wide spectral range. Our stochasticapproach is efficient and easy to implement on top of anyelectronic structure package. The present work opens the way to systematic calculations of optical spectra ofsemiconductors and insulators at finite temperature.
Methods
The calculations were performed within den-sity functional theory in the local density approximation[36, 37], using planewave basis sets and norm-conservingpseudopotentials [38] as implemented in the
QuantumESPRESSO suite [39]. We obtained vibrational frequenciesand eigenmodes via the frozen-phonon method [40, 41].The optical matrix elements including the non-local com-ponents of the pseudopotential [42] were evaluated us-ing
Yambo [43]. Calculations with/without the nonlocalcomponents of the pseudopotential are compared in Fig-ure S1 [25]. In order to address the band gap problemwe used a scissor correction ∆ = 0 .
75 eV in all calcu-lations, close to the GW value of Ref. 44. The non-locality of the scissor operator was taken into accountin the oscillator strengths [42] via the renormalizationfactors ( ε c − ε v ) / ( ε c − ε v + ∆), thereby ensuring thatthe f -sum rule be correctly fulfilled. A comparison be-tween the absorption spectra calculated with or withoutthe scissor correction is shown in Figure S5 [25]. We av-eraged over the atomic configurations using ImportanceSampling Monte Carlo integration [16]. The estimator[45] of ε ( ω ; T ) in Eq. (2) was obtained using configura-tions generated from a random set of normal coordinates { x ν } , as determined from the quantile function of theGaussian distribution, x ν = (2 (cid:104) x ν (cid:105) T ) / erf − (2 t −
1) [46].The t values (one for each normal coordinate, 0 < t < × × ∗ Current address: Center for Atomic-Scale Materials De-sign (CAMD), Department of Physics, Technical Univer-sity of Denmark.[1] C. Kittel,
Introduction to solid state physics , 5th ed. (Wi-ley, New York, 1976).[2] N. W. Ashcroft and N. D. Mermin,
Solid state physics (Harcourt College Publishers, New York, 1976).[3] L. Hall, J. Bardeen, and F. Blatt, Phys. Rev. , 559(1954).[4] F. Bassani and G. Pastori Parravicini, Electronic Statesand Optical Transition in Solids (Pergamon Press, Ox-ford, 1975).[5] J. Noffsinger, E. Kioupakis, C. G. V. de Walle, S. G.Louie, and M. L. Cohen, Phys. Rev. Lett. , 167402(2012).[6] F. Giustino, M. L. Cohen, and S. G. Louie, Phys. Rev.B , 165108 (2007).[7] N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, andD. Vanderbilt, Rev. Mod. Phys. , 1419 (2012).[8] P. B. Allen and V. Heine, J. Phys. C , 2305 (1976).[9] M. Cardona and M. L. W. Thewalt, Rev. Mod. Phys. ,1173 (2005).[10] A. Marini, Phys. Rev. Lett. , 106405 (2008).[11] F. Giustino, S. G. Louie, and M. L. Cohen, Phys. Rev.Lett. , 265501 (2010).[12] E. Cannuccia and A. Marini, Phys. Rev. Lett. ,255501 (2011).[13] G. Antonius, S. Ponc´e, P. Boulanger, M. Cˆot´e, andX. Gonze, Phys. Rev. Lett. , 215501 (2014).[14] F. E. Williams, Phys. Rev. , 281 (1951).[15] M. Lax, J. Chem. Phys. , 1752 (1952).[16] C. E. Patrick and F. Giustino, Nat. Commun. , 2006(2013). [17] F. Della Sala, R. Rousseau, A. G¨orling, and D. Marx,Phys. Rev. Lett. , 183401 (2004).[18] M. A. Green and M. J. Keevers, Prog. Photovolt.: Res.Appl. , 189 (1995).[19] R. Braunstein, M. R. Arnold, and F. Herman, Phys.Rev. , 695 (1958).[20] H. A. Weakliem and D. Redfield, J. Appl. Phys. , 1491(1979).[21] G. N. Watson, J. London Math. Soc. s1-8 , 194 (1933).[22] C. E. Patrick and F. Giustino, J. Phys.: Condens. Matter , 365503 (2014).[23] P. Y. Yu and M. Cardona, Fundamentals of semiconduc-tors , 4th ed. (Springer, Heidelberg, 2010) pp. 200 – 343.[24] L. X. Benedict, E. L. Shirley, and R. B. Bohn, Phys.Rev. B , R9385 (1998).[25] See Supplemental Material at [URL] for SupplementalFigures S1–S7, and Supplemental Reference [49].[26] V. Alex, S. Finkbeiner, and J. Weber, J. Appl. Phys. , 6943 (1996).[27] A. Marini, S. Ponc´e, and X. Gonze, Phys. Rev. B ,224310 (2015).[28] P. Lautenschlager, M. Garriga, L. Vina, and M. Car-dona, Phys. Rev. B , 4821 (1987).[29] M. Cardona, Solid State Commun. , 3 (2005).[30] I. Errea, M. Calandra, and F. Mauri, Phys. Rev. B ,064302 (2014).[31] S. Baroni, S. de Gironcoli, A. Dal Corso, and P. Gian-nozzi, Rev. Mod. Phys. , 515 (2001).[32] J. P. Perdew, M. Ernzerhof, and K. Burke, J. Chem.Phys. , 9982 (1996).[33] J. Heyd, G. E. Scuseria, and M. Ernzerhof, J. Chem.Phys. , 8207 (2003).[34] G. Onida, L. Reining, and A. Rubio, Rev. Mod. Phys. , 601 (2002).[35] B. Monserrat, N. D. Drummond, and R. J. Needs, Phys.Rev. B , 144302 (2013).[36] D. M. Ceperley and B. J. Alder, Phys. Rev. Lett. , 566(1980).[37] J. P. Perdew and A. Zunger, Phys. Rev. B , 5048(1981).[38] M. Fuchs and M. Scheffler, Comput. Phys. Commun. , 67 (1999).[39] P. Giannozzi et al. , J. Phys.: Condens. Matter , 395502(2009).[40] K. Kunc and R. M. Martin, in Ab initio calculation ofphonon spectra , edited by J. Devreese, V. Van Doren,and P. Van Camp (Plenum, New York, 1983) p. 65.[41] G. J. Ackland, M. C. Warren, and S. J. Clark, J. Phys.:Condens. Matter , 7861 (1997).[42] A. F. Starace, Phys. Rev. A , 1242 (1971).[43] A. Marini, C. Hogan, M. Gr¨uning, and D. Varsano, Com-put. Phys. Commun. , 1392 (2009).[44] H. Lambert and F. Giustino, Phys. Rev. B , 075117(2013).[45] R. E. Caflisch, Acta Numer. , 1 (1998).[46] S. E. Brown, I. Georgescu, and V. A. Mandelshtam, J.Chem. Phys. , 044317 (2013).[47] S. Joe and F. Y. Kuo, SIAM J. Sci. Comput. , 2635(2008).[48] F. Aurenhammer, ACM Comput. Surv. , 345 (1991).[49] Y. Okada and Y. Tokumaru, J. Appl. Phys. , 314(1984). Supplemental Material for“Stochastic approach to phonon-assisted optical absorption”
FIG. S1: (a) Comparison between the optical absorption spectrum of silicon calculated in the present work (blue solid line), the spectrumcalculated in Ref. [5] using the HBB theory (red solid line), and the experimental data from Ref. [19] (gray filled discs). All spectra are for T = 78 K. The contributions to the optical matrix elements arising from the nonlocal components of the pseudopotential were not includedin Ref. [5]. For consistency in this figure we calculated spectra without such contributions. (b) and (c) Imaginary part of the dielectricfunction of silicon calculated in the present work at 78 K, with (blue) or without (gray) including the contributions to the optical matrixelements arising from the nonlocal components of the pseudopotential. We note that the scale in (a) is logarithmic, while in (b) and (c)we use a linear scale for clarity. The inclusion of the nonlocal components of the pseudopotential leads to a modification of the oscillatorstrength in the order of 10–15%. We note that in Figs. 1 and 2 of the main text the non-local contributions are correctly included. Thepseudopotential used for these calculations has the local component in the d channel. FIG. S2: Dependence of the calculated indirect band gap of silicon on the energy cutoff E c used in the linear fits of Fig. 2. The fits wereperformed in the energy range [0 , E c ], and the temperature considered in this example is 300 K. FIG. S3: (a) Calculated second-derivative spectra of the real part of the dielectric function of silicon at T = 0 K within the WL formalism.In experiment the direct band gap of silicon is determined by identifying the E (cid:48) transition with the first dip in this spectrum [28].(b) Calculated temperature dependence of the direct gap of silicon within the WL theory (blue filled discs), using the second-derivativemethod illustrated in panel (a). For comparison the experimental data from Ref. [28] are shown as gray filled discs. A single-oscillator fit(see main text) describes adequately the calculated temperature dependence (thin blue line), and we obtain a B = 44 meV and Θ = 337 . ±
17 meV and 267 ±
123 K [28], respectively. The high-temperature asymptote isshown as a straight blue dashed line. Above 350 K the E (cid:48) and E transitions merge (both in experiments and in our calculated spectra)and the direct band gap cannot be determined using this procedure. FIG. S4: (a) Sensitivity of the calculated imaginary part of the dielectric function of silicon to the number of nuclear configurationsused for evaluating the Monte Carlo estimator. The thin blue line and the thick gray line correspond to averages over 6 and 1 randomconfigurations, respectively. (b) Absorption coefficient of silicon calculated using two different Gaussian broadening parameters, 30 meV(blue line) and 1 meV (black line). These calculations are for a 8 × × FIG. S5: Calculated absorption coefficient of bulk silicon at 300 K with (black sold line) and without (blue solid line) the scissor operator.The experimental data are from Ref. [18] (gray filled discs). The thin vertical lines indicate the indirect and direct band gaps (withoutscissor), with nuclei in their equilibrium positions. It is seen that the scissor yields essentially a rigid blue-shift of the entire spectrum, anddoes not introduce any artifacts. We can rationalize this observation as follows. The scissor correction modifies the imaginary part of thedielectric function in two ways: firstly by rigidly shifting the transition energies, and secondly by renormalizing the oscillator strengths.This latter aspect is important since the scissor operator is a projection over the unoccupied manifold, and its non-locality needs to betaken into account in order to obtain the correct dipole matrix element, see Ref. [42]. Let us call the dielectric function without scissor (cid:15) ( ω ), and that with scissor as (cid:15) s ( ω ). Based on the previous considerations, in the independent-particle approximation the imaginary partsare related as follows: (cid:15) s2 ( ω ) = (1 − ∆ / ¯ hω ) (cid:15) ( ω − ∆ / ¯ h ). From this expression it is immediate to verify that the dielectric function afterscissor correction correctly fulfills the same f -sum rule as the original function. Using the above relation we find that the absorptioncoefficients with/without scissor are related as follows: κ s ( ω ) = [ n ( ω − ∆ / ¯ h ) /n s ( ω )] κ ( ω − ∆ / ¯ h ). In the present case the ratio of therefractive indices with/without scissor is found to be approximately equal to unity in the photon energy range of interest (0.91–1.01 forenergies up to 4.5 eV). This reflects the fact that below the direct gap the refractive index is dominated by the real part of the dielectricfunction, and that there is no sharp structure in (cid:15) ( ω ) below the direct gap. As a result the absorption spectrum in the presence of scissorcorrection appears simply as a blue-shifted version of the un-modified spectrum, as shown in the figure. FIG. S6: (a) Numerical convergence of the calculated absorption coefficient of silicon with respect to the size of the BVK supercell (bluelines), compared with the experimental data from Ref. [18] (gray filled discs). Dotted, dash-dotted, and solid lines refer to 2 × ×
2, 4 × × × × × × k -points. FIG. S7: Calculated indirect band gap of silicon at clamped nuclei as a function of temperature. For each temperature we used theexperimental lattice parameters from Ref. [49]. The discs show the change in band gap relative to the value at 100 K, and the linesare guides to the eye. As expected the change between 100 K and 500 K is very small ( < ∼∼