---

Categories: Physics, Quantum mechanics.

Lindhard function

The Lindhard function describes how an external perturbation affects a gas of free electrons, and can be regarded as the backbone of quantum mechanical models of metals.

Derivation

Consider the following general harmonic perturbation $\hat{H}_1$, where $U(\vb{r})$ is an arbitrary potential, $\hat{n}(\vb{r})$ is the electron density operator, and $\eta \to 0^{+}$ is needed to ensure convergence:

$$\begin{aligned} \hat{H}_{1,S}(t) = e^{i (\omega + i \eta) t} \int_{-\infty}^\infty U(\vb{r}) \: \hat{n}_S(\vb{r}) \dd{\vb{r}} \end{aligned}$$

Here, the subscript $S$ refers to the Schrödinger picture. The first-order effect of this perturbation on the expected density $\Expval{\hat{n}}$ is then given by the Kubo formula below, where the subscript $I$ refers to the interaction picture, and the expectation $\Expval{}_0$ is taken with respect to the thermal equilibrium that existed before the perturbation was applied:

$$\begin{aligned} \delta\!\Expval{\hat{n}(\vb{r}, t)} &= -\frac{i}{\hbar} \int_{-\infty}^\infty \Theta(t - t') \Expval{\Comm{\hat{n}_I(\vb{r}, t)}{\hat{H}_{1,I}(t')}}_0 \dd{t'} \\ &= \iint_{-\infty}^\infty \chi(\vb{r}, \vb{r}'; t, t') \: U(\vb{r}') \: e^{i (\omega + i \eta) t'} \dd{t'} \dd{\vb{r}'} \end{aligned}$$

Where $\Theta$ is a Heaviside step function, and the response $\chi$ is the density-density correlation function (i.e. a two-particle Green’s function):

$$\begin{aligned} \chi(\vb{r}, \vb{r}'; t, t') \equiv - \frac{i}{\hbar} \Theta(t - t') \Expval{\Comm{\hat{n}_I(\vb{r}, t)}{\hat{n}_I(\vb{r}', t')}}_0 \end{aligned}$$

Let us assume that the unperturbed system (i.e. without $U$) was spatially uniform, so that $\chi$ only depends on the difference $\vb{r} - \vb{r}'$. We then take its Fourier transform $\vb{r} - \vb{r}' \to \vb{q}$:

$$\begin{aligned} \chi(\vb{q}; t, t') &= \int_{-\infty}^\infty \chi(\vb{r} - \vb{r}'; t, t') \: e^{- i \vb{q} \cdot (\vb{r} - \vb{r}')} \dd{\vb{r}} \\ &= -\frac{i}{\hbar} \frac{\Theta(t \!-\! t')}{(2 \pi)^{2 D}} \iiint_{-\infty}^\infty \Expval{\Comm{\hat{n}_I(\vb{q}_1, t)}{\hat{n}_I(\vb{q}_2, t')}}_0 \: e^{i \vb{q}_1 \cdot \vb{r}} e^{i \vb{q}_2 \cdot \vb{r}'} e^{- i \vb{q} \cdot (\vb{r} - \vb{r}')} \dd{\vb{q}_1} \dd{\vb{q}_2} \dd{\vb{r}} \end{aligned}$$

Where both $\hat{n}_I$ have been written as inverse Fourier transforms, giving a factor $(2 \pi)^{-D}$ each, where $D$ is the number of spatial dimenions. We rearrange to get a Dirac delta function $\delta$:

$$\begin{aligned} \chi(\vb{q}; t, t') &= -\frac{i}{\hbar} \frac{\Theta(t \!-\! t')}{(2 \pi)^{2 D}} \iiint \Expval{\Comm{\hat{n}_I(\vb{q}_1, t)}{\hat{n}_I(\vb{q}_2, t')}}_0 \: e^{i (\vb{q}_1 - \vb{q}) \cdot \vb{r}} e^{i (\vb{q}_2 + \vb{q}) \cdot \vb{r}'} \dd{\vb{q}_1} \dd{\vb{q}_2} \dd{\vb{r}} \\ &= -\frac{i}{\hbar} \frac{\Theta(t \!-\! t')}{(2 \pi)^D} \iint \Expval{\Comm{\hat{n}_I(\vb{q}_1, t)}{\hat{n}_I(\vb{q}_2, t')}}_0 \: \delta(\vb{q}_1 \!-\! \vb{q}) \: e^{i (\vb{q}_2 + \vb{q}) \cdot \vb{r}'} \dd{\vb{q}_1} \dd{\vb{q}_2} \\ &= -\frac{i}{\hbar} \frac{\Theta(t \!-\! t')}{(2 \pi)^D} \int \Expval{\Comm{\hat{n}_I(\vb{q}, t)}{\hat{n}_I(\vb{q}_2, t')}}_0 \: e^{i (\vb{q}_2 + \vb{q}) \cdot \vb{r}'} \dd{\vb{q}_2} \end{aligned}$$

On the left, $\vb{r}'$ does not appear, so it must also disappear on the right. To handle this properly, we consider an arbitrary finite piece of volume $V$ of this infinite system. Clearly $\int_V \dd{\vb{r}'} = V$, so:

$$\begin{aligned} \chi(\vb{q}; t, t') &= -\frac{i}{\hbar} \frac{\Theta(t \!-\! t')}{(2 \pi)^D} \frac{1}{V} \int_V \int_{-\infty}^\infty \Expval{\Comm{\hat{n}_I(\vb{q}, t)}{\hat{n}_I(\vb{q}_2, t')}}_0 \: e^{i (\vb{q}_2 + \vb{q}) \cdot \vb{r}'} \dd{\vb{q}_2} \dd{\vb{r}'} \end{aligned}$$

This gives a Dirac delta function when we recover the full infinite system by letting $V \to \infty$ (we will take this limit properly later, but beware that some authors set $V = 1$ until then):

$$\begin{aligned} \chi(\vb{q}; t, t') &= -\frac{i}{\hbar V} \Theta(t - t') \int_{-\infty}^\infty \Expval{\Comm{\hat{n}_I(\vb{q}, t)}{\hat{n}_I(\vb{q}_2, t')}}_0 \: \delta(\vb{q}_2 \!+\! \vb{q}) \dd{\vb{q}_2} \\ &= -\frac{i}{\hbar V} \Theta(t - t') \Expval{\Comm{\hat{n}_I(\vb{q}, t)}{\hat{n}_I(-\vb{q}, t')}}_0 \end{aligned}$$

Similarly, if the unperturbed Hamiltonian $\hat{H}_0$ is time-independent, then $\chi$ only depends on the time difference $t - t'$. Note that $\delta{\Expval{\hat{n}}}$ already has the form of a Fourier transform $t - t' \to \omega + i \eta$, which gives us an opportunity to rewrite $\chi$ in the Lehmann representation:

$$\begin{aligned} \chi(\vb{q}, \omega) = \frac{1}{Z V} \sum_{\nu \nu'} \frac{\matrixel{\nu}{\hat{n}_S(\vb{q})}{\nu'} \matrixel{\nu'}{\hat{n}_S(-\vb{q})}{\nu}}{\hbar (\omega + i \eta) + E_\nu - E_{\nu'}} \Big( e^{-\beta E_\nu} - e^{- \beta E_{\nu'}} \Big) \end{aligned}$$

Where $\Ket{\nu}$ and $\Ket{\nu'}$ are many-particle eigenstates of $\hat{H}_0$, and $Z$ is the grand partition function. To get ready for the calculations ahead, we swap $\nu$ and $\nu'$ in the second term, so the response function is as shown below. All operators are in the Schrödinger picture from now on:

$$\begin{aligned} \chi(\vb{q}, \omega) = \frac{1}{Z V} \sum_{\nu \nu'} \bigg( \frac{\matrixel{\nu}{\hat{n}(\vb{q})}{\nu'} \matrixel{\nu'}{\hat{n}(-\vb{q})}{\nu}} {\hbar (\omega + i \eta) + E_\nu - E_{\nu'}} - \frac{\matrixel{\nu}{\hat{n}(-\vb{q})}{\nu'} \matrixel{\nu'}{\hat{n}(\vb{q})}{\nu}} {\hbar (\omega + i \eta) + E_{\nu'} - E_\nu} \bigg) e^{-\beta E_\nu} \end{aligned}$$

Next, we need to rewrite $\hat{n}(\vb{q})$ somehow. To do this, we must neglect electron-electron interactions, leaving us with a simple Fermi gas. In that model, $V$ is a cube with periodic boundary conditions, meaning that $\vb{q}$ has a discrete spectrum of allowed values, and the single-electron states are just plane waves. This lets us write, in the second quantization:

$$\begin{aligned} \hat{n}(\vb{q}) = \sum_{\sigma \vb{k}} \hat{c}_{\sigma,\vb{k}}^\dagger \hat{c}_{\sigma,\vb{k} + \vb{q}} \qquad \qquad \hat{n}(-\vb{q}) = \hat{n}^\dagger(\vb{q}) \end{aligned}$$

When neglecting interactions, it is tradition to rename $\chi$ to $\chi_0$. We insert $\hat{n}$, suppressing spin:

$$\begin{aligned} \chi_0(\vb{q}, \omega) &= \frac{1}{Z V} \sum_{\vb{k} \vb{k}'} \sum_{\nu \nu'} \bigg( \frac{\matrixel{\nu}{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}}{\nu'} \matrixel{\nu'}{\hat{c}_{\vb{k}' + \vb{q}}^\dagger \hat{c}_{\vb{k}'}}{\nu}} {\hbar (\omega + i \eta) + E_\nu - E_{\nu'}} - \frac{\matrixel{\nu}{\hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}}}{\nu'} \matrixel{\nu'}{\hat{c}_{\vb{k}'}^\dagger \hat{c}_{\vb{k}' + \vb{q}}}{\nu}} {\hbar (\omega + i \eta) + E_{\nu'} - E_\nu} \bigg) e^{-\beta E_\nu} \end{aligned}$$

Here, $\matrixel{\nu}{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}}{\nu'}$ is only nonzero if $\Ket{\nu'}$ is contructed from $\Ket{\nu}$ by moving an electron from $\vb{k}$ to $\vb{k} \!+\! \vb{q}$, and analogously for the other inner products. As a result, $\vb{k} = \vb{k}'$ (and $\sigma = \sigma'$).

For the same reason, the energy difference $E_\nu - E_{\nu'}$ can simply be replaced by the cost of the single-particle excitation $\xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}$, where $\xi_{\vb{k}}$ is the energy of a $\vb{k}$-orbital. Therefore:

$$\begin{aligned} \chi_0(\vb{q}, \omega) &= \frac{1}{Z V} \sum_{\vb{k}} \sum_{\nu \nu'} \bigg( \frac{\matrixel{\nu}{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}}{\nu'} \matrixel{\nu'}{\hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}}}{\nu}} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} - \frac{\matrixel{\nu}{\hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}}}{\nu'} \matrixel{\nu'}{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}}{\nu}} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} \bigg) e^{-\beta E_\nu} \end{aligned}$$

Notice that we have eliminated all dependence on $\Ket{\nu'}$, so we remove it by $\sum_{\nu'} \Ket{\nu'} \Bra{\nu'} = 1$:

$$\begin{aligned} \chi_0(\vb{q}, \omega) &= \frac{1}{Z V} \sum_{\vb{k}} \sum_{\nu} \bigg( \frac{\matrixel{\nu}{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}} \hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}}}{\nu}} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} - \frac{\matrixel{\nu}{\hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}} \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}}{\nu}} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} \bigg) e^{-\beta E_\nu} \\ &= \frac{1}{Z V} \sum_{\vb{k}} \sum_{\nu} \frac{\matrixel{\nu}{\comm{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}} {\hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}}} \: e^{- \beta \hat{H}_0}}{\nu}} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} \end{aligned}$$

Where we recognized the commutator, and eliminated $E_\nu$ using $\hat{H}_0 \Ket{n} = E_\nu \Ket{\nu}$. The resulting expression has the form of a matrix trace $\Tr$, and therefore a thermal expectation $\Expval{}_0$:

$$\begin{aligned} \chi_0(\vb{q}, \omega) &= \frac{1}{Z V} \sum_{\vb{k}} \frac{\Tr\!\big(\comm{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}} {\hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}}} \: e^{- \beta \hat{H}_0} \big)} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} \\ &= \frac{1}{V} \sum_{\vb{k}} \frac{\expval{\comm{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}}{\hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}}}}_0} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} \end{aligned}$$

The commutator can be evaluated, and in this particular case it turns out to be:

$$\begin{aligned} \comm{\hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}}{\hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k}}} = \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k}} - \hat{c}_{\vb{k} + \vb{q}}^\dagger \hat{c}_{\vb{k} + \vb{q}} \end{aligned}$$

We substitute this result into $\chi_0$, and reintroduce the spin index $\sigma$ associated with $\vb{k}$:

$$\begin{aligned} \chi_0(\vb{q}, \omega) = \frac{1}{V} \sum_{\sigma \vb{k}} \frac{\expval{\hat{c}_{\sigma,\vb{k}}^\dagger \hat{c}_{\sigma,\vb{k}} - \hat{c}_{\sigma,\vb{k}+\vb{q}}^\dagger \hat{c}_{\sigma,\vb{k}+\vb{q}}}_0} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} \end{aligned}$$

The operator $\hat{c}_{\sigma.\vb{k}}^\dagger \hat{c}_{\sigma.\vb{k}}$ simply counts the number of electrons in state $(\sigma, \vb{k})$, which is given by the Fermi-Dirac distribution $n_F$. This gives us the Lindhard response function:

$$\begin{aligned} \boxed{ \chi_0(\vb{q}, \omega) = \frac{1}{V} \sum_{\sigma \vb{k}} \frac{n_F(\xi_{\vb{k}}) - n_F(\xi_{\vb{k} + \vb{q}})} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} } \end{aligned}$$

This is its most general form, but we can go further by exploiting the Fermi gas’ periodic boundary conditions, which enable us to turn the sum into an integral. Furthermore, electrons are spin-1/2 particles, so $\sum_{\sigma}$ simply gives a constant factor $2$:

$$\begin{aligned} \boxed{ \chi_0(\vb{q}, \omega) = \frac{2}{(2 \pi)^D} \int_{-\infty}^{\infty} \frac{n_F(\xi_{\vb{k}}) - n_F(\xi_{\vb{k} + \vb{q}})} {\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} \dd{\vb{k}} } \end{aligned}$$

Beware that some authors define $\chi_0$ as the response of the charge density to an electric potential, leading to an extra factor $e_0^2$, where $e_0$ is the electron charge. In contrast, our result describes the particle density response due to an energy potential.

Evaluation

Let us consider $D = 3$ (calculations in 1D or 2D are also possible), and rewrite the integrand by substituting $\vb{k}' \equiv - \vb{k} - \vb{q}$ in the second term and using $\xi_{-\vb{k}} = \xi_{\vb{k}}$:

$$\begin{aligned} \chi_0(\vb{q}, \omega) &= \frac{2}{(2 \pi)^3} \int_{-\infty}^{\infty} \bigg( \frac{n_F(\xi_{\vb{k}})}{\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} - \frac{n_F(\xi_{\vb{k} + \vb{q}})}{\hbar (\omega + i \eta) + \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}}} \bigg) \dd{\vb{k}} \\ &= \frac{1}{4 \pi^3} \bigg( \int_{-\infty}^{\infty} \frac{n_F(\xi_{\vb{k}})}{\xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}} + \hbar (\omega + i \eta)} \dd{\vb{k}} - \int_{-\infty}^{\infty} \frac{n_F(\xi_{- \vb{k}'})}{\xi_{- \vb{k}' - \vb{q}} - \xi_{- \vb{k}'} + \hbar (\omega + i \eta)} \dd{\vb{k}'} \bigg) \\ &= \frac{1}{4 \pi^3} \int_{-\infty}^{\infty} \bigg( \frac{1}{\xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}} + \hbar (\omega + i \eta)} + \frac{1}{\xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}} - \hbar (\omega + i \eta)} \bigg) n_F(\xi_{\vb{k}}) \dd{\vb{k}} \end{aligned}$$

Considering the system’s symmetry, it is natural to switch to spherical coordinates $(|\vb{k}|, \theta, \varphi)$:

$$\begin{aligned} \chi_0 &= \frac{1}{4 \pi^3} \int_0^\infty \!\! \int_0^\pi \! \int_{-\pi}^\pi \! \bigg( \frac{1}{\xi_{\vb{k}} \!-\! \xi_{\vb{k} + \vb{q}} + \hbar (\omega \!+\! i \eta)} + \frac{1}{\xi_{\vb{k}} \!-\! \xi_{\vb{k} + \vb{q}} - \hbar (\omega \!+\! i \eta)} \bigg) n_F(\xi_{\vb{k}}) \, |\vb{k}|^2 \sin{\theta} \dd{\varphi} \dd{\theta} \dd{|\vb{k}|} \\ &= \frac{2 \pi}{4 \pi^3} \int_0^\infty \!\! \int_0^\pi \bigg( \frac{1}{\xi_{\vb{k}} \!-\! \xi_{\vb{k} + \vb{q}} + \hbar (\omega \!+\! i \eta)} + \frac{1}{\xi_{\vb{k}} \!-\! \xi_{\vb{k} + \vb{q}} - \hbar (\omega \!+\! i \eta)} \bigg) n_F(\xi_{\vb{k}}) \, |\vb{k}|^2 \sin{\theta} \dd{\theta} \dd{|\vb{k}|} \end{aligned}$$

At absolute zero $T = 0$, the Fermi gas has a well-known ground state, in which the occupied single-particle plane wave states form a sphere in $\vb{k}$-space with radius $k_F$. This is a good approximation for nonzero $T$ too, because in practice, the electrons’ Fermi level $\xi_F \equiv \hbar^2 k_F^2 / (2 m)$ is huge compared to thermal excitations. Formally, $n_F(\xi_\vb{k}) \approx \Theta(k_F - |\vb{k}|)$, so we adjust the integration limit:

$$\begin{aligned} \chi_0 &\approx \frac{1}{2 \pi^2} \int_0^{k_F} \! \int_0^\pi \bigg( \frac{1}{\xi_{\vb{k}} \!-\! \xi_{\vb{k} + \vb{q}} + \hbar (\omega \!+\! i \eta)} + \frac{1}{\xi_{\vb{k}} \!-\! \xi_{\vb{k} + \vb{q}} - \hbar (\omega \!+\! i \eta)} \bigg) |\vb{k}|^2 \sin{\theta} \dd{\theta} \dd{|\vb{k}|} \end{aligned}$$

Next, let us take a closer look at the excitation energy in the denominator. Thanks to symmetry, we can always choose our coordinate system such that $\vb{q}$ points in the $\theta = 0$ direction, yielding:

$$\begin{aligned} \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}} &= \frac{\hbar^2}{2 m} \Big( (\vb{k} \cdot \vb{k}) - (\vb{k} \!+\! \vb{q}) \cdot (\vb{k} \!+\! \vb{q}) \Big) \\ &= - \frac{\hbar^2}{2 m} \Big( |\vb{q}|^2 + 2 |\vb{k}| |\vb{q}| \cos{\theta} \Big) \end{aligned}$$

By defining the dimensionless variables $x \equiv |\vb{q}| / (2 k_F)$ and $\kappa \equiv |\vb{k}| / k_F$, and recognizing the Fermi energy $\xi_F = \hbar^2 k_F^2 / (2 m)$, this turns into:

$$\begin{aligned} \xi_{\vb{k}} - \xi_{\vb{k} + \vb{q}} &= - 4 \xi_F (x^2 + x \kappa \cos{\theta}) \end{aligned}$$

We insert this back into $\chi_0$, and abbreviate $y \equiv \hbar (\omega + i \eta) / (4 \xi_F)$ and $c \equiv \cos{\theta}$, giving:

$$\begin{aligned} \chi_0 &= \frac{k_F^3}{8 \pi^2 \xi_F} \int_0^1 \! \int_{-1}^1 \bigg( \frac{1}{-x^2 - x \kappa c + y} + \frac{1}{-x^2 - x \kappa c - y} \bigg) \kappa^2 \dd{c} \dd{\kappa} \end{aligned}$$

From now on, we need to treat $\chi_0$’s real and imaginary parts separately using the Sokhotski-Plemelj theorem. The following is valid for the real part, but is (for now) still written as a complex expression. The point is that, to the best of my knowledge, it is not practical to extract the imaginary part from the “real part’s complex result” so to speak. So, hiding the principal value:

$$\begin{aligned} \chi_0 &= \frac{-k_F^3}{8 \pi^2 \xi_F} \int_0^1 \! \int_{-1}^1 \bigg( \frac{1}{x \kappa c + x^2 - y} + \frac{1}{x \kappa c + x^2 + y} \bigg) \kappa^2 \dd{c} \dd{\kappa} \\ &= \frac{-m k_F}{4 \pi^2 \hbar^2} \int_0^1 \bigg[ \frac{1}{x \kappa} \ln(x \kappa c + x^2 - y) + \frac{1}{x \kappa} \ln(x \kappa c + x^2 + y) \bigg]_{-1}^1 \kappa^2 \dd{\kappa} \\ &= \frac{-m k_F}{4 \pi^2 \hbar^2} \frac{1}{x} \int_0^1 \bigg( \ln\!\bigg( \frac{\kappa + x - y / x}{-\kappa + x - y / x} \bigg) + \ln\!\bigg( \frac{\kappa + x + y / x}{-\kappa + x + y / x} \bigg) \bigg) \kappa \dd{\kappa} \end{aligned}$$

We abbreviate this further by defining $z_{+} \equiv x + y / x$ and $z_{-} \equiv x - y / x$, yielding:

$$\begin{aligned} \chi_0 &= \frac{-m k_F}{4 \pi^2 \hbar^2} \frac{1}{x} \int_0^1 \ln\!\bigg( \frac{(\kappa + z_{-})(\kappa + z_{+})}{(\kappa - z_{-})(\kappa - z_{+})} \bigg) \kappa \dd{\kappa} \end{aligned}$$

Looking up this integral leads to the expression below. Note that all terms disappear for $\kappa = 0$:

$$\begin{aligned} \chi_0 &= \frac{-m k_F}{4 \pi^2 \hbar^2} \frac{1}{x} \bigg[ \kappa (z_{+} \!+\! z_{-}) - z_{+}^2 \tanh^{-1}\!\bigg(\!\frac{\kappa}{z_{+}}\!\bigg) - z_{-}^2 \tanh^{-1}\!\bigg(\!\frac{\kappa}{z_{-}}\!\bigg) + \frac{1}{2} \kappa^2 \ln\!\bigg( \frac{(\kappa + z_{+})(\kappa + z_{-})}{(\kappa - z_{-})(\kappa - z_{+})} \bigg) \bigg]_0^1 \\ &= \frac{-m k_F}{\pi^2 \hbar^2} \frac{1}{4 x} \bigg( 2 x - z_{+}^2 \tanh^{-1}\!\bigg(\!\frac{1}{z_{+}}\!\bigg) - z_{-}^2 \tanh^{-1}\!\bigg(\!\frac{1}{z_{-}}\!\bigg) + \frac{1}{2} \ln\!\bigg( \frac{(1 + z_{+})(1 + z_{-})}{(1 - z_{-})(1 - z_{+})} \bigg) \bigg) \\ &= \frac{-m k_F}{\pi^2 \hbar^2} \frac{1}{8 x} \bigg( 4 x - z_{+}^2 \ln\!\bigg( \frac{z_{+} + 1}{z_{+} - 1} \bigg) - z_{-}^2 \ln\!\bigg( \frac{z_{-} + 1}{z_{-} - 1} \bigg) + \ln\!\bigg( \frac{(z_{+} + 1)(z_{-} + 1)}{(z_{+} - 1)(z_{-} - 1)} \bigg) \bigg) \\ &= \frac{-m k_F}{\pi^2 \hbar^2} \bigg( \frac{1}{2} + \frac{1}{8 x} (1 - z_{+}^2) \ln\!\bigg( \frac{z_{+} + 1}{z_{+} - 1} \bigg) + \frac{1}{8 x} (1 - z_{-}^2) \ln\!\bigg( \frac{z_{-} + 1}{z_{-} - 1} \bigg) \bigg) \end{aligned}$$

Finally, we also know that the density of states $g$ per unit volume of a Fermi gas is given by:

$$\begin{aligned} g(\xi_F) = \frac{2}{4 \pi^2} \bigg( \frac{2 m}{\hbar^2} \bigg)^{3/2} \xi_F^{1/2} = \frac{m k_F}{\pi^2 \hbar^2} \end{aligned}$$

Now we officially take the real part, and arrive at the following expression for the response function of non-interacting electrons:

$$\begin{aligned} \boxed{ \Real{\chi_0(\vb{q}, \omega)} = - g(\xi_F) \: \mathcal{F}\bigg(\frac{|\vb{q}|}{2 k_F}, \: \frac{\hbar \omega}{4 \xi_F}\bigg) } \end{aligned}$$

Where $\mathcal{F}$ (which is sometimes also called the Lindhard function), is defined like so:

$$\begin{aligned} \boxed{ \mathcal{F}(x, y) \equiv \frac{1}{2} + \frac{1}{8 x} \bigg( 1 - \Big(x + \frac{y}{x}\Big)^{\!2} \bigg) \ln\!\bigg| \frac{x + y/x + 1}{x + y/x - 1} \bigg| + \frac{1}{8 x} \bigg( 1 - \Big(x - \frac{y}{x}\Big)^{\!2} \bigg) \ln\!\bigg| \frac{x - y/x + 1}{x - y/x - 1} \bigg| } \end{aligned}$$

Now we move on to $\chi_0$’s imaginary part. Recall the following expression from earlier:

$$\begin{aligned} \chi_0 &= \frac{- k_F^3}{8 \pi^2 \xi_F} \frac{1}{x} \int_0^1 \! \int_{-1}^1 \bigg( \frac{1}{x + y/x + \kappa c} + \frac{1}{x - y/x + \kappa c} \bigg) \kappa^2 \dd{c} \dd{\kappa} \end{aligned}$$

The Sokhotski-Plemelj theorem then states that the imaginary part is as shown below. There is a nasty detail here regarding the positive infinitesimal $\eta \to 0^{+}$: the first input term’s denominator contains $+ i \eta$, so the theorem can be applied normally, but the second term has $-i \eta$, meaning the theorem’s output must be negated!

$$\begin{aligned} \Imag{\chi_0} &= \frac{-2 m k_F^3}{8 \pi^2 \hbar^2 k_F^2} \frac{1}{x} \int_0^1 \! \int_{-1}^1 \bigg( \!-\! \pi \delta(z_{+} + \kappa c) + \pi \delta(z_{-} + \kappa c) \bigg) \kappa^2 \dd{c} \dd{\kappa} \\ &= \frac{\pi m k_F}{4 \pi^2 \hbar^2} \frac{1}{x} \int_0^1 \bigg[ \frac{1}{\kappa} \Theta(z_{+} + \kappa c) - \frac{1}{\kappa} \Theta(z_{-} + \kappa c) \bigg]_{-1}^1 \kappa^2 \dd{\kappa} \\ &= \frac{\pi}{4 x} \: g(\xi_F) \int_0^1 \bigg( \Theta(z_{+} + \kappa) - \Theta(z_{+} - \kappa) - \Theta(z_{-} + \kappa) + \Theta(z_{-} - \kappa) \bigg) \kappa \dd{\kappa} \end{aligned}$$

These integrals are easily evaluated according to the following straightforward relations:

$$\begin{aligned} \int_0^1 \Theta(\alpha - \kappa) \: \kappa \dd{\kappa} &= \frac{1}{2} \big( \alpha^2 - (\alpha^2 - 1) \: \Theta(\alpha - 1) \big) \: \Theta(\alpha) \\ \int_0^1 \Theta(\alpha + \kappa) \: \kappa \dd{\kappa} &= \frac{1}{2} \big( \alpha^2 \: \Theta(\alpha) - \alpha^2 + 1 \big) \: \Theta(\alpha + 1) \end{aligned}$$

Let us assume that $y > 0$ for now, such that $z_{+} > 0$ and $z_{+} > z_{-}$. We can then naively distinguish several options to be considered:

1. $z_{+} > z_{-} > 1 > 0 > -1$.
2. $z_{+} > 1 > z_{-} > 0 > -1$.
3. $z_{+} > 1 > 0 > z_{-} > -1$.
4. $z_{+} > 1 > 0 > -1 > z_{-}$.
5. $1 > z_{+} > z_{-} > 0 > -1$.
6. $1 > z_{+} > 0 > z_{-} > -1$.
7. $1 > z_{+} > 0 > -1 > z_{-}$.

It can be shown that #7 is impossible due to the definitions of $z_{+}$ and $z_{-}$. A tedious investigation of the six remaining cases then yields:

1. In cases #1 or #4, or, equivalently, if both $|z_{+}| > 1$ and $|z_{-}| > 1$, the result vanishes:

   $$\begin{aligned} \Imag{\chi_0} = 0 \end{aligned}$$

   An alternative statement of those conditions can be found, namely:

   $$\begin{aligned} \begin{cases} y > x^2 + x &\mathrm{if} \: x < 1 \\ x^2 - x > y > 0 &\mathrm{if} \: x > 1 \end{cases} \end{aligned}$$

2. In cases #2 or #3, or, equivalently, if both $|z_{+}| > 1$ and $|z_{-}| < 1$, the result is parabolic in $y$. These conditions can be restated as $x^2 + x > y > |x^2 - x|$.

   $$\begin{aligned} \Imag{\chi_0} &= -\frac{\pi}{8 x} \: g(\xi_F) \: (1 - z_{-}^2) \end{aligned}$$

3. In cases #5 or #6, or, equivalently, if both $|z_{+}| < 1$ and $|z_{-}| < 1$ the result is linear in $y$. These conditions can be restated as $-(x^2 - x) > y > 0$.

   $$\begin{aligned} \Imag{\chi_0} &= \frac{\pi}{8 x} \: g(\xi_F) \: (z_{-}^2 - z_{+}^2) = - \frac{\pi}{2 x} \: g(\xi_F) \: y \end{aligned}$$

Repeating the analysis for $y < 0$ and putting everything together, we arrive at the following formula for the imaginary part of the non-interacting response function:

$$\begin{gathered} \boxed{ \Imag{\chi_0(\vb{q}, \omega)} = - g(\xi_F) \: \mathcal{G}\bigg( \frac{|\vb{q}|}{2 k_F}, \: \frac{\hbar \omega}{4 \xi_F} \bigg) } \\ \boxed{ \mathcal{G}(x, y) \equiv \begin{cases} 0 &\mathrm{if} \: |y| > x^2 + x > 0 \\ \sgn(y) \displaystyle\frac{\pi}{8 x} \bigg( 1 - \Big( x - \displaystyle\frac{|y|}{x} \Big)^2 \bigg) &\mathrm{if} \: x^2 + x > |y| > |x^2 - x| \\ \displaystyle\frac{\pi}{2 x} y &\mathrm{if} \: \!-\!(x^2 - x) > |y| > 0 \\ 0 &\mathrm{if} \: x^2 - x > |y| > 0 \end{cases} } \end{gathered}$$

We are being a bit sloppy here by ignoring the boundaries, e.g. $|y| = x^2 - x$. Fortunately, they are all naturally continuous, so either of the two neighboring formulas can be used.

Interpretation

From this, we would like to get the dielectric function $\varepsilon_r$. When an external electric field is applied, the electrons respond and thereby modify the net field inside the material. We include this effect in our energy potential $U$, such that the net electric potential $\Phi_\mathrm{tot} = U / e_0$ (where $e_0 < 0$). This is not the same as including direct electron-electron interactions!

We thus have a linear relation for the induced particle density $\delta\!\Expval{\hat{n}(\vb{q}, \omega)} = \chi_0(\vb{q}, \omega) \: U(\vb{q})$ thanks to the convolution theorem. The corresponding induced charge density is given by $\rho_\mathrm{ind} = e_0^2 \chi_0 \Phi_\mathrm{tot}$, so the standard formula for $\varepsilon_r$ yields:

$$\begin{aligned} \boxed{ \varepsilon_r(\vb{q}, \omega) = 1 - U_{ee}(\vb{q}) \: \chi_0(\vb{q}, \omega) } \end{aligned}$$

Where $U_{ee}(\vb{q}) = e_0^2 / (\varepsilon_0 |\vb{q}|^2)$ is Coulomb repulsion. This is the Lindhard dielectric function of a free non-interacting electron gas, at any temperature and for any dimensionality.

Lindhard theory can reproduce several notable discoveries in various limits. For long wavelengths $\vb{q} \to 0$, i.e. $x \to 0$, it can be shown that $\mathcal{F}(x, y) \approx - x^2 / (3 y^2)$, and therefore:

$$\begin{aligned} \chi_0(\vb{q}, \omega) &\approx \frac{m k_F}{\pi^2 \hbar^2} \frac{|\vb{q}|^2}{4 k_F^2} \frac{16 \xi_F^2}{3 \hbar^2 \omega^2} = \frac{k_F^3}{3 \pi^2 m} \frac{|\vb{q}|^2}{\omega^2} \end{aligned}$$

For a Fermi gas, we know that $k_F^3 = 3 \pi^2 n$, where $n$ is the unperturbed density. Upon insertion into $\varepsilon_r$, this gives the same result as the Drude model, predicting plasma oscillations, also known as plasmons (when treated as quantized particles) or cold Langmuir waves:

$$\begin{aligned} \boxed{ \varepsilon_r(0, \omega) = 1 - \frac{\omega_p^2}{\omega^2} } \qquad \qquad \boxed{ \omega_p \equiv \sqrt{\frac{e_0^2 n}{\varepsilon_0 m}} } \end{aligned}$$

Meanwhile, in the static limit $\omega \to 0$, i.e. $y \to 0$, we arrive at a function that is also encountered when using the Hartree-Fock method to approximate the ground state of jellium:

$$\begin{aligned} \boxed{ \mathcal{F}(x, 0) = \frac{1}{2} + \frac{1}{4 x} (1 - x^2) \ln\!\bigg| \frac{x + 1}{x - 1} \bigg| } \end{aligned}$$

If we assume small $|\vb{q}|$ and $x$ again, it can be shown that $\mathcal{F}(x, 0) \approx 1$, in which case:

$$\begin{aligned} \chi_0(\vb{q}, 0) = - \frac{m k_F^3}{\pi^2 \hbar^2 k_F^2} = - \frac{3 \pi^2 m n}{\pi^2 \hbar^2 k_F^2} = - \frac{3 m n}{\hbar^2 k_F^2} = - \frac{3 n}{2 \xi_F} \end{aligned}$$

The resulting dielectric function also appears in the derivation of Thomas-Fermi screening, which can be regarded as the zero-temperature limit of Debye screening:

$$\begin{aligned} \boxed{ \varepsilon_r(\vb{q}, 0) = 1 + \frac{k_\mathrm{TF}^2}{|\vb{q}|^2} } \qquad \qquad \boxed{ k_\mathrm{TF} \equiv \sqrt{\frac{3 e_0^2 n}{2 \varepsilon_0 \xi_F}} } \end{aligned}$$

But probably the greatest achievement of Lindhard theory is its appearance in the random phase approximation (where $\chi_0$ is sometimes known as the retarded pair-bubble) for treating electron-electron Coulomb interactions in the Fermi gas model. So although we neglected interactions in our derivation, the result is an essential building block of interacting models.

References

1. K.S. Thygesen, Advanced solid state physics: linear response theory, 2013, unpublished.
2. H. Bruus, K. Flensberg, Many-body quantum theory in condensed matter physics, 2016, Oxford.
3. G. Grosso, G.P. Parravicini, Solid state physics, 2nd edition, Elsevier.
4. G.F. Giuliani, G. Vignale, Quantum theory of the electron liquid, 2005, Cambridge.
5. B. Mihaila, Lindhard function of a $d$-dimensional Fermi gas, 2011, arXiv.

---