*2.1. Semiclassical Boltzmann Transport*

To calculate macroscopic transport coefficients, such as electronic or thermal conductivity, one needs to analyze the microscopic processes happening when electrons are scattered by phonons or impurities in a metal or semiconductor. The idea is to treat the electronic excitations as particles and

follow their motion over time. For example, an electron in a metal exposed to an electric field **E** gets accelerated to a velocity **v**, and over time gains the kinetic energy

$$
\Delta \epsilon = \mathbf{v} \mathbf{E} t.\tag{1}
$$

However, the velocity of the electron cannot grow forever as the electron will be scattered sooner or later with impurities or by phonons, changing its direction and losing part of its kinetic energy. If we do not know the exact microscopic details of how this scattering happens, we might just consider a macroscopic relaxation time *τ*. This is the average time an electron gets accelerated by the electric field before being scattered, i.e., the time between electron–phonon transitions. After the collision takes place the particle starts again with its acceleration parallel to the electric field until further scattering. The averaged velocity (over time and many different scattering events) due to acceleration and scattering is responsible for the finite electrical conductivity, as the current density is **J** = *ne***v**, where *n* is the electron concentration and *e* its charge.

Let us now consider the conductivity of a metal or semiconductor in more details. An electron in a crystal occupies states according to its distribution function *f*(**<sup>r</sup>**, **k**, *t*), giving the probability of finding an electron at position **r** and time *t* with the crystal momentum **k**. At equilibrium, when no fields or temperature gradients are applied, this is given by Fermi function,

$$f\_0(\mathbf{k}, \mu, T) = \frac{1}{e^{(\varepsilon\_n(\mathbf{k}) - \mu)/k\_B T} + 1},\tag{2}$$

where *μ* is the total chemical potential and  *n*(**k**) the energy of an electron in the n-th band with momentum **k**, *kB* the Boltzmann's constant and *T* the temperature. However, if we now apply, for example, a potential bias, electrons will ge<sup>t</sup> excited and the actual distribution functions shifts to higher energies. The equal interplay between acceleration on the one hand and collisions and scattering on the other, will possibly set up a steady-state condition. When we switch off the electrical field again, the system relaxes back to its equilibrium state.

In 1872, Boltzmann laid down an equation for *f* connecting thermodynamics with non-equilibrium kinetics [8]. Although this was long before the birth of modern quantum mechanics, his transport theory consists of a probabilistic description for the one-particle distribution function *f*(**<sup>r</sup>**, **k**, *t*). This equation for out-of-equilibrium situations is called the semiclassical Boltzmann transport equation (BTE). Here, we will lay down the basic equation and see how to solve it for systems which are close to equilibrium, while in Section 2.2 we review a fully quantum-mechanical treatment for transport properties. As mentioned before, we will treat crystal electrons as semiclassical particles fulfilling the equations

$$\frac{d\mathbf{k}}{dt}\_{\prime} = \frac{1}{\hbar} \mathbf{F}\_{\prime \prime} \tag{3}$$

$$\mathbf{v}\_{\parallel} = \frac{1}{\hbar} \nabla\_{\mathbf{k}} \epsilon\_{\parallel}(\mathbf{k}).\tag{4}$$

Here, **F***e* is the force acting on the electron and **v** its velocity. Considering the total change in time of the distribution function

$$\frac{\mathrm{d}f(\mathbf{r},\mathbf{k},t)}{\mathrm{d}t} = \frac{\partial f}{\partial t} + \frac{\partial \mathbf{r}}{\partial t} \cdot \frac{\partial f}{\partial \mathbf{r}} + \frac{\partial \mathbf{k}}{\partial t} \frac{\partial f}{\partial \mathbf{k}} = \frac{\partial f}{\partial t} + \mathbf{v} \cdot \nabla\_{\mathbf{r}} f + \frac{1}{\hbar} \mathbf{F}\_{\mathbf{r}} \cdot \nabla\_{\mathbf{k}} f,\tag{5}$$

and setting it equal to the scattering rate, we arrive at the BTE,

$$\frac{\partial f}{\partial t} + \mathbf{v} \cdot \nabla\_{\mathbf{r}} f + \frac{1}{\hbar} \mathbf{F}\_{\mathbf{r}} \cdot \nabla\_{\mathbf{k}} f = \left. \frac{\partial f}{\partial t} \right|\_{\text{scatt}}.\tag{6}$$

The BTE states that the electron distribution changes due to in- and out-going scattering events which together form *∂ f ∂t* |scatt. Although the Boltzmann transport theory is probabilistic, it is still classical as in quantum mechanics one cannot specify the canonical variables **p** and **r** simultaneously, due to the Heisenberg uncertainty principle, Δ*p*Δ*x* ≥ *h*¯. The assumption here is that the uncertainty in space and momentum is small enough compared to the system size that the electrons can be treated as particles. The scattering term, *∂ f ∂t* |scatt, is the most interesting and complicated part of the BTE which accounts for the change of the electron probability due to electron–electron or electron–phonon scattering events, which we can treat quantum mechanically.

Without taking care of the exact form of the underlying scattering mechanisms, the constant relaxation time approximation (CRTA)assumes that the system seeks to return to the equilibrium configuration *f*0 in a time *τ* after being disturbed by external electric or temperature fields. In general, *τ*, known as relaxation time, should depend on the direction of the scattering events, energy, and on the exact scattering mechanisms. However, the easiest version of the CRTA assumes one constant momentum relaxation time for all modes, direction, and scattering processes and the collision term can be written in the form [9,10]

$$\frac{\partial f(\mathbf{r}, \mathbf{k}, t)}{\partial t}\Big|\_{\text{scatt}} \simeq -\frac{\delta f(\mathbf{r}, \mathbf{k}, t)}{\tau} \tag{7}$$

Here, the rate of change of *f*(**<sup>r</sup>**, **k**, *t*) due to collisions is assumed to be proportional to the deviation from the equilibrium distribution, *δ f* = *f* − *f*0. To better understand the meaning of the CTRA, we assume no electric field present (no **F***e*) and spatial uniformity (gradient with respect to **r** vanishes). Then the BTE, Equation (6), in the CRTA takes the form

$$\frac{\partial f}{\partial t} = -\frac{\delta f}{\pi},\tag{8}$$

which has the familiar solution for the deviation from equilibrium

$$
\delta f(t) = \delta f(0)e^{-t/\tau}.\tag{9}
$$

This means that the perturbed system relaxes with a typical timescale *τ* to its equilibrium distribution. The CRTA is a crude approximation that fails for certain systems as we will see later in this section.

When going beyond the simple CRTA and concentrating on the phonon-limited carrier mobilities in semiconductors, it is common to consider the electron–phonon scattering as the dominant mechanism. Furthermore, we adopt the notation for crystals in k-space, as such, the distribution function depends on the band index *n* and the momentum **k**, i.e., *fn***<sup>k</sup>**. Therefore, for time-independent and homogenous fields, the BTE in k-space reads

$$\frac{\partial \mathbf{k}}{\partial t} \frac{\partial f\_{n\mathbf{k}}}{\partial \mathbf{k}} = \left. \frac{\partial f\_{n\mathbf{k}}}{\partial t} \right|\_{\text{scatt}}.\tag{10}$$

Neglecting magnetic fields and considering only the presence of electric fields, the BTE for phonon-limited carrier mobilities is given by [11,12]

$$\begin{array}{rcl} \epsilon \mathbf{E} \frac{\partial f\_{\mathsf{k}\mathsf{k}}}{\partial \mathbf{k}} &=& \frac{\Omega}{(2\pi)^{3\mathsf{k}}} \sum\_{\mathsf{n},\mathsf{p}} \int \mathbf{d} \mathbf{q} |g\_{\mathsf{m}\mathsf{n}\mathsf{p}}(\mathbf{k}, \mathbf{q})|^{2} \\ & \times \left\{ f\_{\mathsf{n}\mathsf{k}} \left[1 - f\_{\mathsf{m}\mathsf{k}+\mathsf{q}} \right] \left[ (n\_{\mathsf{q}\mathsf{p}} + 1) \delta(\varepsilon\_{\mathsf{n}\mathsf{k}} - \varepsilon\_{\mathsf{n}\mathsf{k}+\mathsf{q}} - \hbar \omega\_{\mathsf{q}\mathsf{p}}) + n\_{\mathsf{q}\mathsf{p}} \delta(\varepsilon\_{\mathsf{n}\mathsf{k}} - \varepsilon\_{\mathsf{n}\mathsf{k}+\mathsf{q}} + \hbar \omega\_{\mathsf{q}\mathsf{p}}) \right] \\ & - f\_{\mathsf{m}\mathsf{k}+\mathsf{q}} \left[1 - f\_{\mathsf{n}\mathsf{k}} \right] \left[ (n\_{\mathsf{q}\mathsf{p}} + 1) \delta(\varepsilon\_{\mathsf{n}\mathsf{k}} - \varepsilon\_{\mathsf{n}\mathsf{k}+\mathsf{q}} + \hbar \omega\_{\mathsf{q}\mathsf{p}}) + n\_{\mathsf{q}\mathsf{p}} \delta(\varepsilon\_{\mathsf{n}\mathsf{k}} - \varepsilon\_{\mathsf{n}\mathsf{k}+\mathsf{q}} - \hbar \omega\_{\mathsf{q}\mathsf{p}}) \right] \end{array} \tag{11}$$

The first two terms in the second line represent scattering of an electron from energy band *n* with momentum **k** into the state *m*, **k** + **q** by either the emission (*<sup>n</sup>***q***v* + 1) or the absorption (*<sup>n</sup>***q***v*) of a phonon with frequency *<sup>ω</sup>***q***<sup>v</sup>*. The factor *fn***k**[<sup>1</sup> − *fm***<sup>k</sup>**+**q**] make sure that the initial electronic state is occupied while the final state is unoccupied. The last two terms represent the backscattering events from state *m*, **k** + **q** towards the state in the energy band *n* with momentum **k**. Here, Ω is the volume of the primitive unit cell, the *δ*-functions ensure energy conservation during the scattering process and *<sup>n</sup>***q***v* is the Bose-Einstein distribution function, giving the probability that a crystal phonon with moment **q** is present in branch *v*. The term *gmnv*(**<sup>k</sup>**, **q**) is the electron–phonon coupling matrix element, which is normally calculated within density functional perturbation theory [13–15]. Taking the derivative with respect to the *i*-componen<sup>t</sup> of the electric field *Ei*, and linearizing the equation around the equilibrium distribution *f* 0, we ge<sup>t</sup> a direct equation for *<sup>∂</sup>Eifn***k** [12]

$$\frac{\partial f\_{n\mathbf{k}}}{\partial E\_i} = \varepsilon \frac{\partial f\_{n\mathbf{k}}^0}{\partial \varepsilon\_{n\mathbf{k}}} v\_{n\mathbf{k},i} \tau\_{n\mathbf{k}} + \frac{\Omega \tau\_{n\mathbf{k}}}{(2\pi)^2 \hbar} \sum\_{m,v} \int d\mathbf{q} \frac{\partial f\_{m\mathbf{k}+\mathbf{q}}}{\partial E\_i} |g\_{mvv}(\mathbf{k}, \mathbf{q})|^2 \tag{12}$$
 
$$\times \left\{ \left[1 + n\_{q\mathbf{v}} - f\_{n\mathbf{k}}^0\right] \delta(\varepsilon\_{n\mathbf{k}} - \varepsilon\_{m\mathbf{k}+\mathbf{q}} + \hbar \omega\_{q\mathbf{k}}) + \left[n\_{q\mathbf{v}} + f\_{n\mathbf{k}}^0\right] \delta(\varepsilon\_{n\mathbf{k}} - \varepsilon\_{m\mathbf{k}+\mathbf{q}} - \hbar \omega\_{q\mathbf{k}}) \right\},$$

where the inverse of the electron energy relaxation time due to the electron–phonon interaction is

$$\tau\_{n\mathbf{k}}^{-1}(T,\mu) = \frac{\Omega}{(2\pi)^{2}\hbar} \sum\_{m,\nu} \int d\mathbf{q} |g\_{m\mathbf{n}\mathbf{r}}(\mathbf{k},\mathbf{q})|^{2} \left\{ \left[ n\_{\mathbf{r}\mathbf{q}} + f\_{m\mathbf{k}+\mathbf{q}}^{0} \right] \delta(\varepsilon\_{n\mathbf{k}} - \varepsilon\_{m\mathbf{k}+\mathbf{q}} + \hbar\omega\_{\mathbf{q}\nu}) \right. $$

$$+ \left[ n\_{\mathbf{r}\mathbf{q}} + 1 - f\_{m\mathbf{k}+\mathbf{q}}^{0} \right] \delta(\varepsilon\_{n\mathbf{k}} - \varepsilon\_{m\mathbf{k}+\mathbf{q}} - \hbar\omega\_{\mathbf{q}\nu}) \right\}. \tag{13}$$

Note here that the Fermi distribution for the electrons is a function of both temperature and chemical potential, while the phonon distribution function depends on the temperature. Equation (12) is the linearized BTE (LBTE) for phonon-limited carrier transport and is valid for most semiconductors where the acceleration of the free carriers is smaller than the thermal energy, *e***Ev***τ kBT*. We would like to point out that in contrast to the CRTA, here the relaxation time depends on the energy and momentum of the electron. The LBTE needs to be solved self-consistently for the variation of the electron distribution function with respect to the applied field, therefore also called iterative BTE. By neglecting the integral part of the LBTE, one obtains an equation for the distribution function which can be solved directly,

$$\frac{\partial f\_{n\mathbf{k}}}{\partial E\_i} = \varepsilon \frac{\partial f\_{n\mathbf{k}}^0}{\partial \varepsilon\_{n\mathbf{k}}} v\_{n\mathbf{k},i} \tau\_{n\mathbf{k}} \tag{14}$$

and is called the self-energy relaxation time approximation (SERTA). This is due to the analogy that the relaxation time *τ* in Equation (13) is related to the Fan-Migdal electron self-energy by *τ* −1 *n***k** = 2ImΣFM*n***<sup>k</sup>** [13]. In the SERTA, the electron mobility takes the form

$$\mu\_{a\S}(T,\mu) = -\frac{e\Omega}{n\_{\varepsilon}(2\pi)^{3}}\sum\_{n\in\text{CB}}\int d\mathbf{k}\frac{\partial f\_{n\mathbf{k}}^{0}(T,\mu)}{\partial \varepsilon\_{n\mathbf{k}}}v\_{n\mathbf{k},\mu}v\_{n\mathbf{k},\emptyset}\tau\_{n\mathbf{k}}(T,\mu). \tag{15}$$

When calculating the mobility within SERTA or LBTE, the relaxation time in Equation (13) is evaluated on a very fine **k**-grid by integrating over a dense **q**-grid for the phonons. This direct Brillouin zone sampling is computationally very demanding and hence cannot be applied for complex material screening. For example, the EPW code [16] calculates the electron–phonon coupling on a coarse grid in the BZ and maps it onto a fine grid by using Wannier's function interpolation. This method is however still very costly and therefore other approaches are needed to tackle the calculations of thermoelectric transport properties. One approach solves the BTE within the CRTA, Equation (7): This leads to good results for electrical conductors where the energy relaxation time depends weakly on the electron energy, , [10,17] and allows for a single constant relaxation time. However, when performing the CRTA one needs to estimate the relaxation time within simplified models, such as the deformation potential (DP) approximation [18–20] or Allen's formalism [21,22]. For most materials, *τ* is strongly anisotropic, depending on energy and carrier concentration and the CRTA cannot be applied and we need to revert to first-principles computations for predicting *τn***k**. The so-called electron–phonon averaged (EPA) approximation [23] turns the demanding integral over momentum in Equation (13) into an integration over energy. This is done by replacing quantities that depend on momentum (|*gmnv*(**<sup>k</sup>**, **q**)|2, *<sup>ω</sup>v***q**) by their energy-dependent averages (*g*<sup>2</sup>*v*( *<sup>n</sup>***k**,  *m***k**+**q**) , *<sup>ω</sup>v*). This allows a much coarser grid in the electron energies therefore reducing computational cost drastically. Technical details and derivations can be found in [23] and the final relaxation time within the EPA is given by

$$\pi^{-1}(\varepsilon,\mu,T) = \frac{2\pi\Omega}{2\hbar} \sum\_{\upsilon} \left\{ g\_{\upsilon}^{2}(\varepsilon,\varepsilon + \overline{\omega\_{\upsilon}}) \left[ n(\overline{\omega}\_{\upsilon},T) + f^{0}(\varepsilon + \overline{\omega}\_{\upsilon}) \right] \rho(\varepsilon + \overline{\omega}\_{\upsilon}) \right. \tag{16}$$

$$+ g\_{\upsilon}^{2}(\varepsilon,\varepsilon - \overline{\omega\_{\upsilon}}) \left[ n(\overline{\omega}\_{\upsilon},T) + 1 - f^{0}(\varepsilon - \overline{\omega}\_{\upsilon}) \right] \rho(\varepsilon - \overline{\omega}\_{\upsilon}) \right\},\tag{16}$$

where *ωv* is the average phonon mode energy and *ρ*() is the electron density of states per unit energy and unit volume. The EPA is implemented in the Quantum Espresso suite [24] and Boltztrap code [25].

Figure 1 compares the relaxation time calculated within Equation (13) by the EPW code and within the EPA, Equation (16) for (a) HfCoSb and (b) HfNiSn. We see that the relaxation times in both approaches are in a good agreement, and furthermore that approximating the relaxation times with a constant, as done in CRTA, might fail for predicting the electron conductivity of real materials.

**Figure 1.** Relaxation time for (**a**) HfCoSb and (**b**) HfNiSn calculated within the EPA (Equation (16)) compared to the complete sampling of the Brillion zone (Equation (13)). Both approaches are in good agreemen<sup>t</sup> and the CRTA (dashed blue line) might fail to accurately calculate the conductivity. Reprinted with permission from [23]. Copyright (2018) John Wiley and Sons.
