A brief and bleak historical picture

Considering that detailed studies of the charge transport properties of Si have been published in the 1950s [1] [2] . it may come as a surprise that quantitative understanding of even the gross-scale features of electron scattering in Si has been embarrassingly incomplete until very recently. Let's consider the ideal simple case of an undoped, infinitely long Si crystal. In this case, electrons are subject to two main types of collisions: With the thermal vibrations of the Si ions away from their equilibrium positions (phonon scattering), and with electrons in the valence band. The latter process, known as "impact ionization" or "pair production", consists of an energetic electron "hitting" an electron in the valence band, exciting this electron across the band-gap (about 1.1 eV in Si at 300K), thus generating an electron-hole pair and recoiling to a lower energy state. Quantitative understanding of electron transport demands that we know how many collisions per unit time an electron of energy E will suffer. A look at various results reported in the recent past, shown in Fig. 1, illustrates the bleak state of affairs: Various groups disagree on the magnitude of the electron-phonon scattering rates (top frame) as soon as energies above a few tenth of electron-Volt (eV) are considered, while published pair production rates (bottom frame) vary over 3 orders of magnitude of more.

Lack of interest in high-energy transport in the early days of electronics explains only partially the situation: The evaluation of the scattering rates is rendered difficult by genuine, daunting physical and numerical problems.

Electron-phonon scattering rate

As far as the electron-phonon scattering rates are concerned, we must evaluate how the (Coulomb) potential felt by an electron in the conduction band is modified whenever an ion is displaced from its equilibrium, zero-temperature position. Several approximations must be made at the outset of the calculation:

• We must assume that the ions do not move 'too much' from their equilibrium position. Thus, we can expand the lattice potential in a power series on the displaced coordinates of the ion, note that the linear term vanishes, as we expand around an energy minimum (the equilibrium position), and we can retain only the next term in the expansion. This is a quadratic term, yielding the potential of the harmonic oscillator. This is called the "harmonic approximation". Higher-order terms, anharmonic by nature, lead to multi-phonon processes.
• We must assume that the electrons surrounding the 'displaced' ion follow its motion adiabatically, remaining always in equilibrium. This is known as the "Born-Oppenheimer" approximation.
• More questionably, we must employ ionic "pseudopotentials", instead of the more accurate screened potentials.
• Even more questionably, we must assume that the ionic (pseudo)potential remains unchanged, as the ion is displaced, just shifting rigidly. This is called the "rigid-(pseudo)-ion" approximation.
• Finally, we shall assume that the electrons couple weakly (linearly) with the phonons, so that we can use first-order perturbation theory and the famous "Fermi Golden Rule".

Despite these simplifying approximations, the problem is still very complex. Let's write the rate at which an electron of quasi-momentum k in band n and of energy $$E_n(k) emits/absorbs a phonon as:

where the electron-phonon matrix element $$M_ep(k' n' ; k n ) is given by:

where M is the mass of the unit cell, N the number of unit cells, and $$N_q is the phonon occupation number. As a rule, we have highlighted in red 'troublesome' quantities, hard to evaluate. In the equation above, we have defined in red a quantity which we call "deformation potential", defined in terms of the (pseudo)potential $$V(q), of the pseudo-wavefunctions (Bloch wavefunctions), $$u_n,G(k), and of the phonon polarization vectors, all shown in red below:

Band-structure calculations help us only partially in determining the pseudopotentials, since they can only provide us with the values, $$V(G), of the pseudopotential evaluated at the vectors G of the reciprocal lattice (see the tutorials on pseudopotentials from the University of Illinois at Urbana-Champaign). Lattice-dynamics models usually follow a similar empirical strategy, by fitting adjustable parameters to the phonon spectra, but yielding eigenvectors (the polarization vectors) of dubious quality.

At low electron energies, the Bardeen-Shockley [1] deformation-potential concept makes it possible to connect the electron-phonon matrix element to the dilatation (d) and uniaxial-shear, (u) deformation potentials in terms of the angle between the wavevector of the emitted or absorbed phonon and the longitudinal axis of the ellipsoidal equi-energy surfaces near the X-valley minima:

as shown by Herring and Vogt [2] . These expressions, or their slightly simpler isotropic versions, have been used by the Modena group[3], in their pioneering work: The drift velocity of electrons in Si was measured experimentally beyond the non-Ohmic regime -- data which still represent an unsurpassed standard of accuracy -- and Monte Carlo techniques were used for the first time to investigate electron transport in Si. This work resulted in an accurate calibration of the low-energy intravalley deformation potentials -- highlighted in blue above and also discussed in the page dealing with the electron mobility in inversion layers -- as well as the intervalley deformation potentials

The Modena work still represents the definite study of electron transport in Si up to fields of a few tens of kV/cm and average energies of a few tenths of an eV. But at higher electron energies, the measurable deformation potentials cease to be related with the electron-phonon matrix elements, direct experimental data are lacking, and the seriousness of the problem emerges clearly.

Pair production rate

Considering now the pair-production process, we must consider the Coulomb interaction between the high-energy and the valence electrons. The Born approximation can be used, because of the large energies involved. Therefore, we can write for the rate at which an electron of quasi-momentum k in band n and of energy $$E_n(k) generates a pair as:

where the Coulomb matrix element, $$M_GG', results from the sum (and interference) of the 'direct' ($$M_d) and 'exchange' ($$M_x) matrix elements [6] :

where the quantities $$I_G(k n; k' n_c) and $$I_G'(p n_v; p' n''_c), etc., are the overlap integrals between the Bloch waves of the final and initial states. As usual in solid-state physics, the formulation of the problem is simple, once the potential is known. In this case, the major complication arises from the fact that the Coulomb interaction is screened by valence and conduction electrons. The dielectric screening 'matrix' (in the indices of the reciprocal-lattice vectors G and G') is indicated in red, as it constitutes our 'troublesome' quantity. In the equation above, the quasi-momentum and energy exchanged in the direct (d) and exchange (x) processes are:

Usually, the dielectric function is evaluated in the Random Phase Approximation (RPA),

which is still quite complicated: Even if the equilibrium distribution function, $$f(k) = f_eq(k), is assumed, the dependence of screening on momentum and energy transfer (the latter being known as 'dynamic screening effects') should be accounted for. Some 'reasonable' simplifications are obtained by considering only the longitudinal response, i.e., by retaining only the diagonal terms, G=G', and by considering only the 'smallest' G vectors, since the strength of the interaction will decrease as $$G^-4 for large G's. Finally, the contributions of the valence and conduction electrons are treated separately, factoring the dielectric function as follows:

where the term in square brackets represents the response of the free carriers. The response of the valence electrons can be treated statically, to a good approximation. Several 'model dielectric functions' are available in the literature in this case. The free carriers will contribute mostly at long wavelength, in which case the 'screening parameter' ß(q) becomes the inverse of the 'usual' Debye-Hückel length. So, one can handle screening, somehow, but the problem is far from being trivial and fully settled!

Assuming now that the screening problem is solved, additional common simplifications to the problem consist in approximating the matrix element, or ignoring momentum conservation. The first approximation, which works quite well in a variety of cases, is the so-called "Constant Matrix Element" (CME) approximation: One replaces the Coulomb matrix element with its average. Since this can be known only after the full-blown calculations are performed, one treats this average matrix element as a 'fitting' quantity whose magnitude is roughly given by:

A much stronger approximation works remarkably well in our case: Kane[6] has shown that for electron-initiated processes in Si so many Umklapp processes take part in the pair production process, even at small electron energies, that ignoring momentum conservation does not result in any major errors. Thus, one replaces the many 'spikes' in the momentum-conserving delta-function -- caused by the sum over G-vectors -- with smooth integrals over the Brillouin Zone, ignoring momentum conservation. The evaluation of the ionization rate is thus dramatically simplified, being reduced to a one-dimensional integral over the density of states in the valence and conduction bands, 'joined' by energy conservation:

The results obtained by using this approximation, usually called the random-k approximation, are identical to the results obtained using more sophisticated approaches, as we shall see later in Fig. 3, bottom frame. Yet, despite the validity of these approximations, numerical calculations are still tough: The full band-structure of Si must be used, the (pseudo)wavefunctions tabulated to evaluate overlap factors, the screening matrix must also be evaluated. These numerical difficulties have represented an unsurmountable obstacle for many years, Kane's[6] effort having remained isolated and somewhat ignored.

The old 'fitting' strategy

One now 'fixes' the electron-phonon scattering rates, thanks to a remarkable Ansatz coming to the rescue. If we assume, as suggested initially by Karl Hess, that the density of states (DOS) controls the energy dependence of the electron-phonon scattering rate, as explained elsewhere in these pages , then we can assume that the acoustic and optical matrix elements take the simplest possible forms:

(Note that for the s-type wavefunctions of Si near the X-valley minima, symmetry considerations impose the selection rule $$(DK)_op=0. Since one is now interested in high-energy transport, away from the valley minima, it is reasonable to disregard this selection rule). The 'fitting parameters' highlighted in blue are to be determined by fitting the drift velocity vs. field characteristics, shown in the top frame of Fig. 2. This determines transport at low energy, since the velocity-field data shown in the figure extend only to fields up to a few kV/cm, at electron energies not exceeding a few tenths of an eV. At higher energies, one can only extrapolate the low-energy deformation potentials, relying on the Ansatz above. (That this was the right thing to do has been proven only recently, as a look at the top frame of Fig. 5 below clearly shows). Finally, one fits the ionization coefficient by playing with the ionization rate: Typically, a simple, Keldysh-like expression has been very popular: Ignoring band-structure complications, one assumes that Keldysh results, obtained for the direct-gap and parabolic-band cases, can be extended to more general situations. The ionization rate then becomes simply:

where $$E_th is some threshold for the ionization process (=(3/2)$$E_gap for direct-gap, parabolic band materials) and P is an empirical dimensionless pre-factor, divided by the electron-phonon scattering (or relaxation) time at the threshold energy. While the exponent b should be close to 2 (the parabolic-band value), one can hope to obtain a better fit of the ionization coefficient letting it vary.

The main problem of this procedure is that there is no guarantee that we shall hit the 'correct' solution: At high electron energies, we have, so to speak, two unknows (the electron-phonon and the ionization rates) in one equation (the experimentally known ionization coefficient). Starting from a weak electron-phonon interaction and a 'hard' ionization process (large pre-factor in the Keldysh formula, large threshold), by increasing the electron-phonon scattering rates at high energies, we can obtain an equally satisfactory fit of the ionization coefficient using a 'softer' ionization rate (small pre-factor, small threshold). Even more serious difficulty is caused by the fact that the ionization coefficient depends sub-linearly on the strength of the ionization rate because of a sort of feedback effect: A high ionization rate causes a reduced high-energy tail in the electron energy-distribution. This, in turn, results in a marginally higher ionization coefficient compared to what results from the higher high-energy tails allowed by a lower ionization rate. Thus, even using the 'right' electron-phonon scattering rates, one could not be sure of having pinned down the 'right' ionization rate. This was the situation DAMOCLES found itself in at the early stages of its history (ca. 1988).

A final complication, which has affected historically the development of a sensible model for the scattering processes at high electron energies, is caused by the widespread use of approximated band-structure models, in the attempt of reducing the large computational requirements posed by full-band Monte Carlo models, such as DAMOCLES.

Considering the use of approximated band-structure models, of arbitrary extensions of low-energy electron-phonon matrix elements to high energies, and of oversimplified pair-production rates, the uncertainty illustrated by the range spanned by the scattering rates shown in Fig. 1, should not constitute a big surprise. Considering in particular that different kinematics (band-structure) and electron-phonon scattering models affect the high-energy tails of the distribution functions, and so the ionization coefficient, exponentially, the huge range covered in the bottom frame of Fig 1 appears expected, as 'exponentially different' ionization rates are required to compensate for the 'linearly different' electron-phonons and band-structure models.

The new 'standard model'

• Experimentally a way was devised by E. A. Cartier and F. R. McFeely[11] to measure directly the electron energy relaxation at high energy. As illustrated in Fig. 3, soft X-rays impinging on a Si wafer in vacuum excite the 2p spin-split core states up into the conduction band. If the energy of the X-rays is large enough to excite electron at energies well above the ionization threshold, analyzing the energy spectrum of the electrons which escape into the vacuum (from the cesiated Si surface) will show attenuated but unbroadened lines: Indeed, the very efficient deep-inelastic ionization process will quickly remove electrons from the energy-window of the analyzer, while phonons will not have sufficient time to broaden the lines via small-energy-loss processes. As the X-ray energy is lowered and the excited electrons approach the ionization threshold, phonon broadening will be increasingly apparent, while the line attenuation will fade with the fading ionization. Using DAMOCLES to simulate the experimental results, the ionization rate was varied, while keeping the phonon scattering rate fixed, in order to reproduce the observed attenuation and broadening with varying photon energy.

  ---

  444×666 jpeg
  
  Figure 3.Schematic diagram illustrating the XPS experimental situation: The 2p core-levels in Si are excited into the conduction band by soft X-rays (100eV of energy or more). Impact ionization is effective in attenuating the lines of the electrons emitted into vacuum, while line-broadening caused by scattering with phonons takes over at smaller X-ray energies, when the excited electrons have energies close to the ionization threshold.

  ---

  Since, even with a cesiated surface, only conduction electrons with kinetic energies exceeding the Si work-function (the energy level of the vacuum, about 1.5 eV under optimal cesiating conditions), a complementary experiment, performed much earlier by DiMaria[12], was considered: Looking at Fig. 4, electrons are now injected from the gate of a p-channel MOSFET into the channel across a thin Si$O$₂ layer. Electron-hole pairs generated by the electrons upon emerging in the Si substrate can be counted, since electrons and holes are collected separately at the source/drain and substrate contacts, respectively. By counting the pairs, and thus measuring the pair-generation-yield (called "quantum yield") as a function of gate bias (and so, of electron energy in the Si substrate), one has a direct measurement of the relative probability that an electron loses energy to phonons rather to electron-hole pairs. Simulating once more with DAMOCLES the experimental situation, the ionization rate could be determined down the energies approaching the ionization threshold.

  ---

  572×642 jpeg
  
  Figure 4.Top: Schematic diagram illustrating the "carrier separation" technique used to measure the ionization yield as a function of kinetic energy of the primary electrons. Bottom: Quantum yield obtained from the carrier-separation technique, compared to DAMOCLES simulation results. (Adapted from Appl. Phys. Lett. 62, 3339 (1993))

  ---

  
  The net result of the XPS and quantum-yield experiments, coupled to DAMOCLES simulations, is shown by the solid black line in the bottom frame of Fig. 5.

• Theoretically, the XPS/quantum-yield results triggered a flurry of activity which resulted in refined calculation, following Kane [6] , confirming the semi-experimental results. In Fig. 5, bottom frame, we show the 'ab initio' results by the group at Osaka University[13], at NTT Labs[14], at the Georgia Institute of Technology[15], as well as the rate obtained by Bude and Mastrapasqua[16], by rescaling calculated rates in order to fit substrate currents in sub-tenth-micron MOSFETs.
  
  

  ---

  419×678 jpeg
  
  Figure 5.Top: Electron-phonon scattering rates determined empirically (IBM[5] and Urbana (1991)[18]) or calculated 'ab-initio' (better: as ab-initio as possible!) by groups at IBM-Urbana[17], Osaka[19] and Urbana (1993-94)[20]. Bottom: Impact-ionization rates from the calculations by the Osaka[13], NTT[14], Georgia Tech[15], and Bell-Laboratory[16], compared to the semi-experimental 'Cartier' rate[11], obtained from XPS and quantum-yield data simulated with DAMOCLES. (Adapted from J. Appl. Phys. 78, 1058 (1995))

  ---

  Note how the orders-of-magnitude uncertainty shown in Fig. 1 has completely disappeared. Even more remarkable is the fact that the random-k approximation comes amazingly close in predicting the right result. This is 'amazing', in view of the simplicity of the procedure and one can only wonder why Kane's suggestion has been ignored for so many years. Another consequences of the validity of the random-k approximation is that the density of states appears to control not only the energy-dependence of the electron-phonon scattering rates, as anticipated by Hess, but also the ionization process. Band-structure effects, via the kinematics (group velocity) and dynamics (scattering rates via the density-of-states) are the main factor controlling electron transport in Si (and other semiconductors as well, as explained in our page dealing with the performance of sub-micron devices based on various cubic semiconductors).

  The availability of faster computers and refined numerical schemes has also allowed more detailed calculations of the electron-phonon scattering rates: Using the rigid-pseudoion, empirical-pseudopotential framework, coupled to empirical lattice-dynamics models to obtain phonon spectra and polarization vectors, an early IBM-Urbana cooperation[17] has demonstrated how well Hess' Ansatz really worked, by showing that the early 'guesswork' of DAMOCLES[5] and Urbana[18], thanks to a sensible band-structure model and Hess' Ansatz, ended up in the 'correct' energy-dependence and magnitude of the scattering rates. Work performed at Osaka University[19] has confirmed these results. But the pinnacle of these efforts is arguably represented by the work done by Doug Yoder and collaborators[20]: Using Harris potentials, they have not only evaluated the electron-phonon scattering rates 'ab-initio' (or very close to it), but they have also implemented the full dependence of the matrix element on the initial and final electronic states into a full-band Monte Carlo simulation.

Looking at Fig. 5, noting the excellent agreement among various 'ab-initio' calculations of the electron-phonon and impact-ionization scattering rates, and noting in Fig. 2 how well the resulting transport simulations reproduce the experimental data, we can say that we have almost reached the stage at which the knowledge of the Si ionic (pseudo)potential is the only infomation required to simulate quite accurately electron transport in Si up to very high (5eV or so) kinetic energies.

Acknowledgments

The 'standard model' just presented is the results of the combined efforts of several groups worldwide, in particular prof. Karl Hess' group at the University of Illinois at Urbana-Champaign, prof. Chihiro Hamaguchi's group at Osaka University, in addition to, of course, the IBM Research Division. The amazing 'consensus' reached by these various 'full-band-efforts', illustrated in Fig. 5, came somewhat as a pleasant surprise as a result of a comparison of various models of electron transport in Si, involving 47 authors worldwide [21]. This comparison was ideated and implemented by Dr. Jack Higman, coordinated by prof. Robert Dutton at Stanford University, and supported by the National Center for Computational Electronics, NCCE, co-directed by the Computational Electronics Group, at the University of Illinois at Urbana-Champaign and by Stanford University. Prof. Nobuyuki Sano and prof. K. Natori of the University of Tsukuba have actively contributed to the clarification of many problems concerning the calculation of impact-ionization rates and to the extension of the model to Ge and GaAs.

Abstract and post-script version of an article on which this page is based are available from the IBM Research Division CyberDigest.

References

1. W. Shockley and J. Bardeen, Phys. Rev. 77, 497 (1950); J. Bardeen and W. Shockley, Phys. Rev. 80, 72 (1950).
2. C. Herring and E. Vogt, Phys. Rev. 101, 944 (1956).
3. C. Canali, C. Jacoboni, F. Nava, G. Ottaviani, and A. Alberigi-Quaranta, Phys. Rev. B 12, 2265 (1975). [Modena 1975]
4. J. Y. Tang and K. Hess, J. Appl. Phys. 54, 5139 (1983). [Urbana 1983]
5. M. V. Fischetti and S. E. Laux, Phys. Rev. B 38, 9721 (1988). [IBM 1988]
6. E. O. Kane, Phys. Rev. 159, 624 (1967).
7. Th. Vogelsang and W. Hänsch, J. Appl. Phys. 70, 1493 (1991). [Siemens 1991]
8. R. Thoma, H. J. Pfeifer, W. L. Engl, W. Quade, R. Brunetti, and C. Jacoboni, J. Appl. Phys. 69, 2300 (1991). [Aachen 1991]
9. J. Bude, K. Hess, and G. J. Iafrate, Phys. Rev. B 45, 10958 (1992). [Urbana 1992]
10. N. Sano and A. Yoshii, Phys. Rev. B 45, 4171 (1992). [NTT 1992]
11. E. A. Cartier, M. V. Fischetti, E. A. Eklund, and F. R. McFeely, Appl. Phys. Lett. 62, 3339 (1993). [IBM 1993]
12. D. J. DiMaria, T. N. Theis, J. R. Kirtley, F. L. Pesavento, and D. W. Dong, J. Appl. Phys. 57, 1214 (1985).
13. Y. Kamakura, H. Mizuno, M. Yamaji, M. Morifuji, K. Taniguchi, C. Hamaguchi, T. Kunikiyo, and M. Takenaka, J. Appl. Phys. 75, 3500 (1994). [Osaka 1994]
14. N. Sano and A. Yoshii, J. Appl. Phys. 75, 5102 (1994). [NTT 1994]
15. J. Kolnik, Y. Wang, I. H. Oguzman, and K. F. Brennan, J. Appl. Phys. 76, 3542 (1994). [Georgia Tech 1994]
16. J. Bude and M. Mastrapasqua, IEEE Electron Device Lett. 16, 439 (1995). [Bell Labs 1995]
17. M. V. Fischetti and J. M. Higman, in "Monte Carlo Device Simulation: Full Band and Beyond", K. Hess Ed., (Kluwer Academic Press, Boston, Massachusetts, 1991), p. 123. [IBM-Urbana 1991]
18. H. Shichijo, J. Y. Tang, J. Bude, and P. D. Yoder, in "Monte Carlo Device Simulation: Full Band and Beyond", K. Hess Ed., (Kluwer Academic Press, Boston, Massachusetts, 1991), p. 285. [Urbana 1991]
19. T. Kunikiyo, M. Takenaka, Y. Kamakura, M. Yamaji, H. Mizuno, M. Morifuji, K. Taniguchi, and C. Hamaguchi, J. Appl. Phys.75 297 (1994). [Osaka 1994]
20. P. D. Yoder, J. M. Higman, J. Bude, and K. Hess, Semicond. Sci. Technol. 7, B357 (1992); P. D. Yoder, V. D. Natoli, and R. M. Martin, J. Appl. Phys. 73, 4378 (1993); P. D. Yoder and K. Hess, Semicond. Sci. Technol. 9, 852 (1994). [Urbana 1993-1994]
21. A. Abramo, L. Baudry, R. Brunetti, R. Castagne, M. Charef, F. Dessenne, P. Dollfus, R. Dutton, W. Engl, R. Fauquembergue, C. Fiegna, M. Fischetti, S. Galdin, N. Goldsman, M. Hackel, C. Hamaguchi, K. Hess, K. Hennacy, P. Hesto, J. Higman, T. Iizuka, C. Jungemann, Y. Kamakura, H. Kosina, T. Kunikiyo, S. Laux, H. Lin, C. Maziar, H. Mizuno, H. Peifer, S. Ramaswamy, N. Sano, P. Scrobohaci, S. Selberherr, M. Takenaka, T. Tang, K. Taniguchi, J. Thobel, R. Thoma, K. Tomizawa, M. Tomizawa, T. Vogelsang, S. Wang, X. Wang, C. Yao, P. Yoder and A. Yoshii, IEEE Transactions of Electron Devices 41, 1646 (1994).

Return to the top of the page

Details about numerical algorithms and parameters used in DAMOCLES to evaluate scattering rates

APPENDIX: List of mathematical symbols

	damoclesNO-SPAM@watson.ibm.com (last updated: January 26, 1999)