Categories: Physics, Quantum mechanics.

# Lindhard function

The Lindhard function describes the response of jellium (i.e. a free electron gas) to an external perturbation, and is a quantum-mechanical alternative to the Drude model.

We start from the Kubo formula for the electron density operator $\hat{n}$, which describes the change in $\Expval{\hat{n}}$ due to a time-dependent perturbation $\hat{H}_1$:

\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'} \end{aligned}

Where the subscript $I$ refers to the interaction picture, and the expectation $\Expval{}_0$ is for a thermal equilibrium before the perturbation was applied. Now consider a harmonic $\hat{H}_1$:

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

Where $S$ is the SchrĂ¶dinger picture, $\eta$ is a positive infinitesimal to ensure convergence later, and $U(\vb{r})$ is an arbitrary potential function. The Kubo formula becomes:

\begin{aligned} \delta\!\Expval{ {\hat{n}}}\!(\vb{r}, 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}

Here, $\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$) is 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)^{2D}} \iiint \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)^{-2 D}$, with $D$ being the number of spatial dimensions. 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)^{2D}} \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. If we choose an arbitrary (hyper)cube of volume $V$ in real space, then clearly $\int_V \dd{\vb{r}'} = V$. Therefore:

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

For $V \to \infty$ we get a Dirac delta function, but in fact the conclusion holds for finite $V$ too:

\begin{aligned} \chi(\vb{q}; t, t') &= -\frac{i}{\hbar} \Theta(t \!-\! t') \frac{1}{V} \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} \Theta(t \!-\! t') \frac{1}{V} \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, $\chi$ only depends on the time difference $t - t'$. Note that $\delta{\Expval{\hat{n}}}$ already has the form of a Fourier transform, 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-electron eigenstates of $\hat{H}_0$, and $Z$ is the grand partition function. According to the convolution theorem $\delta{\Expval{\hat{n}}}(\vb{q}, \omega) = \chi(\vb{q}, \omega) \: U(\vb{q})$. In anticipation, we swap $\nu$ and $\nu''$ in the second term, so the general response function is written as:

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

All operators are in the SchrĂ¶dinger picture from now on, hence we dropped the subscript $S$.

To proceed, we need to rewrite $\hat{n}(\vb{q})$ somehow. If we neglect electron-electron interactions, the single-particle states are simply plane waves, in which case:

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

Starting from the general definition of $\hat{n}$, 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}$:

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

If we impose periodic boundary conditions on our $D$-dimensional hypercube of volume $V$, then $\vb{k}$ becomes discrete, with per-value spacing $2 \pi / V^{1/D}$ along each axis.

Consequently, each orbital $\psi_\vb{k}$ uniquely occupies a volume $(2 \pi)^D / V$ in $\vb{k}$-space, so we make the approximation $\sum_{\vb{k}} \approx V / (2 \pi)^D \int_{-\infty}^\infty \dd{\vb{k}}$. This becomes exact for $V \to \infty$, in which case $\vb{k}$ also becomes continuous again, which is what we want for jellium.

We apply this standard trick from condensed matter physics to $\hat{n}$, and $V$ cancels out:

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

The summation variable $\vb{k}$ has an associated spin $\sigma$, and $\hat{n}$ does not carry any spin.

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

\begin{aligned} \chi_0 &= \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 &= \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 &= \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 a thermal expectation $\Expval{}_0$:

\begin{aligned} \chi_0 &= \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}

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

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.

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}

From this, we would like to get the dielectric function $\varepsilon_r$. Recall its definition, where $U_\mathrm{tot}$, $U_\mathrm{ext}$, and $U_\mathrm{ind}$ are the total, external and induced potentials, respectively:

\begin{aligned} U_\mathrm{tot} = U_\mathrm{ext} + U_\mathrm{ind} = \frac{U_\mathrm{ext}}{\varepsilon_r} \end{aligned}

Note that these are all energy potentials: this choice is justified because all energy potentials are caused by electric fields in this case. The electric potential is recoverable as $\Phi_\mathrm{tot} = q_e U_\mathrm{tot}$, where $q_e < 0$ is the charge of an electron.

From the Lindhard response function $\chi_0$, we get the induced particle density offset $\delta{\Expval{\hat{n}}}$ caused by a potential $U$. The density $\delta{\Expval{\hat{n}}}$ should be self-consistent, implying $U = U_\mathrm{tot}$. In other words, we have a linear relation $\delta{\Expval{\hat{n}}} = \chi_0 U_\mathrm{tot}$, so the standard formula for $\varepsilon_r$ gives:

\begin{aligned} \boxed{ \varepsilon_r(\vb{q}, \omega) = 1 - \frac{U_{ee}(\vb{q})}{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}

Where $U_{ee}(\vb{q}) = q_e^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.

## 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.