Introduction to the Jellium model

The jellium model consists of free electrons interacting with a uniform background positive charge and with themselves.
 
The Hamiltonian, in a plane-wave basis, is \[\hat{H}_{\text{Jellium}}=\sum_{k\sigma}\frac{\hbar^2k^2}{2m}c_{k\sigma}^\dagger c_{k\sigma}^{\phantom{\dagger}}+\frac{4\pi e^2}{2V}\sum_{kk'\sigma\sigma', \ q\ne 0}\frac{1}{q^2}c_{k+q\sigma}^\dagger c_{k'-q\sigma'}^\dagger c_{k'\sigma'}^{\phantom{\dagger}}c_{k\sigma}^{\phantom{\dagger}}\] \[V(q)=\begin{cases}\frac{4\pi e^2}{2Vq^2} & q\ne0 \\ 0 & q=0.\end{cases}\] The potential vanishes when \(q=0\) due to the cancellation by the uniform positive background.

Variational principle

Study first with the variational principle. Consider \[E_{\text{trial}}=\frac{\left\langle\psi\middle|\hat{H}\middle|\psi\right\rangle}{\left\langle\psi\middle|\psi\right\rangle}\] We always have \(E_{\text{trial}}\ge E_{gs}\).
 
Proof: Start with \(\left|\psi\right\rangle=\sum_nc_n\left|n\right\rangle\) where \(\hat{H}\left|n\right\rangle=E_n\left|n\right\rangle\) \[E_{\text{trial}}=\frac{\sum_n|c_n|^2E_n}{\sum_n|c_n|^2}\ge\frac{\sum_n|c_n|^2E_0}{\sum_n|c_n|^2}\] For our trial state, we will look at a single Slater determinant. Define \[\alpha_{p\sigma}^\dagger=\sum_la_l^pc_{l\sigma}^\dagger\] and let \[\left|\psi\right\rangle=\prod_{p,\sigma}\alpha_{p\sigma}^\dagger\left|0\right\rangle\] Note that \[\{\alpha_{p\sigma}^\dagger,\alpha_{p'\sigma'}^{\phantom{\dagger}}\}=\sum_{ll'}a_l^pa_{l'}^{p'*}\{c_{l\sigma}^\dagger,c_{l'\sigma'}^{\phantom{\dagger}}\}=\sum_la_l^pa_{l}^{p'*}\delta_{\sigma\sigma'}\] If we choose the \(a_l^p\) vectors to be an orthonormal set, then \[\sum_la_l^pa_l^{p'*}=\delta_{pp'}\] so \(\{\alpha_{p\sigma}^\dagger,\alpha_{p'\sigma'}^{\phantom{\dagger}}\}=\delta_{pp'}\). Next, define a set of vectors \(b_l^p\) such that \[\sum_pb_l^pa_{l'}^p=\delta_{ll'}\] Then, \[\sum_pb_l^p\alpha_{p\sigma}^\dagger=\sum_p\sum_{l'}b_l^pa_{l'}^pc_{l'\sigma}^\dagger=c_{l\sigma}^\dagger.\] So, we can now compute averages. \[\begin{aligned} \left\langle\psi\middle|c_{l\sigma}^\dagger c_{l\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle=& \ \sum_{pp'}b_l^pb_l^{p'*}\left\langle\psi\middle|\alpha_{p\sigma}^\dagger\alpha_{p'\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle \\ =& \ \sum_{pp'}b_l^pb_l^{p'*}\delta_{pp'}\delta(p,\sigma\in\psi) \\ =& \ \sum_p|b_l^p|^2\delta(p,\sigma\in\psi) \end{aligned}\] and \[\begin{aligned} \left\langle\psi\middle|c_{l_1\sigma}^\dagger c_{l_2\sigma'}^\dagger c_{l_3\sigma'}^{\phantom{\dagger}}c_{l_4\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle=& \ \sum_{p_1p_2p_3p_4}b_{l_1}^{p_1}b_{l_2}^{p_2}b_{l_3}^{p_3*}b_{l_4}^{p_4*}\left\langle\psi\middle|\alpha_{p_1\sigma}^\dagger\alpha_{p_2\sigma'}^\dagger\alpha_{p_3\sigma'}^{\phantom{\dagger}}\alpha_{p_4\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle \\ =& \ \sum_{p_1p_2p_3p_4}b_{l_1}^{p_1}b_{l_2}^{p_2}b_{l_3}^{p_3*}b_{l_4}^{p_4*}\left[\delta_{p_1p_4}\delta_{p_2p_3}\left\langle\psi\middle|\alpha_{p_1\sigma}^\dag\alpha_{p_2\sigma'}^\dag\alpha_{p_2\sigma'}^{\phantom{\dagger}}\alpha_{p_1\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle \right.\\ &\left. +\delta_{p_1p_3}\delta_{p_2p_4}\delta_{\sigma\sigma}\left\langle\psi\middle|\alpha_{p_1\sigma}^\dag\alpha_{p_2\sigma}^\dag\alpha_{p_2\sigma}^{\phantom{\dagger}}\alpha_{p_1\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle\right] \end{aligned}\] Note when \(p_1=p_2\) and \(\sigma=\sigma'\) we get zero in both cases since \(\alpha_{p\sigma}^2=0\). So, we can assume \(p_1\ne p_2\). \[\begin{gathered} =\sum_{p_1p_2p_3p_4}b_{l_1}^{p_1}b_{l_2}^{p_2}b_{l_3}^{p_3*}b_{l_4}^{p_4*}\left[\delta_{p_1p_4}\delta_{p_2p_3}\left\langle\psi\middle|\alpha_{p_1\sigma}^\dag\alpha_{p_1\sigma}^{\phantom{\dagger}}\alpha_{p_2\sigma'}^\dag\alpha_{p_2\sigma'}^{\phantom{\dagger}}\middle|\psi\right\rangle \right. \\ \left. -\delta_{p_1p_3}\delta_{p_2p_4}\delta_{\sigma\sigma}\left\langle\psi\middle|\alpha_{p_1\sigma}^\dag\alpha_{p_1\sigma}^{\phantom{\dagger}}\alpha_{p_2\sigma}^\dag\alpha_{p_2\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle\right] \end{gathered}\] Since in this expression, the term with \(p_1=p_2=p_3=p_4\) and \(\sigma=\sigma'\) cancels from both terms, we do not need to worry about the restriction \(p_1\ne p_2\) anymore since the total terms vanish when \(p_1=p_2\). \[=\sum_{p_1}b_{l_1}^{p_1}b_{l_4}^{p_1*}\delta(p_1\sigma\in\psi)\sum_{p_2}b_{l_2}^{p_2}b_{l_3}^{p_2*}\delta(p_2{\sigma'}\in\psi)-\sum_{p_1}b_{l_1}^{p_1}b_{l_3}^{p_1*}\delta(p_1\sigma\in\psi)\sum_{p_2}b_{l_2}^{p_2}b_{l_2}^{p_2*}\delta(p_2{\sigma'}\in\psi)\] Note that we have just shown \[\left\langle\psi\middle|c_{l_1\sigma}^\dag c_{l_2\sigma'}^\dag c_{l_3\sigma'}^{\phantom{\dagger}}c_{l_4\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle=\left\langle\psi\middle|c_{l_1\sigma}^\dag c_{l_4\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle\left\langle\psi\middle|c_{l_2\sigma'}^\dag c_{l_3\sigma'}^{\phantom{\dagger}}\middle|\psi\right\rangle-\left\langle\psi\middle|c_{l_1\sigma}^\dag c_{l_3\sigma'}^{\phantom{\dagger}}\middle|\psi\right\rangle\left\langle\psi\middle|c_{l_2\sigma'}^\dag c_{l_4\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle\] which is called Wick’s theorem – replace averages of products of operators by products of averages of pairs of operators. Holds only for noninteracting particles/Slater determinants.
 
Now we choose the \(a_k^p\) values. Suppose we take \[\left|\psi_{trial}\right\rangle=\prod_{k,\sigma}c_{k\sigma}^\dag\left|0\right\rangle\] \[\hat{H}=\sum_{k\sigma}\frac{\hbar^2k^2}{2m}c_{k\sigma}^\dag c_{k\sigma}^{\phantom{\dagger}}+\frac{1}{2}\sum_{kk'\sigma\sigma',q\ne 0}\frac{4\pi e^2}{Vq^2}c_{k+q\sigma}^\dag c_{k'-q\sigma'}^\dag c_{k'\sigma'}^{\phantom{\dagger}}c_{k\sigma}^{\phantom{\dagger}}\] so \[\left\langle\psi\middle|c_{k+q\sigma}^\dag c_{k'-q\sigma'}^\dag c_{k'\sigma'}^{\phantom{\dagger}}c_{k\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle=\left\langle\psi\middle|c_{k+q\sigma}^\dag c_{k\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle\left\langle\psi\middle|c_{k'-q\sigma'}^\dag c_{k'\sigma'}^{\phantom{\dagger}}\middle|\psi\right\rangle-\left\langle\psi\middle|c_{k+q\sigma}^\dag c_{k'\sigma'}^{\phantom{\dagger}}\middle|\psi\right\rangle\left\langle\psi\middle|c_{k'-q\sigma'}^\dag c_{k\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle\] The first term vanishes unless \(q=0\). But when \(q=0\), \(V(q)=0\) so we neglect those terms for jellium. The second requires \(k'=k+q\), \(\sigma=\sigma'\), or \(q=k'-k\). So we find \[\left\langle\psi\middle|H\middle|\psi\right\rangle\sum_{k\sigma}\frac{\hbar^2k^2}{2m}\left\langle\psi\middle|c_{k\sigma}^\dag c_{k\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle-\sum_{k\ne k',\sigma}\frac{4\pi e^2}{|k'-k|^2V}\left\langle\psi\middle|c_{k'\sigma}^\dag c_{k'\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle\left\langle\psi\middle|c_{k\sigma}^\dag c_{k\sigma}^{\phantom{\dagger}}\middle|\psi\right\rangle\] Since \(\frac{\hbar^2k^2}{2m}\) is an increasing function of \(k\), we minimize the kinetic energy by filling the lowest \(k\) levels first. This is called the “bathtub principle”. We use the state with minimum kinetic energy to estimate the ground state energy of jellium.

Estimating the ground state

If we have \(N\) particles the density is \(N/V\). Then we get \[\sum_{k<k_F,\sigma}1=N\implies 2\frac{V}{(2\pi)^3}\int_0^{k_F}4\pi k^2 \ dk=N\] So, \[n=\frac{N}{V}=\frac{2}{2\pi^2}\int_0^{k_F}k^2dk=\frac{k_F^3}{3\pi^2}\implies k_F=(3\pi^2n)^{1/3}\] The kinetic energy becomes \[2\sum_{k<k_F}\frac{\hbar^2k^2}{2m}=\frac{2V}{(2\pi)^3}4\pi\frac{\hbar^2}{2m}\int_0^{k_F}k^4 \ dk=\frac{\hbar^2}{2m}\frac{k_F^5}{5\pi^2}V=V\frac{(3\pi^2n)^{5/3}}{5\pi^2}\frac{\hbar^2}{2m}\] The potential energy term becomes \[-\frac{1}{2}2\frac{V^2}{(2\pi)^6}\int_{k<k_F}d^3k\int_{k'<k_F}\frac{4\pi e^2}{V|k'-k|^2}\] Do the \(k'\) integration first. Choose \(z\)-axis along \(\mathbf{k}\) direction. \[|\mathbf{k}'-\mathbf{k}|^2=k'^2-\mathbf{k}'\cdot\mathbf{k}+k^2=k'^2-2k'k\cos{\theta'}+k^2\] \[\begin{aligned} \int_0^{k_F}k'^2dk'\int_{-1}^1d\cos{\theta'}\int_0^{2\pi}d\phi'\frac{1}{k^2-2kk'\cos{\theta'}+k'^2}=& \ 2\pi\int_0^{k_F}k'^2dk'\left(-\frac{1}{2kk'}\right)\ln(k^2-2kk'\cos{\theta}+k'^2)\biggr\rvert_{-1}^1\\=& \ -\frac{\pi}{k}\int_0^{k_F}k'dk'\ln\left(\frac{k^2-2kk'+k'^2}{k^2+2kk'+k'^2}\right) \\ =& \ \frac{2\pi}{k}\int_0^{k_F}dk'k'\ln\left \vert\frac{k+k'}{k-k'}\right\vert \\ =& \ 2\pi\left(\frac{k_F^2-k^2}{2k}\ln\left\vert\frac{k_F+k}{k_F-k}\right\vert+k_F\right) \end{aligned}\] Now do the \(k\) integral. \[=\frac{V^2}{(2\pi)^6}\frac{4\pi e^2}{V}\cdot2\pi\cdot4\pi\int_0^{k_F}dk\left(\frac{1}{2}(kk_F^2-k^3)\ln\left\vert\frac{k_F+k}{k_F-k}\right\vert+k^2k_F\right)=\frac{Ve^2}{2\pi^3}\frac{k_F^4}{2}.\] So \[\frac{E_{trial}}{V}=\frac{\hbar^2}{2m}\frac{k_F^5}{5}-\frac{e^2}{4\pi^3}k_F^4\] Recall \(n=\frac{k_F^3}{3\pi^2}\). Express energy in Rydbergs\(=e^2/2a_0\). \[E_N=\frac{E}{Vk_F^3/3\pi^2}=\frac{e^2}{2a_0}\left(\frac{3\hbar^2k_F^2a_0}{5me^2}-\frac{3k_Fa_0}{2\pi}\right).\] Define \(r_s=r/a_0\) where \(r\) is the radius of a sphere with a single electron. \[4\pi r^3=\frac{1}{n}, \ r=\left(\frac{3}{4\pi n}\right)^{1/3}\] \[r_s=\left(\frac{3}{4\pi n}\right)^{1/3}\frac{1}{a_0}=\left(\frac{9\pi}{4}\right)^{1/3}\frac{1}{k_Fa_0}\] \[\frac{E}{N}=\frac{e^2}{2a_0}\left(\frac{3}{5}(k_Fa_0)^2-\frac{3}{2\pi}(k_Fa_0)\right)=\left(\frac{2.210}{r_S^2}-\frac{0.916}{r_S}\right)\text{ Ry}.\]