3  Energy bands

This chapter explores how periodicity leads to the emergence of energy bands and band gaps for states living in a crystal.

3.1 Band fundamentals¹

The key result underlying all band theory is Bloch’s theorem: the eigenstates of a single electron in a periodic potential can always be written as Bloch states. First things first, what is a Bloch state?

Definition 3.1: Bloch state

A Bloch state is a simultaneous eigenstate of all lattice translation operators and satisfies

\[\psi_{\mathbf{k}}(\mathbf{r}+\mathbf{R})=e^{i\mathbf{k}\cdot\mathbf{R}}\,\psi_{\mathbf{k}}(\mathbf{r}) \tag{3.1}\]

for every \(\mathbf{R}\in\mathcal{BL}\). The wavevector \(\mathbf{k}\) is called quasi-momentum and is customarily used as a subscript, \(\psi_{\mathbf{k}}\), to label the state. Equivalently, every Bloch state can be factored as a phase factor times a periodic atomic factor

\[\psi_{\mathbf{k}}(\mathbf{r})=e^{i\mathbf{k}\cdot\mathbf{r}}\,u_{\mathbf{k}}(\mathbf{r}) \tag{3.2}\]

where \(u_{\mathbf{k}}(\mathbf{r}+\mathbf{R})=u_{\mathbf{k}}(\mathbf{r})\) for any \(\mathbf{R}\in\mathcal{BL}\). The equivalence of the two formulations is left to the reader.

CautionNote: \(\psi_\mathbf{k}(\mathbf{r})\) and \(\psi_{\mathbf{k}+\mathbf{g}}(\mathbf{r})\) with \(\mathbf{g}\in\mathcal{RL}\) describe just the same state!

Note here that the factoring \(e^{i\mathbf{k}\cdot\mathbf{r}}\,u_{\mathbf{k}}(\mathbf{r})\) is somewhat arbitrary, it can be done in many equivalent ways. For any \(\mathbf{g}\in\mathcal{RL}\) one has \(e^{i\mathbf{g}\cdot\mathbf{R}}=1\) for all \(\mathbf{R}\in\mathcal{BL}\), thus the translation eigenvalue is identical for \(\mathbf{k}\) and \(\mathbf{k}+\mathbf{g}\):

\[e^{i(\mathbf{k}+\mathbf{g})\cdot\mathbf{R}}=e^{i\mathbf{k}\cdot\mathbf{R}}\cdot e^{i\mathbf{g}\cdot\mathbf{R}}=e^{i\mathbf{k}\cdot\mathbf{R}}.\]

Does it mean the state is the same? Yes, Equation 3.2 admits alternative decompositions

\[\psi_{\mathbf{k}}(\mathbf{r}) =e^{i\mathbf{k}\cdot\mathbf{r}}\,u_{\mathbf{k}}(\mathbf{r}) =e^{i(\mathbf{k}+\mathbf{g})\cdot\mathbf{r}}\,\underbrace{e^{-i\mathbf{g}\cdot\mathbf{r}}\,u_{\mathbf{k}}(\mathbf{r})}_{\text{still lattice-periodic}},\]

where the new function \(e^{-i\mathbf{g}\cdot\mathbf{r}}u_\mathbf{k}(\mathbf{r})\), being a product of two periodic functions, is still periodic on the Bravais lattice. In other words \(\mathbf{k}\sim\mathbf{k}+\mathbf{g}\) are equivalent and label the same wavefunction and physical state, only the partition between the phase factor and the atomic factor is different. Among all the equivalent \(\mathbf{k}\) we choose the smallest one, i.e. the one in the WS cell of the \(\Gamma\) point, aka first Brillouin zone (see Section 2.3).

Periodicity of \(\mathcal{H}\) implies the problem is a bit simpler than it seems…

Theorem 3.1: Bloch’s theorem

Let \(\mathcal{H}\) be a single-particle Hamiltonian with \(V\) periodic on the Bravais lattice,

\[\mathcal{H} = \frac{\mathbf{p}^2}{2m} + V(\mathbf{r})\]

where \(V(\mathbf{r}+\mathbf{R}) = V(\mathbf{r})\) for any \(\mathbf{R}\in\mathcal{BL}\). Then \(\mathcal{H}\) can be diagonalised in a basis of Bloch states,

\[\mathcal{H}\,\psi_{n\mathbf{k}}(\mathbf{r}) = E_n(\mathbf{k})\,\psi_{n\mathbf{k}}(\mathbf{r}), \tag{3.3}\]

where \(n\) is the band index and \(\mathbf{k}\) the quasi-momentum. In addition, when two Bloch states have inequivalent \(\mathbf{k}_1\) and \(\mathbf{k}_2\), then we always have \(\langle\psi_{\mathbf{k}_1}|\mathcal{H}|\psi_{\mathbf{k}_2}\rangle=0\), i.e. zero coupling even when \(\psi_{\mathbf{k}_1}\) and \(\psi_{\mathbf{k}_2}\) are not eigenstates.

Show proof

First proof — translation symmetry. Let \(\mathcal{T}_{\mathbf{\Delta r}}\) be the translation operator by a vector \(\mathbf{\Delta r}\). By hypothesis the Hamiltonian is translationally invariant, \([\mathcal{H},\mathcal{T}_{\mathbf{R}}]=0\) for all \(\mathbf{R}\in\mathcal{BL}\), so \(\mathcal{H}\) and all the \(\mathcal{T}_{\mathbf{R}}\) admit a common eigenbasis. Let \(\psi\) be a simultaneous eigenstate of such a basis:

\[\mathcal{H}\psi = E\psi, \qquad \mathcal{T}_{\mathbf{R}}\psi = \lambda(\mathbf{R})\psi.\]

The group property \(\mathcal{T}_{\mathbf{R}_1}\mathcal{T}_{\mathbf{R}_2}=\mathcal{T}_{\mathbf{R}_1+\mathbf{R}_2}\) implies \(\lambda(\mathbf{R}_1+\mathbf{R}_2)=\lambda(\mathbf{R}_1)\,\lambda(\mathbf{R}_2)\). It is therefore enough to know the phases induced by the primitive translations, \(\lambda(\mathbf{t}_j)=e^{i\varphi_j}\), to obtain it for any \(\mathbf{R}=n_1\mathbf{t}_1+n_2\mathbf{t}_2+n_3\mathbf{t}_3\)…

\[ \begin{aligned} \lambda(\mathbf{R}) &= \lambda\!\left(\sum_{j=1}^3 n_j\,\mathbf{t}_j\right) = \prod_{j=1}^3 \lambda(\mathbf{t}_j)^{n_j} = \prod_{j=1}^3 \exp(in_j\varphi_j)\\ &= \exp\!\left(i\sum_j n_j\,\varphi_j\right) = \exp\!\left(i\sum_j \frac{\mathbf{g}_j\cdot\mathbf{R}}{2\pi}\,\varphi_j\right) = e^{i\mathbf{k}\cdot\mathbf{R}}, \end{aligned} \]

with \(\mathbf{k}=\sum_j(\varphi_j/2\pi)\,\mathbf{g}_j\). Thus the state satisfies Equation 3.1. \(\blacksquare\)

Second proof — plane-wave basis. Consider the plane-wave basis \(w_{\mathbf{k}}(\mathbf{r})=e^{i\mathbf{k}\cdot\mathbf{r}}/\sqrt{V}\), with \(V\) the crystal volume. Expanding the periodic potential on the reciprocal lattice, \(V(\mathbf{r})=\sum_{\mathbf{q}\in\mathcal{RL}}V_{\mathbf{q}}\,e^{i\mathbf{q}\cdot\mathbf{r}}\), the Hamiltonian acts on a plane wave as

\[\mathcal{H}\,w_{\mathbf{k}}(\mathbf{r}) = \varepsilon_{\mathbf{k}}\,w_{\mathbf{k}}(\mathbf{r}) + \sum_{\mathbf{q}\in\mathcal{RL}} V_{\mathbf{q}}\,w_{\mathbf{k}+\mathbf{q}}(\mathbf{r}),\]

where \(\varepsilon_{\mathbf{k}}=\hbar^2\mathbf{k}^2/2m\). The result is a superposition of plane waves whose wavevectors all differ from \(\mathbf{k}\) by reciprocal lattice vectors. It therefore belongs to a subspace \(S_{\mathbf{k}}\) generated by the superposition of plane waves \(w_{\mathbf{k}+\mathbf{q}}(\mathbf{r})\) with \(\mathbf{q}\in\mathcal{RL}\):

\[\mathcal{H}w_\mathbf{k}(\mathbf{r}) \in S_{\mathbf{k}} \equiv \mathrm{span}\left\{\,w_{\mathbf{k}+\mathbf{q}}(\mathbf{r})\,\bigl|\,\mathbf{q}\in\mathcal{RL}\,\right\},\]

and more generally \(\mathcal{H}\,S_{\mathbf{k}}\subseteq S_{\mathbf{k}}\). The \(S_{\mathbf{k}}\) provide a decomposition of the Hilbert space:

• \(S_{\mathbf{k}_1}=S_{\mathbf{k}_2}\iff\mathbf{k}_1-\mathbf{k}_2\in\mathcal{RL}\): two subspaces are therefore either coincident or orthogonal.
• A generic plane wave \(w_\mathbf{k}(\mathbf{r})\) belongs to \(S_\mathbf{k}\) and the Hilbert space thus decomposes as \(\bigoplus_{\mathbf{k}\in\mathrm{BZ}} S_{\mathbf{k}}.\)

Therefore \(\mathcal{H}\) is already block diagonal and eigenstates can be searched inside each single \(S_{\mathbf{k}}\). Any element of \(S_{\mathbf{k}}\) is a Bloch state with quasi-momentum \(\mathbf{k}\), in fact

\[\psi(\mathbf{r}) = \sum_{\mathbf{q}\in\mathcal{RL}} c_{\mathbf{q}}\,e^{i(\mathbf{k}+\mathbf{q})\cdot\mathbf{r}} = e^{i\mathbf{k}\cdot\mathbf{r}}\,\underbrace{\sum_{\mathbf{q}\in\mathcal{RL}} c_{\mathbf{q}}\,e^{i\mathbf{q}\cdot\mathbf{r}}}_{u_{\mathbf{k}}(\mathbf{r})},\]

and since \(\mathbf{q}\in\mathcal{RL}\), the function \(u_{\mathbf{k}}(\mathbf{r})\) is lattice-periodic. Given \(\mathcal{H}\) is basically block-diagonal over the \(S_\mathbf{k}\) spaces, it does not mix any pair of Bloch states with inequivalent quasi-momenta. \(\blacksquare\)

CautionParametric Hamiltonian formulation: \(\mathcal{H}(\mathbf{k})\)

Note that, given \(\left[\mathbf{p},e^{i\mathbf{k}\cdot\mathbf{r}}\right] = \hbar\mathbf{k}\,e^{i\mathbf{k}\cdot\mathbf{r}}\) and the consequent

\[\mathbf{p}\,\psi_\mathbf{k}(\mathbf{r})=e^{i\mathbf{k}\cdot\mathbf{r}}(\mathbf{p}+\hbar\mathbf{k})u_\mathbf{k}(\mathbf{r}),\]

any operator equation for \(\psi_{\mathbf{k}}(\mathbf{r})\) can be easily rewritten as one for \(u_{\mathbf{k}}(\mathbf{r})\) by substituting \(\mathbf{p}\) with \(\mathbf{p}+\hbar\mathbf{k}\). In particular, the eigenvalue problem of Equation 3.3 for \(\psi_{n\mathbf{k}}(\mathbf{r})\) can be easily rewritten as

\[\mathcal{H}(\mathbf{k})\,u_{n\mathbf{k}}(\mathbf{r}) = E_n(\mathbf{k})\,u_{n\mathbf{k}}(\mathbf{r}), \tag{3.4}\]

for the periodic functions \(u_{n\mathbf{k}}(\mathbf{r})\), where we introduced a parametric Hamiltonian

\[\mathcal{H}(\mathbf{k}) = \frac{(\mathbf{p}+\hbar\mathbf{k})^{2}}{2m} + V(\mathbf{r}). \tag{3.5}\]

Now \(\mathbf{k}\) is just a parameter of the problem and a discrete energy spectrum has to be found for eigenfunctions \(u_\mathbf{k}(\mathbf{r})\) defined on a single cell with periodic boundaries.

3.1.1 Band properties: Time-reversal symmetry

A fundamental property of the band structure follows from time-reversal symmetry: in the absence of magnetic fields, the energy bands are symmetric under inversion of the wavevector.

Theorem 3.2: Kramers degeneracy

For a Hamiltonian with a real potential \(V(\mathbf{r})\) and no spin–orbit coupling, the energy bands satisfy

\[E_n(\mathbf{k}) = E_n(-\mathbf{k}). \tag{3.6}\]

That is, the band structure is symmetric about the origin in \(\mathbf{k}\)-space.

Show proof (spinless)

In short, the time-reversal operator on a spinless state acts on eigenstates as a simple conjugation

\[\psi(\mathbf{r}) \to \psi^*(\mathbf{r}),\]

thus the phase factor goes from \(e^{+i\mathbf{k}\cdot\mathbf{r}}\) to \(e^{-i\mathbf{k}\cdot\mathbf{r}}\), and the time-reversed eigenstate is still a Bloch state with a reversed quasi-momentum and identical energy. The inclusion of spin in the case of an \(\mathcal{H}\) with spin-orbit coupling requires taking into account how spin states transform under time reversal, but this leads to similar conclusions: basically, \(E_{n\uparrow}(\mathbf{k})=E_{n\downarrow}(-\mathbf{k})\). \(\blacksquare\)

This symmetry implies that the band structure is mirror symmetric about \(\mathbf{k} = 0\) and, by the periodicity of the reciprocal lattice, also about the Brillouin zone boundaries.

3.1.2 Band properties: Smoothness

The energy bands \(E_n(\mathbf{k})\) are smooth functions of \(\mathbf{k}\) away from degeneracy points. This can be derived by applying Theorem 3.3 to the parametric family \(\mathcal{H}(\mathbf{k})\) of Equation 3.5,

\[\nabla_{\mathbf{k}}E_n(\mathbf{k}) = \langle u_{n\mathbf{k}}|\,\nabla_{\mathbf{k}}\mathcal{H}(\mathbf{k})\,|u_{n\mathbf{k}}\rangle = \frac{\hbar}{m}\,\langle u_{n\mathbf{k}}|\,\mathbf{p}+\hbar\mathbf{k}\,|u_{n\mathbf{k}}\rangle = \hbar\,\langle\psi_{n\mathbf{k}}|\,\mathbf{v}\,|\psi_{n\mathbf{k}}\rangle,\]

where \(\mathbf{v}=\mathbf{p}/m\) is the velocity operator. This implies that the band (group) velocity equals

\[\mathbf{v}_n(\mathbf{k}) = \frac{1}{\hbar}\nabla_{\mathbf{k}}E_n(\mathbf{k}), \tag{3.7}\]

a relation that will play a central role in semiclassical transport (Section 8.1).

Theorem 3.3: Hellmann-Feynman

Let \(\mathcal{H}(\lambda)\) be a Hamiltonian depending smoothly on a real parameter \(\lambda\), with a non-degenerate normalised eigenstate \(|\psi(\lambda)\rangle\) and eigenvalue \(E(\lambda)\). Then

\[\frac{dE}{d\lambda}=\langle\psi|\frac{\partial\mathcal{H}}{\partial\lambda}|\psi\rangle.\]

Show proof

We start from the equation \(E(\lambda)=\langle\psi|\mathcal{H}|\psi\rangle\) and proceed with a differentiation with respect to \(\lambda\)

\[\frac{dE}{d\lambda} = \left\langle\partial_\lambda\psi\right|\mathcal{H}|\psi\rangle + \langle\psi|\frac{\partial\mathcal{H}}{\partial\lambda}|\psi\rangle + \langle\psi|\mathcal{H}\left|\partial_\lambda\psi\right\rangle.\]

Acting with \(\mathcal{H}\) on the bra and ket pulls out a factor \(E(\lambda)\) from the first and third terms, giving

\[\frac{dE}{d\lambda} = E\,\frac{d}{d\lambda}\langle\psi|\psi\rangle + \langle\psi|\frac{\partial\mathcal{H}}{\partial\lambda}|\psi\rangle.\]

The first term vanishes by normalisation, in fact \(\langle\psi|\psi\rangle=1\) for all \(\lambda\). We thus obtain

\[\frac{dE}{d\lambda} = \langle\psi(\lambda)|\,\frac{\partial\mathcal{H}}{\partial\lambda}\,|\psi(\lambda)\rangle.\]

The generalisation to vector parameters such as \(\mathbf{k}\) is obvious. \(\blacksquare\)

Where does the non-degeneracy hypothesis enter? In practice, it guarantees that \(\psi(\lambda)\) is differentiable. More clues will be given when the more general degenerate case is discussed (see Berry phase section).

3.1.3 Reduced, extended and repeated schemes

The intrinsic ambiguity in the definition of \(\mathbf{k}\) leads to equivalent, yet visually very different, “schemes” for displaying the band structure \(E_n(\mathbf{k})\):

• Extended zone scheme. Assign each band \(n\) to the \(n\)-th Brillouin zone by “unfolding” it by \(\mathbf{g}\in\mathcal{RL}\). This representation makes the link to the free-electron parabola most transparent.

• Reduced zone scheme. Restrict \(\mathbf{k}\) to the first Brillouin zone. Each band index \(n\) gives a separate branch \(E_n(\mathbf{k})\). This is the most common representation in solid-state physics.

• Repeated zone scheme. Plot \(E_n(\mathbf{k})\) for all \(n\) over all \(\mathbf{k}\)-space, exploiting the periodicity. Visually redundant and slightly confusing, but useful for seeing how bands connect across zone boundaries.

import { renderZoneSchemes } from "./_interactive/ch03-zone-schemes.js"
{
  const container = document.createElement("div");
  container.style.cssText = "width:100%;";
  yield container;
  renderZoneSchemes(container);
}

Figure 3.1: Zone schemes for a 1D periodic potential \(V(x)=2V_1\cos(2\pi x/a)\). Left: extended/reduced zone schemes: the dashed grey parabola is the free-electron dispersion \(E = \hbar^2 k^2/(2m)\). Right: repeated zone scheme. Use the slider to control the gap size (strength of the periodic potential).

CautionIs this really a “folding”?!

The term “folding” is very very common but somewhat unfortunate, even if visually and intuitively compelling. In reality, the “folding” operation corresponds to a translation by a reciprocal lattice vector, which looks like a folding only due to the presence of Kramers degeneracy, which is responsible for the symmetry \(+\mathbf{k}\leftrightarrows-\mathbf{k}\).

This is fine as long as one is well aware that the “folded” version of a Bloch state at \(k=\pi/a+\varepsilon\) is not a Bloch state at \(k=\pi/a-\varepsilon\), but rather one at \(k=-\pi/a+\varepsilon\). The dispersion looks folded, but equivalent Bloch states do not map according to a folding but according to a translation.

CautionWhat is the n-th Brillouin zone?! ²

This is not a very frequently used concept, but the fact that the FBZ is named “first” clearly suggests that a second, third and so on Brillouin zone can be defined as well.

Recall that the FBZ is the WS cell of the \(\Gamma\) point: it contains all the \(\mathbf{k}\) in reciprocal space that have \(\Gamma\) as their closest point in the \(\mathcal{RL}\), and each such \(\mathbf{k}\) is the smallest representative of its Bloch state. The n-th Brillouin zone is the natural generalisation: it is the set of \(\mathbf{k}\) that have \(\Gamma\) as their \(n\)-th closest point in the \(\mathcal{RL}\). Equivalently, the \(n\)-th BZ is the region one reaches from \(\Gamma\) by crossing exactly \(n-1\) perpendicular bisectors of reciprocal lattice vectors.

In 1D the successive Brillouin zones are just the intervals \([(n-1)\pi/a,\,n\pi/a]\) and \([-n\pi/a,\,-(n-1)\pi/a]\). In higher dimensions the shapes become quickly non-trivial: see Figure 3.2 for a 2D square reciprocal lattice, where the slider progressively fills the first seven zones with shaded colors.

import { renderBrillouinZones } from "./_interactive/ch03-brillouin-zones.js"
{
  const container = document.createElement("div");
  container.style.cssText = "width:100%;";
  yield container;
  renderBrillouinZones(container);
}

Figure 3.2: First seven Brillouin zones of a 2D square reciprocal lattice. Drag the slider to progressively reveal zones \(1,2,\dots,n\). Note how each successive zone, though disconnected into several pieces, can be refolded onto the FBZ.

3.2 How do energy bands emerge?³

Here we want to build intuition on how energy bands emerge from periodicity. We will look at the same problem from complementary angles, using simple 1D models where computations are trivial.

TheoremNo band overlaps in one dimension

1D is somewhat special and simpler:

1. \(\mathcal{H}\psi = E\psi\) is a second-order ODE and admits no more than two linearly independent solutions;
2. these are indeed typically two \(E=E_n(k)=E_n(-k)\) (Equation 3.6).

Consequently, bands in 1D cannot overlap: at every energy there are at most two states, and no additional degeneracies are possible. This is very untrue in 2D and 3D, where band crossings and degeneracies are common.

3.2.1 The transfer matrix method

Perspective n.1. Here we want to see how energy bands emerge from constructive and destructive interference from periodically arranged scattering centers. Even in a simplified 1D case, taking into account all the possible scattering paths is not obvious: considering all the combinations of transmission and reflection events over each barrier leads to a rapidly divergent set of options.

import { renderScatteringVsTM } from "./_interactive/ch03-scattering-vs-tm.js"
{
  const container = document.createElement("div");
  container.style.cssText = "width:100%;position:relative;margin:0 auto;";
  yield container;
  renderScatteringVsTM(container);
}

Figure 3.3: Left: Scattering picture, given the incoming particle beams we calculate the amplitudes of outgoing ones. Right: transfer picture, we describe what happens on the right side of the barrier as a function of what happens on the left side. Below: the key advantage is in the calculation of the cumulative effects of multiple barriers on particle transmission and reflection.

A powerful and general approach for 1D scattering problems is the transfer matrix formalism. Consider a region where the potential has compact support (i.e. \(V(x)=0\) outside \([0,a]\)). A plane wave of energy \(E\) incident from the left is partially reflected and partially transmitted. We define the scattering amplitudes \(r\) (reflection) and \(t\) (transmission):

\[\psi(x) = \begin{cases} e^{iqx} + r\,e^{-iqx} & x < 0 \\ t\,e^{iqx} & x > a \end{cases} \tag{3.8}\]

where \(E(q)=\hbar^2q^2/2m\). As mentioned, this approach becomes very hard to use in the presence of multiple barriers since the number of possible transmission and reflection paths becomes infinite and quite complicated to analyse (see Figure 3.3: one barrier is obvious, two feasible, beyond two… good luck). The transfer matrix \(\mathcal{M}\) takes a different approach and relates the coefficients of the plane waves on the two sides of the barrier:

\[\begin{pmatrix} A_R \\ B_R \end{pmatrix} = \mathcal{M} \begin{pmatrix} A_L \\ B_L \end{pmatrix} = \begin{pmatrix} m_{11} & m_{12} \\ m_{21} & m_{22} \end{pmatrix} \begin{pmatrix} A_L \\ B_L \end{pmatrix}, \tag{3.9}\]

where \(A\) and \(B\) are the right-moving and left-moving amplitudes, respectively, i.e. we assume

\[\psi(x) = \begin{cases} A_L\,e^{iqx} + B_L\,e^{-iqx} & x < 0 \\ A_R\,e^{iqx} + B_R\,e^{-iqx} & x > a. \end{cases}\]

The key advantage of the transfer matrix is its multiplicative composition: the cumulative effect of \(n\) barriers described by the matrices \(\mathcal{M}_i\) can be simply calculated by a matrix product

\[\mathcal{M}_\text{tot} = \mathcal{M}_n\,\mathcal{M}_{n-1} \cdots \mathcal{M}_2\mathcal{M}_1.\]

TheoremKey properties of transfer matrices

Transfer matrices are unimodular \(\det\mathcal{M}=1\) and their matrix elements are linked by conjugation

\[\begin{cases} m_{22} &= m_{11}^*\\ m_{21} &= m_{12}^* \end{cases} \tag{3.10}\]

This is required by time reversal symmetry and probability conservation.

Show proof

Time reversal symmetry. If the system is symmetric under this operation and \(\psi(x)\) is a good solution, then also \(\psi^*\) is a good solution. This operation maps \(A_L\to B_L^*\), and similarly all the other coefficients, thus leading to

\[\begin{pmatrix} B_R^* \\ A_R^* \end{pmatrix} = \begin{pmatrix} m_{11} & m_{12} \\ m_{21} & m_{22} \end{pmatrix} \begin{pmatrix} B_L^* \\ A_L^* \end{pmatrix},\]

that has to be compared with the conjugate of the original Equation 3.9

\[\begin{pmatrix} A_R^* \\ B_R^* \end{pmatrix} = \begin{pmatrix} m_{11}^* & m_{12}^* \\ m_{21}^* & m_{22}^* \end{pmatrix} \begin{pmatrix} A_L^* \\ B_L^* \end{pmatrix},\]

leading finally to \(m_{22}=m_{11}^*\) and \(m_{21}=m_{12}^*\).

Probability current density. Probability current densities on the two sides

\[\begin{aligned} J_L &= (|A_L|^2-|B_L|^2)v_g\\ J_R &= (|A_R|^2-|B_R|^2)v_g \end{aligned}\]

where \(v_g\) is the group velocity, have to be the same; this implies

\[\begin{aligned} |A_R|^2-|B_R|^2 &= |m_{11}A_L+m_{12}B_L|^2-|m_{12}^*A_L+m_{11}^*B_L|^2\\ &= |m_{11}|^2(|A_L|^2-|B_L|^2)-|m_{12}|^2(|A_L|^2-|B_L|^2)\\ &= \det\mathcal{M}\,(|A_L|^2-|B_L|^2), \end{aligned}\]

thus the matrix has to be unimodular: \(\det\mathcal{M}=1.\) \(\blacksquare\)

Starting from the scattering setup in Equation 3.8, the amplitudes on the left side are \((A_L,B_L) = (1,r)\), while on the right side \((A_R,B_R) = (t,0)\) (no wave coming from the right). From the definition Equation 3.9 we have

\[\begin{cases} t &= m_{11}+m_{12}r\\ 0 &= m_{21}+m_{22}r \end{cases}\]

and thus \(r = -m_{21}/m_{22}\), while the first equation yields

\[t = m_{11} + m_{12}\,r = m_{11} - m_{12}\,\frac{m_{21}}{m_{22}} = \frac{\det\mathcal{M}}{m_{22}}\]

i.e., using \(\det\mathcal{M} = 1\), we have \(m_{22}=m^*_{11} = 1/t\) and finally

\[\mathcal{M} = \begin{pmatrix} 1/t^* & -r^*/t^* \\ -r/t & 1/t \end{pmatrix}. \tag{3.11}\]

Note that the formalism becomes singular in the limit of zero transmission, which makes sense (how can you tell what happens beyond an infinite barrier?). In the following we explore a few relevant limit cases.

CautionCase study - Resonant tunnelling

Assume we have one barrier with \(t = |t|\,e^{i\phi_t}\) and \(r\), with \(T = |t|^2\) and \(R = |r|^2 = 1 - T\). See Figure 3.4, we now add a second identical barrier. Its transfer matrix is easy to calculate: an incoming wave will have to travel a further \(a\) to reach the barrier and \(a\) again after reflection. Differently, any phase delay before the barrier becomes a phase anticipation in the transmission beyond the barrier. In short,

\[ \begin{cases} r &\to re^{2iqa}\\ t &\to t. \end{cases} \]

In terms of transfer matrix formalism, the matrix translated by \(a\) becomes

\[\mathcal{M}(a) = \begin{pmatrix} m_{11} & m_{12}e^{-2iqa} \\ m_{21}e^{+2iqa} & m_{22} \end{pmatrix}.\]

The total transfer matrix is \(\mathcal{M}_\text{DB} = \mathcal{M}(a)\mathcal{M}(0)\), but we only need to calculate one element

\[m_{\text{DB},11} = \frac{1}{t_\text{DB}^*} = m_{11}^2+m_{12}m_{21}e^{-2iqa} = \frac{e^{2i\phi_t}+Re^{-2iqa}}{T}\]

to finally get the transmission

\[T_\text{DB} = \frac{T^2}{|1+Re^{-2i\varphi}|^2} \tag{3.12}\]

where \(\varphi = \phi_t+qa\) is the total phase accumulated in a trip between the barriers (propagation phase \(qa\) plus the scattering phase \(\phi_t\)). The denominator oscillates between \((1-R)^2\) and \((1+R)^2\) as \(\varphi\) varies:

\[\begin{cases} T_\text{DB}^\text{max} = \cfrac{T^2}{(1-R)^2} = 1 \quad &\text{when } 2\varphi = (2n+1)\pi, \\ T_\text{DB}^\text{min} = \cfrac{T^2}{(1+R)^2} \quad &\text{when } 2\varphi = 2n\pi. \end{cases}\]

At resonance (\(T_\text{DB} = 1\)), one can show that the wave function between the barriers is strongly amplified due to trapping in a quasi-bound state. Off resonance, reflection dominates and \(T_\text{DB} \ll T\).

How can two barriers transmit more than one?

Beyond math, one can intuitively wonder how it is possible that two barriers are much less opaque than a single one. The reason might become less obscure by looking at numerical (and interactive/didactic: do not always expect perfect transmission due to the limitations of the reactive code behind the figure!) results in Figure 3.4: select the first resonance in the double barrier and observe what happens in the following two cases.

• We use a finite wave packet using the “pulse” button. Most of the packet is still reflected, since its energy distribution is way too broad to be transmitted through the sharp resonance of the double barrier;
• We use an infinite train of waves using the “on” button. Here resonant tunnelling emerges, but it clearly relies on “pumping” a resonant state trapped between the barriers. The state is “leaky” and continuously decays by tunnelling out on the left and right sides of the double barrier structure. At resonance, conditions are such that the “left leak” coherently cancels the reflection from the first barrier, while building full transmission on the right side. The resonant state decay becomes evident once the train is turned “off”.

Exercise. In the limit of \(T\to 0\), assuming \(\varphi\approx qa\) and free propagation between the barriers, it is possible to show that \(T_\text{DB}(E)\) tends to a Lorentzian peak, with a width \(\Delta=\hbar/\tau\) where \(\tau\) is the lifetime of the quasi-bound state bouncing between the two barriers. The demonstration of this fact is left to the reader.

import { renderSplitStep } from "./_interactive/ch03-splitstep.js"
{
  const W = Math.min(width, 960);
  const container = document.createElement("div");
  container.style.cssText = "width:100%; max-width:" + W + "px; margin:0 auto;";
  yield container;
  renderSplitStep(container);
}

Figure 3.4: Time-dependent split-step Fourier simulation of a 1D Schrödinger wavepacket scattering off a selectable potential. Left panel: the transmission probability \(T(E)\) obtained by transfer matrix on the discretised \(V(x)\). Right panel: \(|\psi|\) (yellow band), \(\operatorname{Re}\psi\) (blue), \(\operatorname{Im}\psi\) (red); the violet band is \(V(x)\).

CautionCase study - \(N\)-barrier limit

The case of \(N\) identical barriers separated by equal distances \(a\) is also quite interesting and directly related to interesting devices such as Bragg reflectors. The problem can be analysed using the transfer matrix formalism, but here we take a different route. In the region \(x\in[0,Na]\), we have a periodic potential and whenever the energy \(E\) belongs to a band, a generic state can be described as

\[\psi(x) = \alpha\psi_{+k}(x)+\beta\psi_{-k}(x)\]

with \(E(k)=E\) and we are excluding the band edges, where only a single state is available. At the \(x=Na\) interface, two states are sufficient to match with any wavefunction present at \(x>Na\), including in particular a transmitted wavefunction \(\psi_R(x)=e^{iq(x-Na)}\). The matching requires two conditions

\[\begin{cases} \psi(Na^+) = 1 &= \alpha\psi_{+k}(Na^{-}) + \beta\psi_{-k}(Na^{-}) = \alpha\psi_{+k}(0)e^{+ikNa} + \beta\psi_{-k}(0)e^{-ikNa}\\ \psi'(Na^+) = iq &= \alpha\psi'_{+k}(Na^{-}) + \beta\psi'_{-k}(Na^{-}) = \alpha\psi'_{+k}(0)e^{+ikNa} + \beta\psi'_{-k}(0)e^{-ikNa} \end{cases}\]

The same matching condition is possible for \(\psi_L(x)=e^{iqx+i\varphi}\) at \(x=0\), as long as \(e^{+ikNa}=e^{-ikNa}=e^{i\varphi}\): under these conditions, no reflection is needed and perfect transmission is achieved. This simply requires \(2kNa\) to be a multiple of \(2\pi\): we have a number of options but we want to exclude: (1) \(k\) values outside the FBZ since they would be clearly redundant; (2) negative \(k\) values since this would be a clear double counting; (3) \(k=0\) and \(k=\pi/a\) since we want to exclude band edges where we most likely have a single solution (ruling out zero-gap limit cases). In conclusion, valid \(k\) values and resonance energies are

\[k_j = \frac{j\pi}{Na}\quad\implies\quad E_{nj}=E_n(k_j)\]

for \(j\) going from \(1\) to \(N-1\). Thus, for every band (as labeled by index \(n\)) we expect to see \(N-1\) resonances at energies \(E_n(k_j)\). See interactive numerical results in Figure 3.5.

import { renderNBarriers } from "./_interactive/ch03-n-barriers.js"
{
  const container = document.createElement("div");
  container.style.cssText = "width:100%;";
  yield container;
  renderNBarriers(container);
}

Figure 3.5: Transmission probability through \(N\) identical square barriers as a function of energy. As \(N\) grows, the discrete resonances merge into the allowed-band continuum of the infinite-periodic limit. Adjust the barrier height and number of barriers with the sliders. Activate the wavefunction visualization and resonance snap using the checkboxes.

CautionCase study - Infinite barriers

For a genuinely infinite periodic potential \(V(x+a) = V(x)\), there is no “incident” wave and no “transmitted” wave and every plane-wave component is simultaneously the outgoing amplitude of one unit cell and the incoming amplitude of the next. The transfer-matrix formalism still applies, but it is now combined with the Bloch condition to select the physically admissible solutions.

Here we consider a 1D potential \(V(x)\) of periodicity \(a\), and call \(V_0(x)\) the potential restricted to \(x\in[0,a]\), see Figure 3.6. For simplicity, the minima of \(V(x)\) are assumed to be zero and located at positions \(t_n=na\); our transfer matrix formalism can then be applied directly to describe states before and after the barrier

\[\begin{cases} \psi(x) = \psi_L(x) = A_Le^{iqx}+B_Le^{-iqx}\, \quad &x \le 0,\\ \psi(x) = \psi_R(x) = A_Re^{iqx}+B_Re^{-iqx}\, \quad &x \ge a, \end{cases}\]

import { renderPeriodicScattering } from "./_interactive/ch03-periodic-scattering.js"
{
  const container = document.createElement("div");
  container.style.cssText = "width:100%;";
  yield container;
  renderPeriodicScattering(container);
}

Figure 3.6: (a) Periodic potential \(V(x+a) = V(x)\) extending over the whole real line. (b) Restriction to one unit cell: \(V_0(x)\) scatters the plane-wave components \(A_L\,e^{iqx}\), \(B_L\,e^{-iqx}\) on its left into \(A_R\,e^{iqx}\), \(B_R\,e^{-iqx}\) on its right. Solid arrows mark right-moving waves, dashed arrows left-moving ones.

where as usual \(E=\hbar^2q^2/2m\). Coefficients \((A_R,B_R)\) can be calculated using the matrix \(\mathcal{M}\) associated to the potential \(V_0(x)\). On the other hand — thanks to Bloch’s theorem and to the fact that we can look for a solution which is also a Bloch state \(\psi_k(x)\) with some yet unknown quasi-momentum \(k\) — they can also be calculated using the Bloch condition \(\psi(a^+)=\psi(0^-)\,e^{ika}\), see Figure 3.6. In matrix format

\[ \begin{aligned} \begin{pmatrix} A_R\,e^{iqa} \\ B_R\,e^{-iqa} \end{pmatrix} &= \mathcal{P}(a)\mathcal{M} \begin{pmatrix} A_L \\ B_L \end{pmatrix} = \begin{pmatrix} m_{11}\,e^{iqa} & m_{12}\,e^{iqa} \\ m_{21}\,e^{-iqa} & m_{22}\,e^{-iqa} \end{pmatrix} \begin{pmatrix} A_L \\ B_L \end{pmatrix}\\ &= e^{ika} \begin{pmatrix} A_L \\ B_L \end{pmatrix} \end{aligned} \tag{3.13}\]

where \(\mathcal{P}(a)\) is a propagation matrix by a distance \(a\), redefining the origin. The matrix \(\mathcal{P}(a)\mathcal{M}\) in the equation has eigenvalue \(\lambda_1=e^{ika}\). Since it is unimodular (and we expect time-reversal symmetry) it must also have an eigenvalue \(\lambda_2=e^{-ika}\), thus the trace \(\lambda_1+\lambda_2=2\cos(ka)\) has to satisfy

\[\mathrm{Tr}\left[\mathcal{P}(a)\mathcal{M}\right] = 2\,\mathrm{Re}\left[\frac{e^{iqa}}{t^*}\right] = 2\cos(ka)\]

and finally

\[F(qa) = \cos(ka). \tag{3.14}\]

where we defined \(F(qa)=\cos(\phi_t+qa)/\sqrt{T}\) and we encounter again \(\varphi(E)=\phi_t(E)+q(E)a\), as in the double-barrier problem. This is the implicit equation for \(E_n(k)\): since \(T(E)\le 1\) the function \(F(qa)\) is an oscillation of amplitude \(\ge 1\) that has to match \(\cos(ka)\in[-1,+1]\). The equation admits a solution only where \(|F(qa)|\le 1\), and the inversion yields a single \(E_n(k)\) band dispersion. See Figure 3.7 for an interactive view of all four pieces — \(F(qa)\), \(\cos(ka)\), the free-electron parabola \(E(q)\), and the resulting band structure \(E_n(k)\) — for a square barrier of width \(b\) and height \(V_0\).

import { renderKronigPenney } from "./_interactive/ch03-kronig-penney.js"
{
  const container = document.createElement("div");
  container.style.cssText = "width:100%;";
  yield container;
  renderKronigPenney(container);
}

Figure 3.7: Kronig–Penney for a square barrier of width \(b\) and height \(V_0\) (dimensionless units \(\hbar^2/(2ma^2) = 1\), \(a = 1\)). Top-left: \(F(qa)\) vs \(qa\), with the allowed-energy strip \(|F|\le 1\) in grey and the band edges \(F=\pm 1\) as dashed dark-grey lines. Top-right: \(\cos(ka)\) vs \(ka\) in the reduced zone \([-\pi/a,\pi/a]\). Bottom-left: the free-electron parabola \(E = q^2\) with allowed bands shaded grey. Bottom-right: the resulting band structure \(E_n(k)\). Move the two sliders to change \(b/a\) and \(V_0\).

CautionAn historical example: the Kronig–Penney model

The Kronig–Penney model (R. de L. Kronig & W. G. Penney, 1931) was the first exactly-solvable 1D periodic model and historically the cleanest demonstration that allowed bands and forbidden gaps are an inevitable consequence of periodicity. Each unit cell consists of a free well of width \(w\) followed by a rectangular barrier of width \(b\) and height \(V_0\), with the periodicity

\[a = b + w.\]

Take a single cell with the well on \([-w, 0]\) and the barrier on \([0, b]\). Inside the well \(V = 0\) and the wavefunction is a plane wave; inside the barrier it decays exponentially for \(E < V_0\) (or oscillates for \(E > V_0\)):

\[\psi(x) = \begin{cases} A\,e^{iqx} + B\,e^{-iqx}, & -w \le x \le 0, \quad q = \sqrt{2mE}/\hbar, \\[4pt] C\,e^{\beta x} + D\,e^{-\beta x}, & 0 \le x \le b, \quad \beta = \sqrt{2m(V_0-E)}/\hbar. \end{cases}\]

Continuity of \(\psi\) and \(\psi'\) at \(x = 0\) together with the Bloch conditions at the cell edges, \(\psi(b) = e^{ika}\,\psi(-w)\) and \(\psi'(b) = e^{ika}\,\psi'(-w)\), give a homogeneous \(4\times 4\) linear system for \((A, B, C, D)\) that can be cast in matrix form

\[\begin{pmatrix} 1 & 1 & -1 & -1 \\[2pt] iq & -iq & -\beta & \beta \\[2pt] e^{-iqw}\,e^{ika} & e^{+iqw}\,e^{ika} & -e^{+\beta b} & -e^{-\beta b} \\[2pt] iq\,e^{-iqw}\,e^{ika} & -iq\,e^{+iqw}\,e^{ika} & -\beta\,e^{+\beta b} & \beta\,e^{-\beta b} \end{pmatrix} \begin{pmatrix} A \\ B \\ C \\ D \end{pmatrix} = \begin{pmatrix} 0 \\ 0 \\ 0 \\ 0 \end{pmatrix}.\]

A non-trivial solution requires the determinant of the matrix to vanish: after some algebra

\[F(qa) \;\equiv\; \cos(qw)\cosh(\beta b) + \frac{\beta^2 - q^2}{2q\beta}\sinh(\beta b)\sin(qw) \;=\; \cos(ka), \tag{3.15}\]

identical in structure to Equation 3.14 and admitting the same band/gap interpretation \(|F(qa)|\le 1\). See Figure 3.7 for an interactive view: it implements exactly this expression — the slider for \(b\) controls the barrier width and the slider for \(V_0\) its height.

A particularly elegant limit is \(b \to 0\), \(V_0 \to \infty\) with the product \(V_0 b \equiv P\hbar^2/(ma)\) kept fixed (Dirac-delta barriers). The Kronig–Penney relation then collapses to

\[\cos(ka) = \cos(qa) + \frac{P}{qa}\sin(qa), \tag{3.16}\]

controlled by the single dimensionless parameter \(P\): for \(P = 0\) one recovers free electrons, while \(P \to \infty\) gives isolated quantum wells (flat bands).

3.2.2 Nearly-free electron model

Perspective n.2. We consider how a weak periodic potential \(V(\mathbf{r})\) can perturb the parabolic dispersion of a free electron. This is called the nearly-free electron (NFE) approach and starts from an “empty lattice” limit and a folded free-electron parabola in the FBZ. We know that any Bloch state with a given \(\mathbf{k}\) can be expanded as a linear combination of plane waves in the \(S_\mathbf{k}\) subspace

\[ w_{\mathbf{k}+\mathbf{g}}(\mathbf{r}) = \frac{1}{\sqrt{V}}e^{i(\mathbf{k}+\mathbf{g})\cdot\mathbf{r}} \]

with energy \(E_{\mathbf{k}+\mathbf{g}}=\hbar^2(\mathbf{k}+\mathbf{g})^2/2m\); here \(V\) is the normalisation volume and \(\mathbf{g}\in\mathcal{RL}\). The generic matrix element in the \(S_\mathbf{k}\) subspace is

\[ \langle w_{\mathbf{k}+\mathbf{h}}|\mathcal{H}|w_{\mathbf{k}+\mathbf{g}}\rangle = E_{\mathbf{k}+\mathbf{g}}\delta_{\mathbf{g},\mathbf{h}} + \frac{1}{V}\int V(\mathbf{r})e^{i(\mathbf{g}-\mathbf{h})\cdot\mathbf{r}}d^3\mathbf{r} = E_{\mathbf{k}+\mathbf{g}}\delta_{\mathbf{g},\mathbf{h}} + V_{\mathbf{h}-\mathbf{g}} \]

where the off-diagonal terms are the Fourier components of \(V(\mathbf{r})\). The strongest, qualitative changes will occur when \(E_{\mathbf{k}+\mathbf{g}} \approx E_{\mathbf{k}+\mathbf{h}}\) and degenerate perturbation theory is required. For simplicity, here we move to 1D at \(k\approx\pi/a\), where the \(g=0\) parabola and the folded \(h=-2\pi/a\) collide. The projection on the degenerate subspace is

\[\mathcal{H}_ ext{deg}=\begin{pmatrix} E_1 & V_1 \\ V_1^* & E_2 \end{pmatrix}, \tag{3.17}\]

where we use \(k=\pi/a+\Delta k\) and \(E_{1/2} = \hbar^2(\pi/a\pm\Delta k)^2/2m\). Setting \(\det(\mathcal{H}_ ext{deg}-E)=0\) leads to

\[ \begin{aligned} E_\pm(k) &= \frac{E_1+E_2}{2} \pm \sqrt{\frac{(E_1-E_2)^2}{4} + |V_1|^2}\\ &= E_0+\frac{\hbar^2\Delta k^2}{2m} \pm \sqrt{4E_0\frac{\hbar^2\Delta k^2}{2m} + |V_1|^2}\\ &\approx E_0\pm|V_1|+\frac{\hbar^2\Delta k^2}{2m}\left(1\pm\frac{2E_0}{|V_1|}\right) \end{aligned} \tag{3.18}\]

where \(E_0=\hbar^2(\pi/a)^2/2m\) is the crossing energy. Thus we get a gap \(2|V_1|\) with an approximately hyperbolic dependence on \(\Delta k\). Note that this is just a perturbation: at higher order \(V_1\) will open a gap also at the next crossing, where the leading term is controlled by \(V_2\) (see Figure 3.8, for an unrestricted numerical diagonalisation).

import { renderNFEbands } from "./_interactive/ch03-nfe-bands.js"
{
  const container = document.createElement("div");
  container.style.cssText = "width:100%; position:relative;";
  yield container;
  renderNFEbands(container);
}

Figure 3.8: Nearly-free electron model in 1D. Left: band structure \(E_n(k)\) in the reduced zone scheme (solid blue), with the folded free-electron parabolas shown dashed for reference. Right: the periodic potential \(V(x) = \sum_{n=1}^{2} 2|V_n|\cos(2\pi n\,x/a + \varphi_n)\) over \(10\) periods.

3.3 Filling the bands

As discussed later, the behavior of a crystal is profoundly impacted by the total or partial filling of its bands. It is thus very important to tell how many states are available in a band.

3.3.1 Number of states⁴

The standard counting method consists in setting periodic conditions at the crystal boundaries. While it might sound more natural to set \(\psi=0\) at the physical edge of a finite crystal, this would require mixing Bloch states with different \(\mathbf{k}\), while fundamentally not changing the final conclusions.

Definition 3.4: Born–von Kármán conditions and number of states in a band⁵

For a finite crystal of \(N_i\) cells along direction \(i\) the Born-von Kármán (BvK) conditions consist in setting

\[\psi(\mathbf{r}+N_i\mathbf{t}_i) = \psi(\mathbf{r}),\]

which, combined with the Bloch condition Equation 3.1, requires \(e^{iN_i\mathbf{k}\cdot\mathbf{t}_i} = 1\), so

\[\mathbf{k} = \sum_i \frac{m_i}{N_i}\,\mathbf{g}_i, \qquad m_i = 0, 1, \dots, N_i - 1,\]

where \(\mathbf{g}_i\) are the reciprocal-lattice primitive vectors and \(m_i\) is restricted to non-equivalent \(\mathbf{k}\) values. Finally there are \(N = N_1 N_2 N_3\) inequivalent \(\mathbf{k}\) points, one per cell. Each band therefore contributes \(N\) states, to be multiplied by any additional degeneracy (typically, we need at least a factor \(2\) because of spin degeneracy).

CautionI am not convinced… what if I still want to set \(\psi=0\) at the edges?

This is a reasonable doubt given that crystals do not actually have periodic boundary conditions: an electron hitting one end of a solid will certainly not disappear and pop up at the other end. So, yes, it is possible to set \(\psi=0\) at the edges, but this yields essentially the same state count, with a more complicate basis.

Let us work in 1D and impose \(\psi(0)=\psi(Na)=0\) for a given energy \(E\). The general solution will be a superposition of Bloch states

\[\psi(x) = \alpha\psi_k(x)+\beta\psi_{-k}(x)\]

where \(E=E(k)\) and we use \(\psi_{-k}(x)=\psi_k^*(x)\). The two coefficients \(\alpha\) and \(\beta\) clearly need to have the same modulus and a suitable relative phase so that \(\alpha\psi_k(0)+\beta\psi_k^*(0)=0\). We now require the same to happen at \(x=Na\), but using the Bloch relation

\[\psi(Na) = \alpha e^{ikNa}\psi_k(0)+\beta e^{-ikNa}\psi_k^*(0)\]

it is easy to show that this requires \(e^{ikNa}=e^{-ikNa}\), i.e. \(e^{2ikNa}=1\), and thus

\[k=\frac{m}{2N}\frac{2\pi}{a}\]

where \(m=1, \dots 2N-1\) is an integer and where we need to impose \(k>0\) to avoid a double counting. So, finally, we have to restrict the sign of \(k\) and we have twice the density of states in the \(k\)-space: the total number of states in the band is still essentially the same. To be precise, here we actually have \(N-1\) options instead of \(N\), because a still electron at the bottom or top of a band cannot sastify the boundary conditions. This difference is however irrelevant in the \(N\to\infty\) limit of a large crystal.

3.3.2 Density of states

The BvK conditions indicate that the volume for every allowed \(\mathbf{k}\)-value is

\[ \frac{|\mathbf{g}_1\cdot(\mathbf{g}_2\times\mathbf{g}_3)|}{N} = \frac{\Omega_k}{N} = \frac{(2\pi)^3}{N\Omega} = \frac{(2\pi)^3}{V} \]

where \(V\) is the crystal volume, and similar expressions can be obtained for lower dimensions using surface \(S\) or length \(L\). Any additional degeneracy, such as spin, will require an additional integer factor often marked as \(g\). Taking the continuous limit for \(N\to\infty\), any sum over the band states can be written as an integral

\[ g\sum_{\mathbf{k}\in \mathrm{FBZ}} \quad\to\quad \frac{gV}{(2\pi)^3}\int d^3k. \]

Whenever the sum only depends on the energy, it is convenient to introduce the density of states (DoS) — called \(D(E)\) or \(g(E)\) depending on books — which counts the states \(dN=D(E)dE\) in a given energy slice \([E, E+dE]\). The number of states up to energy \(E\) for a given \(E(\mathbf{k})\) is

\[N(E) = \frac{gV}{(2\pi)^3}\int \theta(E-E(\mathbf{k}))d^3k, \tag{3.19}\]

and thus the DoS is defined as

\[D(E) = \frac{dN}{dE} = \frac{gV}{(2\pi)^3}\int \delta(E - E(\mathbf{k}))d^3k. \tag{3.20}\]

Reconsider now states in the \([E,E+dE]\) energy range: their number is proportional to the distance between the two isoenergetic surfaces \(E(\mathbf{k})=E\) and \(E(\mathbf{k})=E+dE\). This depends on how “fast” \(E(\mathbf{k})\) changes in \(\mathbf{k}\), which depends on the local gradient: \(dE=\nabla_\mathbf{k}E\cdot d\mathbf{k}\). Thus, the volume spanned by a portion \(dS_k\) of the \(E(\mathbf{k})\) surface will be \(dS_kdk_\perp = dS_kdE/|\nabla_\mathbf{k} E|\), where \(dk_\perp\) is perpendicular to the energy surface and parallel to \(\nabla_\mathbf{k}E\). Finally, we have a new definition

\[D(E) = \frac{gV}{(2\pi)^3}\int_{E(\mathbf{k})=E} \frac{dS_k}{|\nabla_\mathbf{k} E|} \tag{3.21}\]

suggesting that whenever \(\nabla_\mathbf{k}E=\mathbf{0}\) something interesting can occur.

Definition 3.5: Van Hove singularities

At stationary points of the band dispersion with \(\nabla_\mathbf{k}E = \mathbf{0}\), the integrand in Equation 3.21 diverges. This leads to so-called van Hove singularities or critical points in the DoS.

When \(E(\mathbf{k})\) is non-degenerate, it is a smooth function of \(\mathbf{k}\) and critical points can be analysed using a Taylor expansion. Note a singularity will occur at a given energy \(E_0\) and at one or more \(\mathbf{k}_0\), but the \(\mathbf{k}_0\) position is irrelevant in the integral so it is assumed \(\mathbf{k}_0=\mathbf{0}\) and any possible multiplicity is neglected. Stopping at the second order, we have

\[ E(\mathbf{k}) = E_0+\frac{1}{2}\,\frac{\partial^2 E}{\partial k_i\partial k_j}k_ik_j. \]

Qualitatively different singularities can be obtained depending on whether we are looking at a band maximum, minimum or saddle point. This depends on the signs of the eigenvalues of the Hessian matrix. See in Figure 3.9 all the available options for 1D, 2D and 3D dispersions (see below for an analysis of selected cases).

import { renderVanHove } from "./_interactive/ch03-van-hove.js"
{
  const W = Math.min(width, 760);
  const container = document.createElement("div");
  container.style.cssText =
    "width:100%; max-width:" + W + "px; margin:0 auto;";
  yield container;
  renderVanHove(container);
}

Figure 3.9: Critical points can be classified based on the signature of the band Hessian. Left: Select the Hessian signature. Right: Visualization of the resulting singularity in the selected dimension.

CautionA few DoS calculations

In the following calculations, the Hessian will be assumed to be diagonal and parametrized in terms of \(m^*\) parameters, which we will later call band mass or effective mass

\[E(\mathbf{k}) = E(\mathbf{0})\pm\frac{\hbar^2k_x^2}{2m_x^*}\pm\frac{\hbar^2k_y^2}{2m_y^*}\pm\frac{\hbar^2k_z^2}{2m_z^*}.\]

1D case. In the parabolic approximation \(E(k)=E(0)\pm\hbar^2k^2/2m^*\), the critical points \(P_0\) and \(P_1\) give rise to the following DoS

\[ \begin{aligned} D_{\mathrm{P0}}(E) &= gL\frac{\sqrt{m^*}}{\pi\hbar}\,\frac{1}{\sqrt{2(E-E_0)}}\\ D_{\mathrm{P1}}(E) &= gL\frac{\sqrt{m^*}}{\pi\hbar}\,\frac{1}{\sqrt{2(E_1-E)}}. \end{aligned} \]

Note here that \(m^*\) at minimum and maximum can easily be different.

Show proof

The calculation is easy to perform starting from \(N(E)\): \(k(E)=\sqrt{2m^*(E-E_0)}/\hbar\), occupied states cover the \([-k(E),+k(E)]\) interval and thus we have \(N(E)=gLk(E)/\pi\). Finally

\[D(E) = \frac{dN}{dE} = gL\frac{\sqrt{m^*}}{\pi\hbar}\,\frac{1}{\sqrt{2(E-E_0)}}.\]

The DoS for the band maximum can be obtained by symmetry.

2D case. The result for an isotropic minimum \(m^*=m_x^*=m_y^*\) is rather simple to calculate and yields

\[ \begin{aligned} D_{\mathrm{Q0}}(E) &= \frac{gSm^*}{2\pi\hbar^2}\,\theta(E-E_0)\\ D_{\mathrm{Q2}}(E) &= \frac{gSm^*}{2\pi\hbar^2}\,\theta(E_2-E) \end{aligned} \]

while in the most general case we need to replace \(m^*\) with \(\sqrt{m_x^*m_y^*}\) or, an even more generic and coordinate-independent expression, \(\sqrt{\det(m^*)}\) where \(m^*\) is the mass tensor (more details will come in the next chapters).

Show proof

We can again start from \(N(E)\): occupied states cover a circle with radius \(k(E)\), i.e. \(\pi k(E)^2 = 2\pi m^*(E-E_0)/\hbar^2\).

\[D(E) = \frac{dN}{dE} = \frac{gSm^*}{2\pi\hbar^2}\,\theta(E-E_0).\]

In the most general case of an anisotropic mass it is convenient to write \(E(\mathbf{q})=E_0+q_x^2+q_y^2\) where \(q_x=\hbar k_x/\sqrt{2m_x^*}\) and \(q_y=\hbar k_y/\sqrt{2m_y^*}\). The occupied states cover the circle \(|\mathbf{q}|<\sqrt{E-E_0}\), which can be easily converted to the original \(\mathbf{k}\) coordinates, yielding an occupied area

\[\frac{2\pi \sqrt{m_x^*m_y^*}}{\hbar^2}\,(E-E_0).\]

As in the 1D case, the DoS in \(Q_2\) can be derived by symmetry.

The saddle point \(Q_1\) requires a slightly less obvious calculation and leads to

\[ D_{\mathrm{Q1}}(E) = A \log\frac{1}{|E-E_1|} \]

where the factor \(A\) is not universal, even in the case of a parabolic approximation. See Section 2.7 of Grosso and Pastori Parravicini (2014) for a derivation.

3D case. Again the isotropic case \(m^*=m_x^*=m_y^*=m_z^*\) easily leads to

\[N(E) = \frac{gV}{(2\pi)^3}\frac{4\pi k(E)^3}{3} = \frac{gV}{\pi^2}\frac{\sqrt{2m^{*3}(E-E_0)^3}}{3\hbar^3}\]

and finally

\[D(E) = \frac{gV\sqrt{2m^{*3}(E-E_0)}}{2\pi^2\hbar^3}.\]

At the \(M_1\) saddle point we have

\[D(E) = C_1-\frac{gV\sqrt{2m^{*3}(E_1-E)}}{2\pi^2\hbar^3}.\]

where \(C_1\) is a non-universal constant. Points \(M_2\) and \(M_3\) follow by symmetry. See Section 2.7 of Grosso and Pastori Parravicini (2014) for a derivation of the DoS at the saddle point.

General functional form. Neglecting the exact multiplicative factors, the functional form of any DoS at a band minimum is practically fixed by dimensional analysis and equals

\[D(E)\propto \frac{\sqrt{||m^*||}}{\hbar^n}(E-E_0)^{n/2-1}.\]

where \(n\) is the dimension of the \(\mathbf{k}\)-space under consideration and \(||m^*||=\det(m^*)\) is the determinant of the mass tensor. Indeed, we can get an inverse length \(\sqrt{2m^*(E-E_0)}/\hbar\) for every \(m^*\) eigenvalue. Given the DoS is proportional to the real space volume (or surface or length depending on the dimension), we have a general and coordinate-independent dimensionless parameter \(\sqrt{||m^*||(E-E_0)^n}/\hbar^n\). The additional \(1/(E-E_0)\) factor is the only possible one to recover the correct dimensions for \(D(E)\).

The resulting concept of Fermi surface is developed in Section 5.4; how the filling of the bands classifies solids into metals and insulators is discussed in Section 8.6, once the band dynamics has been introduced.

In conclusion…

ImportantTake home messages…

At the end of this chapter you should know…

• Bloch theorem. Definition of Bloch state and Bloch theorem. Definition of quasi-momentum \(\mathbf{k}\) and difference from standard momentum: why and how is \(\mathbf{k}\) redundant? Reduced, extended and repeated schemes.
• Basic band properties. Time-reversal/Kramers symmetry, smoothness and group velocities, etc.
• Transfer matrix. Bands emerge from constructive/destructive interference of barrier scattering. Why use transfer matrix? How are matrices defined and used? Key use cases: resonant tunnelling and the \(N\to\infty\) limit.
• Nearly free electron. Why does a periodic potential open gaps in a folded parabola? What are the leading terms opening the individual gaps, and why?
• Counting the states. How to count the electronic states in a band? Density of states: definition, dimensionality fingerprints and van Hove singularities; why do we care about the density of states, at all?

Ashcroft, Neil W., and N. David Mermin. 1976. Solid State Physics. Harcourt College Publishers.

Grosso, Giuseppe, and Giuseppe Pastori Parravicini. 2014. Solid State Physics. 2nd ed. Academic Press.

---

1. Grosso and Pastori Parravicini (2014, sec. 1).

2. Ashcroft and Mermin (1976, chap. 9 p. 162).

3. Grosso and Pastori Parravicini (2014, sec. 1).

4. Grosso and Pastori Parravicini (2014, sec. 2.7).

5. Grosso and Pastori Parravicini (2014, sec. 2.6.3).