B is for Bunching

Bose or Fermi statistics gives rise to correlations in the positions of particles in a many body system, even in product states. These correlations affect the energy of many particle states, and provide the basis of the simplest approximate theories.

Reading: Nazarov and Danon (2013), Baym (2018)

---

1 Correlation Functions

In Lecture 1 we introduced the pair distribution function

$$\begin{array}{}\text{(1)}& {\rho }_2(x_1,x_2)=N(N-1)\int d{x}_3\dots d{x}_N{|\mathrm{\Psi }(x_1,x_2,\dots ,x_N)|}^2,\end{array}$$

which measures the likelihood of finding a pair of particles at $x_1$ and $x_2$. We also saw that ${\rho }_2(x_1,x_2)$ could be related to the expectation of the product of the density operators $\rho (x)$ at different points.

$$\begin{array}{}\text{(2)}& {\rho }_2(x_1,x_2)=⟨\mathrm{\Psi }|\rho (x_1)\rho (x_2)|\mathrm{\Psi }⟩-{\rho }_1(x_1)\delta (x_1-x_2),\end{array}$$

where ${\rho }_1(x)=⟨\mathrm{\Psi }|\rho (x)|\mathrm{\Psi }⟩$. ${\rho }_2(x,y)$ is the simplest example of a correlation function. Note also that the quantity $⟨\mathrm{\Psi }|{\rho }_q{\rho }_{-q}|\mathrm{\Psi }⟩$, used to quantify the crystalline order in a solid in Lecture 3, is just the Fourier transform of ${\rho }_2(x_1,x_2)$. Higher correlation functions, involving products of more than two density operators, provide more detailed information on the distribution of the particles.

1.1 Correlation Functions in Second Quantization

We evaluated ${\rho }_2(x,y)$ for the 1D Fermi gas:

$$\begin{array}{}\text{(3)}& {\rho }_2(x,y)=n^2[1-{\left(\frac{\sin [k_{\text{F}}(x-y)]}{k_{\text{F}}(x-y)}\right)}^2].\end{array}$$

Let’s see how to reproduce this result using second quantization. ${\rho }_2(x,y)$ can be written in second quantized notation as

$$\begin{array}{}\text{(4)}& {\rho }_2(x,y)=⟨\mathrm{\Psi }|{\psi }^{†}(x){\psi }^{†}(y){\psi }^{\phantom{†}}(y){\psi }^{\phantom{†}}(x)|\mathrm{\Psi }⟩.\end{array}$$

Check

Perhaps the simplest way to get this formula is to use the correspondence that we found in Lecture 5 $$\begin{array}{rl}|\mathrm{\Psi }⟩& ⟷\mathrm{\Psi }(x_1,\dots ,x_N)\\ \psi (X)|\mathrm{\Psi }⟩& ⟷\sqrt{N}\mathrm{\Psi }(X,x_1,\dots ,x_{N-1})\end{array}$$

Show that Equation 1 implies Equation 4.

Equation 4 differs from $⟨\mathrm{\Psi }|\rho (x)\rho (y)|\mathrm{\Psi }⟩$ by the ordering of operators. Using the canonical commutation relations for the fields gives the relation Equation 2 for both bosons and fermions.

Let’s evaluate Equation 4 for the ground state of the Fermi gas. We will find the result for a general product state, using the same strategy as we used in the last lecture for the single particle density matrix. That is, we substitute the expansion of the fields ${\psi }^{\phantom{†}}(x)$ and ${\psi }^{†}(x)$ in terms of the orthonormal single particle states making up the product state

$$\begin{array}{r}{\psi }^{\phantom{†}}(x)=\sum _{\beta }{\phi }_{\beta }^{}(x)a_{\beta }^{\phantom{†}},\\ {\psi }^{†}(x)=\sum _{\beta }{\phi }_{\beta }^{\ast }(x)a_{\beta }^{†}.\end{array}$$

This gives

$${\rho }_2(x,y)=\sum _{\alpha ,\beta ,\gamma ,\delta }{\phi }_{\alpha }^{\ast }(x){\phi }_{\beta }^{\ast }(y){\phi }_{\gamma }^{}(y){\phi }_{\delta }^{}(x)⟨\mathrm{\Psi }|a_{\alpha }^{†}{a}_{\beta }^{†}{a}_{\gamma }^{\phantom{†}}{a}_{\delta }^{\phantom{†}}|\mathrm{\Psi }⟩.$$

If we are considering the expectation in a state of the form $|\mathbf{N}⟩$, we can see that an annihilation operator for a given single particle state must be accompanied by a creation operator for the same state. There are therefore two possibilities

$$\begin{array}{rl}& \alpha =\delta ,\beta =\gamma ,\text{ or }\\ & \alpha =\gamma ,\beta =\delta ,\end{array}$$

which give rise to two groups of terms

$$\begin{array}{rl}⟨\mathbf{N}|a_{\alpha }^{†}{a}_{\gamma }^{†}{a}_{\gamma }^{\phantom{†}}{a}_{\alpha }^{\phantom{†}}|\mathbf{N}⟩& =N_{\alpha }{N}_{\gamma }\\ ⟨\mathbf{N}|a_{\alpha }^{†}{a}_{\gamma }^{†}{a}_{\alpha }^{\phantom{†}}{a}_{\gamma }^{\phantom{†}}|\mathbf{N}⟩& =\pm N_{\alpha }{N}_{\gamma }\text{if }\alpha \ne \gamma ,\end{array}$$

the $\pm$ corresponding to bosons and fermions respectively. Thus we have

$$\begin{array}{}\text{(5)}& \begin{array}{r}{\rho }_2(x,y)=\sum _{\alpha ,\beta }{N}_{\alpha }{N}_{\beta }[|{\phi }_{\alpha }(x){|}^2|{\phi }_{\beta }(y){|}^2\pm {\phi }_{\alpha }^{\ast }(x){\phi }_{\alpha }^{}(y){\phi }_{\beta }^{\ast }(y){\phi }_{\beta }^{}(x)].\end{array}\end{array}$$

You might notice that this expression weights the case $\alpha =\beta =\gamma =\delta$ by a factor $2N_{\alpha }^2$ when it should have $N_{\alpha }(N_{\alpha }-1)$ in the case of bosons, and zero for fermions. On the other hand, such cases amount to a sum over a single index, where the general case is a sum over a double index. It seems reasonable, then, that this accounting error is not important in the thermodynamic limit. For instance, in the case that ${\phi }_{\alpha }(x)$ are plane wave states we have the usual prescription

$$\sum _{\alpha }(\cdots )⟶L\int (\cdots )\frac{dk}{2\pi },$$

assuming the integrand is smooth. In this case the error in Equation 5 is a factor of $L^{-1}$ smaller than what we retain. We do however have to be careful in Bose condensates, where one state has a finite fraction of the particles and this argument does not apply.

We can express the result Equation 5 in terms of the density and density matrix as

$$\begin{array}{}\text{(6)}& {\rho }_2(x,y)={\rho }_1(x){\rho }_1(y)\pm g(x,y)g(y,x),\end{array}$$

which, in the case of the ground state of the Fermi gas, reproduces Equation 3. We see that the correlation function vanishes as the separation $x-y\to 0$, because $g(x,x)=⟨\rho (y)⟩$. This is, of course, another manifestation of the exclusion principle: it is not possible for two fermions to sit on top of each other. The scale of the ‘hole’ in the correlation function is of course set by the mean interparticle separation, or Fermi wavelength. Note also the decaying oscillations, indicating liquid-like correlations in the positions of the particles. These are sometimes known as Friedel oscillations.

Correlation function ${\rho }_2(x,0)$ for the Fermi gas.

For bosons the situation is very different. If $g(x,y)\to 0$ as $|x-y|\to \mathrm{\infty }$, the value of the correlation function as $|x-y|\to 0$ is twice the value at $|x-y|\to \mathrm{\infty }$. This characteristic behavior is often termed bunching: a pair of bosons is more likely to be found at two nearby points than at two distant points.

Nothing about the result Equation 6 is special to 1D of course: one just has to recalculate the density and density matrix. Remember, though, that it does only apply to product states.

1.2 The Hanbury Brown and Twiss Effect

The result for the density correlations Equation 5 or Equation 6 can be viewed as a kind of interference effect that shows up in the correlations of the intensity of a quantum waves, even when there is no interference in the intensity itself. To illustrate this interpretation, we consider a classic experiment from with Bose condensates described in Andrews et al. (1997).

Consider a gas of $N$ noninteracting bosons occupying the lowest energy level of some potential well: a Bose condensate. If the ground state wavefunction is ${\phi }_0(\mathbf{r})$, the $N$-body wavefunction for such a state is

$$\mathrm{\Psi }({\mathbf{r}}_1,{\mathbf{r}}_2,\dots ,{\mathbf{r}}_N)=\prod _i^N{\phi }_0({\mathbf{r}}_i),$$

which we can write in second quantized notation as

$$|\mathrm{\Psi }⟩=\frac{1}{\sqrt{N!}}{\left(a_0^{†}\right)}^N|\text{VAC}⟩,$$

where $a_0^{†}$ creates a particle in the state ${\phi }_0(\mathbf{r})$. Imagine that we took another well, also filled with $N$ bosons, and placed it alongside the first. If we switch off the potentials at some instant, the particles will fly out, with wavefunctions orginating in the two wells overlapping. Precisely this experiment was reported in Andrews et al. (1997). What do we expect to see?

Let us denote by ${\phi }_L(\mathbf{r})$ and ${\phi }_R(\mathbf{r})$ the ground states of two spatially separated potential wells. First, consider a state where each boson is in a superposition of ${\phi }_L(\mathbf{r})$ and ${\phi }_R(\mathbf{r})$. Such a situation could arise by starting from a single well and adiabatically splitting it in two. We can write such a state as

$$|{\overline{N}}_L,{\overline{N}}_R{⟩}_{\theta }\equiv \frac{1}{\sqrt{N!}}{[\sqrt{\frac{{\overline{N}}_L}{N}}{e}^{-i\theta /2}{a}_L^{†}+\sqrt{\frac{{\overline{N}}_R}{N}}{e}^{i\theta /2}{a}_R^{†}]}^N|\text{VAC}⟩,$$

where ${\overline{N}}_{L,R}$ are the expectation values of particle number in each state $N={\overline{N}}_L+{\overline{N}}_R$. We allow the system to evolve for some time $t$, so that the two `clouds’ begin to overlap (typically achieved by allowing free expansion i.e. turning off the confining potentials). If we consider the time evolution in the Heisenberg picture then, as we saw last time, the field operator obeys the free particle Schrödinger equation (ignoring interactions)

$$i\frac{\partial {\psi }^{\phantom{†}}(\mathbf{r},t)}{\partial t}=-\frac{1}{2m}{\mathrm{\nabla }}^2{\psi }^{\phantom{†}}(\mathbf{r},t).$$

We write the field operator as

$${\psi }^{\phantom{†}}(\mathbf{r})={\phi }_L(\mathbf{r},t)a_L^{\phantom{†}}+{\phi }_R(\mathbf{r},t)a_R^{\phantom{†}}+\cdots ,$$

where the wavefunctions ${\phi }_{L/R}(\mathbf{r},t)$ are evolving freely, and the dots denote the other states in some complete orthogonal set that includes ${\phi }_L(\mathbf{r})$ and ${\phi }_R(\mathbf{r})$: we can ignore them because they are empty. A calculation of the expectation value of $\rho (\mathbf{r})={\psi }^{†}(\mathbf{r}){\psi }^{\phantom{†}}(\mathbf{r})$ gives

$$\begin{array}{}\text{(7)}& \begin{array}{r}{\rho }_1(\mathbf{r},t)={\overline{N}}_L|{\phi }_L(\mathbf{r},t){|}^2+{\overline{N}}_R|{\phi }_R(\mathbf{r},t){|}^2+\stackrel{\equiv {\rho }_{\mathrm{int}(\mathbf{r},t)}}{\overbrace{2\sqrt{{\overline{N}}_L{\overline{N}}_R}\mathrm{Re}{e}^{i\theta }{\phi }_L^{\ast }(\mathbf{r},t){\phi }_R(\mathbf{r},t)}}.\end{array}\end{array}$$

If the clouds begin to overlap, the last term in Equation 7 comes into play. Its origin is in quantum interference between the two coherent subsystems, showing that the relative phase has a real physical effect.

Check

Consider a Gaussian wavefunction of width $R_0$ at time $t=0$. Show (by substitution into the Schrödinger equation is fine) that this function evolves as

$$\begin{array}{}\text{(8)}& \phi (\mathbf{r},t)=\frac{1}{{\left(\pi R_t^2\right)}^{3/4}}\exp [-\frac{{\mathbf{r}}^2(1+it/m{R}_0^2))}{2R_t^2}],\end{array}$$

where

$$R_t^2=R_0^2+{\left(\frac{t}{m{R}_0}\right)}^2.$$

Equation 8 illustrates a very important point about the expansion of a gas. After a long period of expansion, the final density distribution is a reflection of the initial momentum distribution. This is simply because faster moving atoms fly further, so after time $t$ an atom with velocity $\mathbf{v}$ will be at position $\mathbf{r}=\mathbf{v}t$ from the center of the trap, provided that this distance is large compared to $R_0$, the initial radius of the gas. The $t\to \mathrm{\infty }$ limit of Equation 8 gives

$$|\phi (\mathbf{r},t\to \mathrm{\infty }){|}^2\propto \exp [-{\left(\frac{m{R}_0\mathbf{r}}{t}\right)}^2],$$

reflecting a Gaussian initial momentum distribution of width $R_0^{-1}$. Imaging the density distribution after expansion is one of the most commonly used experimental techniques in ultracold physics, and yields information about the momentum distribution $n(\mathbf{p})\equiv a_{\mathbf{p}}^{†}{a}_{\mathbf{p}}^{\phantom{†}}$ before expansion.

Consider the evolution of two Gaussian wavepackets with width $R_0$ at $t=0$, separated by a distance $d\gg R_0$

$${\phi }_{L,R}(\mathbf{r},t)=\frac{1}{{\left(\pi R_t^2\right)}^{3/4}}\exp [-\frac{{(\mathbf{r}\pm \mathbf{d}/2)}^2(1+it/m{R}_0^2))}{2R_t^2}],$$

The final term of Equation 7 is then

$$\begin{array}{rl}{\rho }_{\mathrm{int}}(\mathbf{r},t)& =A(\mathbf{r},t)\cos [\theta +\frac{\mathbf{r}\cdot \mathbf{d}}{m{R}_0^2{R}_t^2}t]\\ A(\mathbf{r},t)& =\frac{2\sqrt{{\overline{N}}_L{\overline{N}}_R}}{{\pi }^{3/2}{R}_t^3}\exp (-\frac{{\mathbf{r}}^2+{\mathbf{d}}^2/4}{R_t^2})\end{array}$$

The interference term therefore consists of regularly spaced fringes, with a separation at long times of $2\pi t/md$.

Now we imagine doing the same thing with two condensates of fixed particle number, which bear no phase relation to one another. The system is described by the product state (often called a Fock state in this context)

$$|N_L,N_R⟩\equiv \frac{1}{\sqrt{N_L!N_R!}}{\left(a_L^{†}\right)}^{N_L}{\left(a_R^{†}\right)}^{N_R}|\text{VAC}⟩.$$

Computing the density in the same way yields

$${\rho }_1(\mathbf{r},t)=N_L|{\phi }_L(\mathbf{r},t){|}^2+N_R|{\phi }_R(\mathbf{r},t){|}^2,$$

which differs from the previous result by the absence of the interference term.

This is not the end of the story, however. When we look at an absorption image of the gas, we are not looking at an expectation value of $\rho (\mathbf{r})$ but rather the measured value of some observable(s) $\rho (\mathbf{r})$. The expectation value just tells us the result we would expect to get if we repeated the same experiment many times and averaged the result. We get more information by thinking about the correlation function of the density at two different points.

An application of our result Equation 5 for the density correlations gives

$$\begin{array}{}\text{(9)}& \begin{array}{rl}{\rho }_2(\mathbf{r},{\mathbf{r}}^{\prime })& ={\rho }_1(\mathbf{r}){\rho }_1({\mathbf{r}}^{\prime })+N_L{N}_R{\phi }_L^{\ast }(\mathbf{r}){\phi }_R^{\ast }({\mathbf{r}}^{\prime }){\phi }_L({\mathbf{r}}^{\prime }){\phi }_R(\mathbf{r})\\ & +N_L{N}_R{\phi }_R^{\ast }(\mathbf{r}){\phi }_L^{\ast }({\mathbf{r}}^{\prime }){\phi }_R({\mathbf{r}}^{\prime }){\phi }_L(\mathbf{r}).\end{array}\end{array}$$

We see that the second line contains interference fringes, with the same spacing as before. The correlation function gives the relative probability of finding an atom at ${\mathbf{r}}^{\prime }$ if there is one at $\mathbf{r}$. We conclude that in each measurement of the density, fringes are present but with a phase that varies between measurements, even if the samples are identically prepared.

Fringes observed by interfering two Bose condensates (Source: Andrews et al. (1997))

The rather surprising implication is that predictions for measured quantities for a system in a Fock state are the same as in a relative phase state, but with a subsequent averaging over the phase.

Check

Prove this by showing that the density matrix

$$\rho ={\int }_0^{2\pi }\frac{d\theta }{2\pi }|{\overline{N}}_L,{\overline{N}}_R{⟩}_{\theta }⟨{\overline{N}}_R,{\overline{N}}_L{|}_{\theta }$$

coincides with that of a mixture of Fock states with binomial distribution of atoms into states ${\phi }_L$, ${\phi }_r$. At large $N$ this distribution becomes sharply peaked at occupations ${\overline{N}}_L$, ${\overline{N}}_R$.

The interference of two independent condensates was observed in 1997 in Andrews et al. (1997). The related question of whether two independent light sources give rise to interference was discussed much earlier in Magyar and Mandel (1963). The occurrence of interference fringes in a correlation function does not depend upon Bose condensation, although the phenomenon is very striking in this case because the fluctuations are parametrically as large as the ${\rho }_1(\mathbf{r}){\rho }_1({\mathbf{r}}^{\prime })$ term in Equation 9. The general phenomenon is known as the Hanbury Brown and Twiss effect, which is the work of two people (not three): Robert Hanbury Brown and Richard Q. Twiss. For the history and early applications of this effect, see Baym (1999) and Kleppner (2008).

2 Hartree–Fock Theory

We now apply these ideas to the approximate calculation of the energy of an interacting many body system.

2.1 The Hartree and Fock Potentials

Recall from last time that a two body interaction has the form

$${\hat{H}}_{\text{int.}}=\sum _{j<k}U({\mathbf{r}}_j-{\mathbf{r}}_k)=\frac{1}{2}\int d{\mathbf{r}}_{1}d{\mathbf{r}}_{2}U({\mathbf{r}}_1-{\mathbf{r}}_2){\psi }^{†}({\mathbf{r}}_1){\psi }^{†}({\mathbf{r}}_2){\psi }^{\phantom{†}}({\mathbf{r}}_2){\psi }^{\phantom{†}}({\mathbf{r}}_1).$$

Since

$$\sum _{j<k}U({\mathbf{r}}_j-{\mathbf{r}}_k)=\frac{1}{2}\int \sum _{j\ne k}\delta ({\mathbf{r}}_1-{\mathbf{r}}_j)\delta ({\mathbf{r}}_2-{\mathbf{r}}_k)U({\mathbf{r}}_1-{\mathbf{r}}_2)d{\mathbf{r}}_{1}d{\mathbf{r}}_2,$$

we can immediately write down the expectation value of the interaction energy in a product state

$$\begin{array}{rl}⟨\hat{V}⟩& =\stackrel{\equiv E_{\text{Hartree}}}{\overbrace{\frac{1}{2}\int d\mathbf{r}d{\mathbf{r}}^{\prime }{\rho }_1(\mathbf{r})U(\mathbf{r}-{\mathbf{r}}^{\prime }){\rho }_1({\mathbf{r}}^{\prime })}}\\ & \stackrel{\equiv E_{\text{Fock}}}{\overbrace{\pm \frac{1}{2}\int d\mathbf{r}d{\mathbf{r}}^{\prime }U(\mathbf{r}-{\mathbf{r}}^{\prime })g(\mathbf{r},{\mathbf{r}}^{\prime })g({\mathbf{r}}^{\prime },\mathbf{r})}}.\end{array}$$

The two terms are known as the Hartree and Fock (or exchange) contributions, respectively. This expression lies at the core of the variational Hartree–Fock method for many body systems, which approximates the ground state by a product state. The Hartree term looks completely reasonable, while the Fock potential doesn’t look like a potential at all, and reflects the non-classical correlations.

2.2 Hartree–Fock for the Electron Gas

How does the Hartree–Fock picture change when we have spin? Let’s consider a system of spin-1/2 fermions. We can describe such a system in terms of field operators ${\psi }_{\sigma }^{\phantom{†}}(\mathbf{r})$, ${\psi }_{\sigma }^{†}({\mathbf{r}}^{\prime })$ satisfying the canonical anticommutation relations

$$\begin{array}{c}\{{\psi }_{{\sigma }_1}^{\phantom{†}}({\mathbf{r}}_1),{\psi }_{{\sigma }_2}^{†}({\mathbf{r}}_2)\}={\delta }_{{\sigma }_1{\sigma }_2}\delta ({\mathbf{r}}_1-{\mathbf{r}}_2)\\ \{{\psi }_{{\sigma }_1}^{\phantom{†}}({\mathbf{r}}_1),{\psi }_{{\sigma }_2}^{\phantom{†}}({\mathbf{r}}_2)\}=\{{\psi }_{{\sigma }_1}^{†}({\mathbf{r}}_1),{\psi }_{{\sigma }_2}^{†}({\mathbf{r}}_2)\}=0.\end{array}$$

The density matrix is a matrix in spin space as well as real space

$$g_{{\sigma }_1{\sigma }_2}({\mathbf{r}}_1,{\mathbf{r}}_2)=⟨\mathrm{\Psi }|{\psi }_{{\sigma }_1}^{†}({\mathbf{r}}_1){\psi }_{{\sigma }_2}^{\phantom{†}}({\mathbf{r}}_2)|\mathrm{\Psi }⟩.$$

From $g_{{\sigma }_1{\sigma }_2}({\mathbf{r}}_1,{\mathbf{r}}_2)$ we can extract the spin density as well the density

$$\rho (\mathbf{r})=\mathrm{tr}[g(\mathbf{r},\mathbf{r})],\mathbf{s}(\mathbf{r})=\frac{1}{2}\mathrm{tr}[\mathit{\sigma }g(\mathbf{r},\mathbf{r})].$$

(I’ve dropped the subscript from ${\rho }_1(\mathbf{r})$ here.) An spin-independent interaction potential is described by a Hamiltonian of the form

$${\hat{H}}_{\text{int.}}=\frac{1}{2}\sum _{{\sigma }_1,{\sigma }_2}\int d{\mathbf{r}}_{1}d{\mathbf{r}}_{2}U({\mathbf{r}}_1-{\mathbf{r}}_2){\psi }_{{\sigma }_1}^{†}({\mathbf{r}}_1){\psi }_{{\sigma }_2}^{†}({\mathbf{r}}_2){\psi }_{{\sigma }_2}^{\phantom{†}}({\mathbf{r}}_2){\psi }_{{\sigma }_1}^{\phantom{†}}({\mathbf{r}}_1).$$

The Hartree–Fock energy is then

$$\begin{array}{rl}⟨{\hat{H}}_{\text{int.}}⟩& =\frac{1}{2}\int d\mathbf{r}d{\mathbf{r}}^{\prime }\rho (\mathbf{r})U(\mathbf{r}-{\mathbf{r}}^{\prime })\rho ({\mathbf{r}}^{\prime })\\ & -\frac{1}{2}\int d\mathbf{r}d{\mathbf{r}}^{\prime }U(\mathbf{r}-{\mathbf{r}}^{\prime })\mathrm{tr}[g(\mathbf{r},{\mathbf{r}}^{\prime })g({\mathbf{r}}^{\prime },\mathbf{r})].\end{array}$$

The Fock term can be rewritten in a more useful way using the identity

$$\begin{array}{}\text{(10)}& {\delta }_{ab}{\delta }_{cd}=\frac{1}{2}[{\mathit{\sigma }}_{ac}\cdot {\mathit{\sigma }}_{db}+{\delta }_{ac}{\delta }_{db}].\end{array}$$

Check

One way to understand Equation 10 is to to think of the two Pauli matrices as acting on two spin 1/2s (as we did in Lecture 4, in which case we can work in the basis $|s_1⟩|s_2⟩$ and we have

$${\mathit{\sigma }}_1\cdot {\mathit{\sigma }}_2+{\mathbb{1}}_1{\mathbb{1}}_2=2\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 0& 1& 0\\ 0& 1& 0& 0\\ 0& 0& 0& 1\end{array}\right).$$ The matrix elements of the right hand side are then $2{\delta }_{s_1,s_2^{\prime }}{\delta }_{s_2,s_1^{\prime }}$.

which gives

$$\begin{array}{rl}{E}_{\text{Fock}}& =-\frac{1}{4}\int d\mathbf{r}d{\mathbf{r}}^{\prime }U(\mathbf{r}-{\mathbf{r}}^{\prime })\mathrm{tr}[g(\mathbf{r},{\mathbf{r}}^{\prime })]\mathrm{tr}[g({\mathbf{r}}^{\prime },\mathbf{r})]\\ & -\frac{1}{4}\int d\mathbf{r}d{\mathbf{r}}^{\prime }U(\mathbf{r}-{\mathbf{r}}^{\prime })\mathrm{tr}[\mathit{\sigma }g(\mathbf{r},{\mathbf{r}}^{\prime })]\cdot \mathrm{tr}[\mathit{\sigma }g({\mathbf{r}}^{\prime },\mathbf{r})].\end{array}$$

Suppose we had a $\delta$-function interaction $U(\mathbf{r})=V_0\delta (\mathbf{r})$. Then the Fock energy can be written

$$\begin{array}{r}{E}_{\text{Fock}}=-\frac{V_0}{4}\int d\mathbf{r}\rho (\mathbf{r}{)}^2-V_0\int d\mathbf{r}\mathbf{s}(\mathbf{r})\cdot \mathbf{s}(\mathbf{r})\end{array}$$

The second term favours ferromagnetism for repulsive interactions. The physical origin is the same as the Hund’s rule coupling in atoms: fermions in different spin states can sit at the same spatial location, while those in the same spin state must be in different locations. For repulsive interactions occupying the same spin state is energetically favourable.

Check

This is most succintly put by the formula

$${\rho }_2(\mathbf{r},\mathbf{r})=\frac{1}{2}\rho (\mathbf{r}{)}^2-2\mathbf{s}(\mathbf{r})\cdot \mathbf{s}(\mathbf{r})$$

The Hartree–Fock energy forms the basis of a variational method using product states as variational wavefunctions. For a Hamiltonian with translational invariance, like

$$\begin{array}{}\text{(11)}& H=\int d\mathbf{r}\frac{1}{2m}\mathrm{\nabla }{\psi }^{†}\cdot \mathrm{\nabla }{\psi }^{\phantom{†}}+\frac{1}{2}\int d\mathbf{r}d{\mathbf{r}}^{\prime }U(\mathbf{r}-{\mathbf{r}}^{\prime }){\psi }^{†}(\mathbf{r}){\psi }^{†}({\mathbf{r}}^{\prime }){\psi }^{\phantom{†}}({\mathbf{r}}^{\prime }){\psi }^{\phantom{†}}(\mathbf{r}),\end{array}$$

this is not too bad, as we are guaranteed to be working with plane wave single particle states. Then the only variational parameters are the occupancies of these states: we’ll meet an example in the next section. If translational symmetry is broken by introducing a potential $U(\mathbf{r})$ into the single particle part of Equation 11, say, the situation is more complicated. We may be tempted to use single particle states that are the eigenstates of the single particle part of the Hamiltonian, but there is no real justification for this, and one has to allow the states, as well as the occupancies, to vary.

2.3 Stoner Criterion for Ferromagnetism

Let us try to put a bit more flesh on the idea that repulsive interactions favour ferromagnetism in fermionic systems. We will continue to use the model interaction $U(\mathbf{r})=V_0\delta (\mathbf{r})$ that we introduced in the previous section. Of course, this isn’t a realistic interaction between electrons in a metal, say, but as we’ll see in Lecture 11, the long-ranged Coulomb interaction is screened and becomes finite-ranged. Thus our model is not a bad approximation to the effective interaction between electrons in a metal.

Polarizing the spins in a Fermi gas is not without cost (otherwise everything would be ferromagnetic!): there is a price to pay in increased kinetic energy. To understand why this is so, consider the ground state kinetic energy of $N$ (spinless) fermions in three dimensions, obtained from

$$\begin{array}{rl}N=\sum _{|\mathbf{k}|<k_{\text{F}}}1& ⟶L^3{\int }_{|\mathbf{k}|<k_{\text{F}}}\frac{d\mathbf{k}}{(2\pi {)}^3}=\frac{k_{\text{F}}^3{L}^3}{6{\pi }^2}\\ E_{\text{kin}}=\sum _{|\mathbf{k}|<k_{\text{F}}}\frac{{\mathbf{k}}^2}{2m}& ⟶L^3{\int }_{|\mathbf{k}|<k_{\text{F}}}\frac{d\mathbf{k}}{(2\pi {)}^3}\frac{{\mathbf{k}}^2}{2m}\\ & =\frac{k_{\text{F}}^5{L}^3}{20{\pi }^{2}m}=L^3\frac{3}{10m}(6{\pi }^2{)}^{2/3}{n}^{5/3},\end{array}$$

(The assumption of a quadratic dispersion is not important here. More generally, we fill a band structure.) where $n=\frac{N}{L^3}$ is the mean density. Assuming a system of fermions with spin is now polarized in the $z$-direction, we have differing densities $n_{\uparrow ,\downarrow }$ of spin up and spin down fermions. Their total energy is

$$E_{\text{kin}}(n_{\uparrow },n_{\downarrow })=\frac{c{L}^3}{m}(n_{\uparrow }^{5/3}+n_{\downarrow }^{5/3}),$$

where $c=\frac{3}{10}(6{\pi }^2{)}^{2/3}$. In terms of the overall density $n=n_{\uparrow }+n_{\downarrow }$ and density of spin $\overline{s}=(n_{\uparrow }-n_{\downarrow })/2$, we have

$$E_{\text{kin}}(n,\overline{s})=\frac{c{L}^3}{m}({[n/2+\overline{s}]}^{5/3}+{[n/2-\overline{s}]}^{5/3}).$$

Alternatively, write this in terms of the polarization $P\equiv \frac{n_{\uparrow }-n_{\downarrow }}{n}$ that varies in the range $[-1,1]$ as

$$E_{\text{kin}}(P)=\frac{E_{\text{K}}}{2}[(1+P{)}^{5/3}+(1-P{)}^{5/3}].$$

We see that, on account of the convexity of $x^{5/3}$, $E_{\text{kin}}^{(0)}(n,\overline{s})$ is minimized for $s=0$.

Let’s compare this with the effect of interactions. In the short-ranged model introduced above, the total Hartree–Fock energy is

$$E_{\text{HF}}(n,\overline{s})=\frac{V_0{L}^3}{2}{n}^2-\frac{V_0{L}^3}{2}(n_{\uparrow }^2+n_{\downarrow }^2)=\frac{V_0{L}^3}{2}(\frac{1}{2}{n}^2-2{\overline{s}}^2).$$

We write this in terms of the polarization as

$$E_{\text{HF}}(P)=\frac{E_V}{2}(1-P^2).$$

Check

Minimize the total energy $E(P)=E_{\text{kin}}(P)+E_{\text{HF}}(P)$ to show

1. For $E_V/E_K<10/9$ the ground state is non-magnetic.
2. As $E_V/E_K$ increases past $10/9$ the magnetization begins to increase.
3. At $E_V/E_K>\frac{5}{6}{2}^{2/3}$ is the ground state is fully polarized.

One shouldn’t take these numerical values too seriously given the simplicity of the model, but they illustrate the physical principles at work behind the appearance of ferromagnetism in metals.

2.4 Excited State Energies

So far we have discussed properties of the ground state only. We can, however, evaluate the Hartree–Fock energy in a product state describing an excited state of a noninteracting system. This gives the first order perturbation theory correction to the excited state energy. For concreteness we will stick with fermions for now, though the method is general.

If we work in a translationally invariant system, the appropriate single particle states are plane waves. We write the field operators in the plane wave basis as

$$\begin{array}{r}{\psi }^{\phantom{†}}(\mathbf{r})\equiv \frac{1}{L^{3/2}}\sum _{\mathbf{k}}\exp (i\mathbf{k}\cdot \mathbf{r})a_{\mathbf{k}}^{\phantom{†}},\\ {\psi }^{†}(\mathbf{r})\equiv \frac{1}{L^{3/2}}\sum _{\mathbf{k}}\exp (-i\mathbf{k}\cdot \mathbf{r})a_{\mathbf{k}}^{†},\end{array}$$

and represent the interaction potential in terms of its Fourier components

$$U(\mathbf{r}-{\mathbf{r}}^{\prime })=\frac{1}{L^3}\sum _{\mathbf{q}}\tilde{U}(\mathbf{q})\exp (i\mathbf{q}\cdot [\mathbf{r}-{\mathbf{r}}^{\prime }]).$$

The interaction Hamiltonian for spinless particles can then be written

$$\begin{array}{}\text{(12)}& {\hat{H}}_{\text{int.}}=\frac{1}{2L^3}\sum _{{\mathbf{k}}_1+{\mathbf{k}}_2={\mathbf{k}}_3+{\mathbf{k}}_4}\tilde{U}({\mathbf{k}}_1-{\mathbf{k}}_4)a_{{\mathbf{k}}_1}^{†}{a}_{{\mathbf{k}}_2}^{†}{a}_{{\mathbf{k}}_3}^{\phantom{†}}{a}_{{\mathbf{k}}_4}^{\phantom{†}}.\end{array}$$

When written in this way, interaction Hamiltonians are sometimes associated with the graphical representation below, one of the ingredients of the Feynman diagram technique for performing perturbation theory calculations in field theories. The incoming lines (arrows in) represent particles being removed (in momentum states ${\mathbf{k}}_3$ and ${\mathbf{k}}_4$) and the outgoing lines represent particles added (momenta ${\mathbf{k}}_1$ and ${\mathbf{k}}_2$). The wiggly line represents $\tilde{V}(\mathbf{q})$. The conservation of momentum at the vertices is a consequence of the translational invariance of the problem: the two integrations over $\mathbf{r}$ and ${\mathbf{r}}^{\prime }$ give rise to two $\delta$-functions that perform this function.

Graphical representation of the interaction Equation 12.

As we’ve already discussed, the expectation value of Equation 12 in a product state of momentum eigenstates gives two terms, with different ‘pairings’ of creation operators with annihilation operators. We can represent these two terms graphically as shown below.

Graphical representation of the Hartree and Fock terms.

Evaluating the two contributions in terms of the occupation numbers gives

$$⟨\mathbf{N}|{\hat{H}}_{\text{int.}}|\mathbf{N}⟩=\frac{1}{2V}\tilde{U}(0)\sum _{{\mathbf{k}}_1,{\mathbf{k}}_2}{N}_{{\mathbf{k}}_1}{N}_{{\mathbf{k}}_2}-\frac{1}{2V}\sum _{{\mathbf{k}}_1,{\mathbf{k}}_2}\tilde{U}({\mathbf{k}}_1-{\mathbf{k}}_2)N_{{\mathbf{k}}_1}{N}_{{\mathbf{k}}_2}$$

While the Hartree term just depends on the total number of particles, the Fock term depends on the individual occupations. The interaction energy to add a single particle to state $\mathbf{k}$ is

$$\mathrm{\Delta }{U}_{\mathbf{k}}=\frac{\tilde{U}(0)}{V}\sum _{{\mathbf{k}}^{\prime }}{N}_{{\mathbf{k}}^{\prime }}-\frac{1}{V}\sum _{{\mathbf{k}}^{\prime }}\tilde{U}(\mathbf{k}-{\mathbf{k}}^{\prime })N_{{\mathbf{k}}^{\prime }}$$

References

Andrews, MR, CG Townsend, H-J Miesner, DS Durfee, DM Kurn, and W Ketterle. 1997. “Observation of Interference Between Two Bose Condensates.” Science 275 (5300): 637–41.

Baym, Gordon. 1999. “The Physics of Hanbury Brown-Twiss Interferometry: From Stars to Nuclear Collisions to Atomic and Condensed Matter Systems.” Measuring The Size Of Things In The Universe: Hbt Interferometry And Heavy Ion Physics: Proceedings Of Cris’ 98, 11.

———. 2018. Lectures on Quantum Mechanics. CRC Press.

Kleppner, Daniel. 2008. “Hanbury Brown’s Steamroller.” Physics Today 61 (8): 8–9.

Magyar, G, and L Mandel. 1963. “Interference Fringes Produced by Superposition of Two Independent Maser Light Beams.” Nature 198 (4877): 255–56.

Nazarov, Yuli V, and Jeroen Danon. 2013. Advanced Quantum Mechanics: A Practical Guide. Cambridge University Press.