--- title: "Jellium" sort_title: "Jellium" date: 2021-11-23 categories: - Physics - Quantum mechanics - Perturbation layout: "concept" --- **Jellium**, also called the **uniform** or **homogeneous electron gas**, is a theoretical material where all electrons are free, and the ions' positive charge is smeared into a uniform background "jelly". This simple model lets us study electron interactions easily. ## Without interactions Let us start by neglecting electron-electron interactions. This is clearly a dubious assumption, but we will stick with it for now. For an infinitely large sample of jellium, the single-electron states are simply plane waves. We consider an arbitrary cube of volume $$V$$, and impose periodic boundary conditions on it, such that the single-particle orbitals are (suppressing spin): $$\begin{aligned} \Inprod{\vb{r}}{\psi_{\vb{k}}} = \psi_{\vb{k}}(\vb{r}) = \frac{1}{\sqrt{V}} \exp(i \vb{k} \cdot \vb{r}) \qquad \quad \vb{k} = \frac{2 \pi}{V^{1/3}} (n_x, n_y, n_z) \end{aligned}$$ Where $$n_x, n_y, n_z \in \mathbb{Z}$$. This is a discrete (but infinite) set of independent orbitals, so it is natural to use the [second quantization](/know/concept/second-quantization/) to write the non-interacting Hamiltonian $$\hat{H}_0$$, where $$\hbar^2 |\vb{k}|^2 / (2 m)$$ is the kinetic energy of the orbital with wavevector $$\vb{k}$$, and $$s$$ is the spin: $$\begin{aligned} \hat{H}_0 = \sum_{s} \sum_{\vb{k}} \frac{\hbar^2 |\vb{k}|^2}{2 m} \hat{c}_{s,\vb{k}}^\dagger \hat{c}_{s,\vb{k}} \end{aligned}$$ Assuming that the temperature $$T = 0$$, the $$N$$-electron ground state of this Hamiltonian is known as the **Fermi sea** or **Fermi sphere** $$\Ket{\mathrm{FS}}$$, and is constructed by filling up the single-electron states starting from the lowest energy: $$\begin{aligned} \Ket{\mathrm{FS}} = \prod_{s} \prod_{j = 1}^{N/2} \hat{c}_{s,\vb{k}_j}^\dagger \Ket{0} \end{aligned}$$ Because $$T = 0$$, all the electrons stay in their assigned state. The energy and wavenumber $$|\vb{k}|$$ of the highest filled orbital are called the **Fermi energy** $$\epsilon_F$$ and **Fermi wavenumber** $$k_F$$, and obey the expected kinetic energy relation: $$\begin{aligned} \boxed{ \epsilon_F = \frac{\hbar^2}{2 m} k_F^2 } \end{aligned}$$ The Fermi sea can be visualized in $$\vb{k}$$-space as a sphere with radius $$k_F$$. Because $$\vb{k}$$ is discrete, the sphere's surface is not smooth, but in the limit $$V \to \infty$$ it becomes perfect. Now, we would like a relation between the system's parameters, e.g. $$N$$ and $$V$$, and the resulting values of $$\epsilon_F$$ or $$k_F$$. The total population $$N$$ must be given by: $$\begin{aligned} N = \sum_{s} \sum_{\vb{k}} \matrixel{\mathrm{FS}}{\hat{c}_{s,\vb{k}}^\dagger \hat{c}_{s,\vb{k}}}{\mathrm{FS}} = \sum_{s} \frac{V}{(2 \pi)^3} \int_{-\infty}^\infty \matrixel{\mathrm{FS}}{\hat{c}_{s,\vb{k}}^\dagger \hat{c}_{s,\vb{k}}}{\mathrm{FS}} \dd{\vb{k}} \end{aligned}$$ Where we have turned the sum over $$\vb{k}$$ into an integral with a constant factor, by using that each orbital exclusively occupies a volume $$(2 \pi)^3 / V$$ in $$\vb{k}$$-space. At zero temperature, this inner product can only be $$0$$ or $$1$$, depending on whether $$\vb{k}$$ is outside or inside the Fermi sphere. We can therefore rewrite using a [Heaviside step function](/know/concept/heaviside-step-function/): $$\begin{aligned} N = \sum_{s} \frac{V}{(2 \pi)^3} \int_{-\infty}^\infty \Theta(k_F - |\vb{k}|) \dd{\vb{k}} = 2 \frac{V}{(2 \pi)^3} \int_{-\infty}^\infty \Theta(k_F - |\vb{k}|) \dd{\vb{k}} \end{aligned}$$ Where we realized that spin does not matter, and replaced the sum over $$s$$ by a factor $$2$$. In order to evaluate this 3D integral, we go to [spherical coordinates](/know/concept/spherical-coordinates/) $$(|\vb{k}|, \theta, \varphi)$$: $$\begin{aligned} N &= \frac{V}{4 \pi^3} \int_0^{2 \pi} \int_0^\pi \int_0^\infty \Theta(k_F - |\vb{k}|) |\vb{k}|^2 \sin(\theta) \dd{|\vb{k}|} \dd{\theta} \dd{\varphi} \\ &= \frac{V}{4 \pi^3} 4 \pi \int_0^{k_F} |\vb{k}|^2 \dd{|\vb{k}|} = \frac{V}{\pi^2} \bigg[ \frac{|\vb{k}|^3}{3} \bigg]_0^{k_F} = \frac{V}{3 \pi^2} k_F^3 \end{aligned}$$ Using that the electron density $$n = N/V$$, we thus arrive at the following relation: $$\begin{aligned} \boxed{ k_F^3 = 3 \pi^2 n } \end{aligned}$$ This result also justifies our assumption that $$T = 0$$: we can accurately calculate the density $$n$$ for many conducting materials, and this relation then gives $$k_F$$ and $$\epsilon_F$$. It turns out that $$\epsilon_F$$ is usually very large compared to the thermal energy $$k_B T$$ at reasonable temperatures, so we can conclude that thermal fluctuations are negligible. Now, $$\epsilon_F$$ is the highest single-electron energy, but about the total $$N$$-particle energy $$E^{(0)}$$? $$\begin{aligned} E^{(0)} = \matrixel{\mathrm{FS}}{\hat{H}_0}{\mathrm{FS}} = \sum_{s} \sum_{\vb{k}} \frac{\hbar^2 |\vb{k}|^2}{2 m} \matrixel{\mathrm{FS}}{\hat{c}_{s,\vb{k}}^\dagger \hat{c}_{s,\vb{k}}}{\mathrm{FS}} \end{aligned}$$ Once again, we turn the sum over $$\vb{k}$$ into an integral, and recognize the spin's irrelevance: $$\begin{aligned} E^{(0)} &= \sum_{s} \frac{V}{(2 \pi)^3} \int_{-\infty}^\infty \frac{\hbar^2 |\vb{k}|^2}{2 m} \matrixel{\mathrm{FS}}{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k}}}{\mathrm{FS}} \dd{\vb{k}} \\ &= \frac{\hbar^2 V}{8 \pi^3 m} \int_{-\infty}^\infty |\vb{k}|^2 \: \Theta(k_F - |\vb{k}|) \dd{\vb{k}} \end{aligned}$$ In spherical coordinates, we evaluate the integral and find that $$E^{(0)}$$ is proportional to $$k_F^5$$: $$\begin{aligned} E^{(0)} &= \frac{\hbar^2 V}{8 \pi^3 m} \int_0^{2 \pi} \int_0^\pi \int_0^\infty \Big( |\vb{k}|^2 \: \Theta(k_F - |\vb{k}|) \Big) |\vb{k}|^2 \sin(\theta) \dd{|\vb{k}|} \dd{\theta} \dd{\varphi} \\ &= \frac{\hbar^2 V}{8 \pi^3 m} 4 \pi \int_0^{k_F} |\vb{k}|^4 \dd{|\vb{k}|} = \frac{\hbar^2 V}{2 \pi^2 m} \bigg[ \frac{|\vb{k}|^5}{5} \bigg]_0^{k_F} = \frac{\hbar^2 V}{10 \pi^2 m} k_F^5 \end{aligned}$$ In general, it is more useful to consider the average kinetic energy per electron $$E^{(0)} / N$$, which we find to be as follows, using that $$k_F^3 = 3 \pi^2 n$$: $$\begin{aligned} \boxed{ \frac{E^{(0)}}{N} = \frac{3 \hbar^2}{10 m} k_F^2 = \frac{3}{5} \epsilon_F } \:\sim\: n^{2/3} \end{aligned}$$ Traditionally, this is expressed using a dimensionless parameter $$r_s$$, defined as the radius of a sphere containing a single electron, measured in Bohr radii $$a_0 \equiv 4 \pi \varepsilon_0 \hbar^2 / (e^2 m)$$: $$\begin{aligned} \frac{4 \pi}{3} (a_0 r_s)^3 = \frac{1}{n} = \frac{3 \pi^2}{k_F^3} \quad \implies \quad r_s = \Big( \frac{3}{4 \pi a_0^3 n} \Big)^{1/3} = \Big( \frac{9 \pi}{4} \Big)^{1/3} \frac{1}{a_0 k_F} \end{aligned}$$ Such that the ground state energy can be rewritten in Rydberg units of energy like so: $$\begin{aligned} \frac{E^{(0)}}{N} = \frac{3 \hbar^2}{10 m} \frac{4 \pi \varepsilon_0 e^2}{4 \pi \varepsilon_0 e^2} \frac{a_0^2 k_F^2}{a_0^2} = \frac{3 e^2}{40 \pi \varepsilon_0} \Big( \frac{9 \pi}{4} \Big)^{2/3} \frac{1}{a_0 r_s^2} \approx \frac{2.21}{r_s^2} \; \mathrm{Ry} \end{aligned}$$ ## With interactions To include Coulomb interactions, let us try [time-independent pertubation theory](/know/concept/time-independent-perturbation-theory/). Clearly, this will give better results when the interaction is relatively weak, if ever. The Coulomb potential is proportional to the inverse distance, and the average electron spacing is roughly $$n^{-1/3}$$, so the interaction energy $$E_\mathrm{int}$$ should scale as $$n^{1/3}$$. We already know that the kinetic energy $$E_\mathrm{kin} = E^{(0)}$$ scales as $$n^{2/3}$$, meaning perturbation theory should be reasonable if $$1 \gg E_\mathrm{int} / E_\mathrm{kin} \sim n^{-1/3}$$, so in the limit of high density $$n \to \infty$$. The two-body Coulomb interaction operator $$\hat{W}$$ is as follows in second-quantized form: $$\begin{aligned} \hat{W} = \frac{1}{2 V} \sum_{s_1 s_2} \sum_{\vb{k}_1 \vb{k}_2} \sum_{\vb{q} \neq 0} \frac{e^2}{\varepsilon_0 |\vb{q}|^2} \hat{c}_{s_1, \vb{k}_1 + \vb{q}}^\dagger \hat{c}_{s_2, \vb{k}_2 - \vb{q}}^\dagger \hat{c}_{s_2, \vb{k}_2} \hat{c}_{s_1, \vb{k}_1} \end{aligned}$$ The first-order correction $$E^{(1)}$$ to the ground state (i.e. Fermi sea) energy is then given by: $$\begin{aligned} E^{(1)} = \matrixel{\mathrm{FS}}{\hat{W}}{\mathrm{FS}} = \frac{e^2}{2 \varepsilon_0 V} \sum_{s_1 s_2} \sum_{\vb{k}_1 \vb{k}_2} \sum_{\vb{q} \neq 0} \frac{1}{|\vb{q}|^2} \matrixel{\mathrm{FS}}{ \hat{c}_{s_1, \vb{k}_1 + \vb{q}}^\dagger \hat{c}_{s_2, \vb{k}_2 - \vb{q}}^\dagger \hat{c}_{s_2, \vb{k}_2} \hat{c}_{s_1, \vb{k}_1} }{\mathrm{FS}} \end{aligned}$$ This inner product can only be nonzero if the two creation operators $$\hat{c}^\dagger$$ are for the same orbitals as the two annihilation operators $$\hat{c}$$. Since $$\vb{q} \neq 0$$, this means that $$s_1 = s_2$$, and that momentum is conserved: $$\vb{k}_2 = \vb{k}_1 \!+\! \vb{q}$$. And of course both $$\vb{k}_1$$ and $$\vb{k}_1 \!+\! \vb{q}$$ must be inside the Fermi sphere, to avoid annihilating an empty orbital. Let $$s = s_1$$ and $$\vb{k} = \vb{k}_1$$: $$\begin{aligned} E^{(1)} &= \frac{e^2}{2 \varepsilon_0 V} \sum_{s} \sum_{\vb{k}} \sum_{\vb{q} \neq 0} \frac{1}{|\vb{q}|^2} \matrixel{\mathrm{FS}}{ \hat{c}_{s, \vb{k} + \vb{q}}^\dagger \hat{c}_{s, \vb{k}}^\dagger \hat{c}_{s, \vb{k} + \vb{q}} \hat{c}_{s, \vb{k}} }{\mathrm{FS}} \\ &= \frac{- e^2}{2 \varepsilon_0 V} \sum_{s} \sum_{\vb{k}} \sum_{\vb{q} \neq 0} \frac{1}{|\vb{q}|^2} \matrixel{\mathrm{FS}}{ \big( \hat{c}_{s, \vb{k} + \vb{q}}^\dagger \hat{c}_{s, \vb{k} + \vb{q}}\big) \big(\hat{c}_{s, \vb{k}}^\dagger \hat{c}_{s, \vb{k}}\big) }{\mathrm{FS}} \\ &= \frac{- e^2}{2 \varepsilon_0 V} \sum_{s} \sum_{\vb{k}} \sum_{\vb{q} \neq 0} \frac{1}{|\vb{q}|^2} \Theta(k_F - |\vb{k}|) \:\Theta(k_F - |\vb{k} \!+\! \vb{q}|) \end{aligned}$$ Next, we convert the sum over $$\vb{q}$$ into an integral in spherical coordinates. Clearly, $$\vb{q}$$ is the "jump" made by an electron from one orbital to another, so the largest possible jump goes from a point on the Fermi surface to the opposite point, and thus has length $$2 k_F$$. This yields the integration limit, and therefore leads to: $$\begin{aligned} E^{(1)} &= \frac{- e^2}{(2 \pi)^3 \varepsilon_0} \sum_{\vb{k}} \int_0^{2 \pi} \!\!\int_0^\pi \!\!\int_0^\infty \Theta(k_F \!-\! |\vb{k}|) \: \Theta(k_F \!-\! |\vb{k} \!+\! \vb{q}|) \frac{|\vb{q}|^2}{|\vb{q}|^2} \sin(\theta_q) \dd{|\vb{q}|} \dd{\theta_q} \dd{\varphi_q} \\ &= \frac{- e^2}{2 \pi^2 \varepsilon_0} \sum_{\vb{k}} \int_0^{2 k_F} \Theta(k_F \!-\! |\vb{k}|) \: \Theta(k_F \!-\! |\vb{k} \!+\! \vb{q}|) \dd{|\vb{q}|} \end{aligned}$$ Where we have used that the direction of $$\vb{q}$$, i.e. $$(\theta_q,\varphi_q)$$, is irrelevant, as long as we define $$\theta_k$$ as the angle between $$\vb{q}$$ and $$\vb{k} \!+\! \vb{q}$$ when we go to spherical coordinates $$(|\vb{k}|, \theta_k, \varphi_k)$$ for $$\vb{k}$$: $$\begin{aligned} E^{(1)} &= \frac{- e^2 V}{16 \pi^5 \varepsilon_0} \int_0^{2 k_F} \!\!\!\!\int_0^{2 \pi} \!\!\!\int_0^\pi \!\!\!\int_0^\infty \!\Theta(k_F \!-\! |\vb{k}|) \: \Theta(k_F \!-\! |\vb{k} \!+\! \vb{q}|) \: |\vb{k}|^2 \sin(\theta_k) \dd{|\vb{k}|} \dd{\theta_k} \dd{\varphi_k} \dd{|\vb{q}|} \\ &= \frac{- e^2 V}{16 \pi^5 \varepsilon_0} \int_0^{2 k_F} \!\!\!\!\int_0^{2 \pi} \!\!\!\int_0^\pi \!\!\!\int_0^{k_F} \!\Theta(k_F \!-\! |\vb{k} \!+\! \vb{q}|) \: |\vb{k}|^2 \sin(\theta_k) \dd{|\vb{k}|} \dd{\theta_k} \dd{\varphi_k} \dd{|\vb{q}|} \end{aligned}$$ Unfortunately, this last step function is less easy to translate into integration limits. In effect, we are trying to calculate the intersection volume of two spheres, both with radius $$k_F$$, one centered on the origin (for $$\vb{k}$$), and the other centered on $$\vb{q}$$ (for $$\vb{k} \!+\! \vb{q}$$). Imagine a triangle with side lengths $$|\vb{k}|$$, $$|\vb{q}|$$ and $$|\vb{k} \!+\! \vb{q}|^2$$, where $$\theta_k$$ is the angle between $$|\vb{k}|$$ and $$|\vb{k} \!+\! \vb{q}|$$. The *law of cosines* then gives the following relation: $$\begin{aligned} |\vb{k}|^2 = |\vb{q}|^2 + |\vb{k} \!+\! \vb{q}|^2 - 2 |\vb{q}| |\vb{k} \!+\! \vb{q}| \cos(\theta_k) \end{aligned}$$ We already know that $$|\vb{k}| < k_F$$ and $$0 < |\vb{q}| < 2 k_F$$, so by isolating for $$\cos(\theta_k)$$, we can obtain bounds on $$\theta_k$$ and $$|\vb{k}|$$. Let $$|\vb{k}| \to k_F$$ in both cases, then: $$\begin{aligned} \cos(\theta_k) = \frac{|\vb{k} \!+\! \vb{q}|^2 + |\vb{q}|^2 - |\vb{k}|^2}{2 |\vb{k} \!+\! \vb{q}| |\vb{q}|} &\:\:\underset{|\vb{q}| \to 0}{>}\:\:\: \frac{k_F^2 + |\vb{q}|^2 - k_F^2}{2 k_F |\vb{q}|} = \frac{|\vb{q}|}{2 k_F} \\ &\underset{|\vb{q}| \to 2 k_F}{<}\:\: \frac{k_F^2 + 4 k_F ^2 - k_F^2}{2 k_F 2 k_F} = 1 \end{aligned}$$ Meaning that $$0 < \theta_k < \arccos{|\vb{q}| / (2 k_F)}$$. To get a lower limit for $$|\vb{k}|$$, we "cheat" by artificially demanding that $$\vb{k}$$ does not cross the halfway point between the spheres, with the result that $$|\vb{k}| \cos(\theta_k) > |\vb{q}|/2$$. Then, thanks to symmetry (both spheres have the same radius), we just multiply the integral by $$2$$, for $$\vb{k}$$ on the other side of the halfway point. Armed with these integration limits, we return to calculating $$E^{(1)}$$, substituting $$\xi \equiv \cos(\theta_k)$$: $$\begin{aligned} E^{(1)} &= \frac{- e^2 V}{16 \pi^5 \varepsilon_0} 2 \int_0^{2 k_F} \!\!\!\int_0^{2 \pi} \!\!\int_0^{\arccos{|\vb{q}| / (2 k_F)}} \!\!\int_{|\vb{q}|/(2 \cos{\theta_k})}^{k_F} |\vb{k}|^2 \sin(\theta_k) \dd{|\vb{k}|} \dd{\theta_k} \dd{\varphi_k} \dd{|\vb{q}|} \\ &= \frac{e^2 V}{8 \pi^5 \varepsilon_0} 2 \pi \int_0^{2 k_F} \!\!\!\int_1^{|\vb{q}| / (2 k_F)} \!\!\int_{|\vb{q}|/(2 \xi)}^{k_F} |\vb{k}|^2 \frac{\sin(\theta_k)}{\sin(\theta_k)} \dd{|\vb{k}|} \dd{\xi} \dd{|\vb{q}|} \\ &= \frac{- e^2 V}{4 \pi^4 \varepsilon_0} \int_0^{2 k_F} \!\!\!\int_{|\vb{q}| / (2 k_F)}^1 \!\!\int_{|\vb{q}|/(2 \xi)}^{k_F} |\vb{k}|^2 \dd{|\vb{k}|} \dd{\xi} \dd{|\vb{q}|} \end{aligned}$$ Where we have used that $$\varphi_k$$ does not appear in the integrand. Evaluating these integrals: $$\begin{aligned} E^{(1)} &= \frac{- e^2 V}{4 \pi^4 \varepsilon_0} \int_0^{2 k_F} \!\!\!\int_{|\vb{q}| / (2 k_F)}^1 \bigg[ \frac{|\vb{k}|^3}{3} \bigg]_{|\vb{q}|/(2 \xi)}^{k_F} \dd{\xi} \dd{|\vb{q}|} \\ &= \frac{- e^2 V}{4 \pi^4 \varepsilon_0} \int_0^{2 k_F} \!\!\!\int_{|\vb{q}| / (2 k_F)}^1 \bigg( \frac{k_F^3}{3} - \frac{|\vb{q}|^3}{24 \xi^3} \bigg) \dd{\xi} \dd{|\vb{q}|} \\ &= \frac{- e^2 V}{4 \pi^4 \varepsilon_0} \int_0^{2 k_F} \bigg[ \frac{k_F^3}{3} x + \frac{|\vb{q}|^3}{48 \xi^2} \bigg]_{|\vb{q}| / (2 k_F)}^1 \dd{|\vb{q}|} \\ &= \frac{- e^2 V}{4 \pi^4 \varepsilon_0} \int_0^{2 k_F} \bigg( \frac{k_F^3}{3} + \frac{|\vb{q}|^3}{48} - \frac{k_F^2 |\vb{q}|}{4} \bigg) \dd{|\vb{q}|} \\ &= \frac{- e^2 V}{4 \pi^4 \varepsilon_0} \bigg[ \frac{k_F^3 |\vb{q}|}{3} + \frac{|\vb{q}|^4}{192} - \frac{k_F^2 |\vb{q}|^2}{8} \bigg]_0^{2 k_F} \\ &= \frac{- e^2 V}{16 \pi^4 \varepsilon_0} k_F^4 = \frac{- e^2 N}{16 \pi^4 \varepsilon_0 n} k_F^4 = -\frac{3 e^2 N}{16 \pi^2 \varepsilon_0} k_F \end{aligned}$$ Per particle, the first-order energy correction $$E^{(1)}$$ is therefore found to be as follows: $$\begin{aligned} \boxed{ \frac{E^{(1)}}{N} = -\frac{3 e^2}{16 \pi^2 \varepsilon_0} k_F } \end{aligned}$$ This can also be written using the parameter $$r_s$$ introduced above, leading to: $$\begin{aligned} \frac{E^{(1)}}{N} = -\frac{3 e^2}{16 \pi^2 \varepsilon_0} \frac{a_0 k_F}{a_0} = -\frac{3 e^2}{16 \pi^2 \varepsilon_0} \Big( \frac{9 \pi}{4} \Big)^{1/3} \frac{1}{a_0 r_s} \end{aligned}$$ Consequently, for sufficiently high densities $$n$$, the total energy $$E$$ per particle is given by: $$\begin{aligned} \boxed{ \frac{E}{N} \approx \bigg( \frac{2.21}{r_s^2} - \frac{0.92}{r_s} \bigg) \; \mathrm{Ry} } \end{aligned}$$ Unfortunately, this is as far as we can go. In theory, the second-order energy correction $$E^{(2)}$$ is as shown below, but it turns out that it (and all higher orders) diverge: $$\begin{aligned} E^{(2)} = \sum_{\Psi_n \neq \mathrm{FS}} \frac{\big| \matrixel{\mathrm{FS}}{\hat{W}}{\Psi_n} \big|^2}{E^{(0)} - E_n} \end{aligned}$$ The only cure for this is to go to infinite order, where all the infinities add up to a finite result. ## References 1. H. Bruus, K. Flensberg, *Many-body quantum theory in condensed matter physics*, 2016, Oxford.