path: root/source/know/concept/lindhard-function/index.md

---
title: "Lindhard function"
sort_title: "Lindhard function"
date: 2026-09-24 # Originally 2022-01-24, major rewrite
categories:
- Physics
- Quantum mechanics
layout: "concept"
---

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](/know/concept/kubo-formula/) below,
where the subscript $$I$$ refers to
the [interaction picture](/know/concept/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](/know/concept/heaviside-step-function/),
and the response $$\chi$$ is the density-density correlation function
(i.e. a two-particle [Green's function](/know/concept/greens-functions/)):

$$\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](/know/concept/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](/know/concept/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](/know/concept/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](/know/concept/grand-canonical-ensemble/).
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](/know/concept/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](/know/concept/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}$$


{% include proof/start.html id="proof-density" -%}
Starting from the definition of $$\hat{n}(\vb{r})$$,
we write out the field operators $$\hat{\Psi}(\vb{r})$$,
and insert the known non-interacting single-electron orbitals
$$\psi_\vb{k}(\vb{r}) = e^{i \vb{k} \cdot \vb{r}} / \sqrt{V}$$,
omitting spin:

$$\begin{aligned}
    \hat{n}(\vb{r})
    \equiv \hat{\Psi}{}^\dagger(\vb{r}) \hat{\Psi}(\vb{r})
    = \sum_{\vb{k} \vb{k}'} \psi_{\vb{k}}^*(\vb{r}) \: \psi_{\vb{k}'}(\vb{r})\: \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k}'}
    = \frac{1}{V} \sum_{\vb{k} \vb{k}'} e^{i (\vb{k}' - \vb{k}) \cdot \vb{r}} \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k}'}
\end{aligned}$$

Taking the Fourier transfom yields a Dirac delta function $$\delta$$:

$$\begin{aligned}
    \hat{n}(\vb{q})
    = \frac{1}{V} \int_{-\infty}^\infty
    \sum_{\vb{k} \vb{k}'} \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k}'} \: e^{i (\vb{k}' - \vb{k} - \vb{q})\cdot \vb{r}} \dd{\vb{r}}
    = \frac{(2 \pi)^D}{V} \sum_{\vb{k} \vb{k}'} \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k}'} \: \delta(\vb{k}' \!-\! \vb{k} \!-\! \vb{q})
\end{aligned}$$

Thanks to our periodic boundary conditions,
[this trick](/know/concept/discrete-spectrum-summation/)
lets us turn a sum into an integral:

$$\begin{aligned}
    \hat{n}(\vb{q})
    &= \frac{(2 \pi)^D}{V} \frac{V}{(2 \pi)^D} \sum_{\vb{k}} \int_{-\infty}^\infty
    \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k}'} \: \delta(\vb{k}' \!-\! \vb{k} \!-\! \vb{q}) \dd{\vb{k}'}
    = \sum_{\vb{k}} \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} + \vb{q}}
\end{aligned}$$

For negated arguments, we simply define $$\vb{k}' \equiv \vb{k} - \vb{q}$$
to show that $$\hat{n}(-\vb{q}) = \hat{n}{}^\dagger(\vb{q})$$,
which can also be understood as a consequence of $$\hat{n}(\vb{r})$$ being real:

$$\begin{aligned}
    \hat{n}(-\vb{q})
    = \sum_{\vb{k}} \hat{c}_{\vb{k}}^\dagger \hat{c}_{\vb{k} - \vb{q}}
    = \sum_{\vb{k}'} \hat{c}_{\vb{k}' + \vb{q}}^\dagger \hat{c}_{\vb{k}'}
    = \hat{n}^\dagger(\vb{q})
\end{aligned}$$
{% include proof/end.html id="proof-density" %}


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}$$


{% include proof/start.html id="proof-commutator" -%}
In general, for any single-particle states labeled by $$m$$, $$n$$, $$o$$ and $$p$$, we have:

$$\begin{aligned}
    \comm{\hat{c}_m^\dagger \hat{c}_n}{\hat{c}_o^\dagger \hat{c}_p}
    &= \hat{c}_m^\dagger \hat{c}_n \hat{c}_o^\dagger \hat{c}_p - \hat{c}_o^\dagger \hat{c}_p \hat{c}_m^\dagger \hat{c}_n
    \\
    &= \hat{c}_m^\dagger \big( \acomm{\hat{c}_n}{\hat{c}_o^\dagger} - \hat{c}_o^\dagger \hat{c}_n \big) \hat{c}_p
    - \hat{c}_o^\dagger \big( \acomm{\hat{c}_p}{\hat{c}_m^\dagger} - \hat{c}_m^\dagger \hat{c}_p \big) \hat{c}_n
\end{aligned}$$

Using the standard fermion anticommutation relations, this becomes:

$$\begin{aligned}
    \comm{\hat{c}_m^\dagger \hat{c}_n}{\hat{c}_o^\dagger \hat{c}_p}
    &= \hat{c}_m^\dagger \big( \delta_{no} - \hat{c}_o^\dagger \hat{c}_n \big) \hat{c}_p
    - \hat{c}_o^\dagger \big( \delta_{pm} - \hat{c}_m^\dagger \hat{c}_p \big) \hat{c}_n
    \\
    &= \hat{c}_m^\dagger \hat{c}_p \: \delta_{no} - \hat{c}_m^\dagger \hat{c}_o^\dagger \hat{c}_n \hat{c}_p
    - \hat{c}_o^\dagger \hat{c}_n \: \delta_{pm} + \hat{c}_o^\dagger \hat{c}_m^\dagger \hat{c}_p \hat{c}_n
    \\
    &= \hat{c}_m^\dagger \hat{c}_p \: \delta_{no} - \hat{c}_o^\dagger \hat{c}_n \: \delta_{pm}
\end{aligned}$$

In this case, $$m = p = \vb{k}$$ and $$n = o = \vb{k} \!+\! \vb{q}$$,
so the Kronecker deltas are unnecessary.
{% include proof/end.html id="proof-commutator" %}


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](/know/concept/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](/know/concept/discrete-spectrum-summation/)
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](/know/concept/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](/know/concept/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](/know/concept/cauchy-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](/know/concept/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](/know/concept/dielectric-function/) $$\varepsilon_r$$.
When an external [electric field](/know/concept/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](/know/concept/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](/know/concept/drude-model/),
predicting *plasma oscillations*, also known as
*plasmons* (when treated as quantized particles)
or [*cold Langmuir waves*](/know/concept/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](/know/concept/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](/know/concept/debye-length/):

$$\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](/know/concept/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*](https://doi.org/10.48550/arXiv.1111.5337),
    2011, arXiv.