Week 3: Many-body quantum mechanics continued

Learning goals for this week:

- Learning the idea of one- and two-particle Fock-space operators 
- Understanding the concept of the Fermi sea and learning how to compute some of its properties
- Particle-hole formalism
- Start of discussing an interacting system

One- and two-particle Fock-space operators

One-particle operators

A general 1-particle operator $\hat F^{(1)}$ is a sum of identical operators, each of which operate on one-particle states \[ \hat F^{(1)} = \sum_{j=1}^N \hat f_j, \] where $\hat f_j$ acts solely on the $j$th particle. An example of such a 1-particle operator is the kinetic energy operator $\hat T = \sum_j \frac{{\hat{\vec p}}_j^2}{2m}$ which is diagonal in spin.

We should now find the corresponding operator in the $N$-particle Fock space.

Let us write the operators $\hat f_j$ in a diagonal basis so that $\hat f_j = \sum_{\nu_1,\ldots,\nu_N} f_{\nu_j}\ket{\nu_1} \cdots \ket{\nu_N}\bra{\nu_1}\cdots\bra{\nu_N}$ where $f_{\nu} = \bra{\nu}\hat f \ket{\nu}$ is the expectation value of the single-particle operator. Next, we want to consider the action of $F^{(1)}$ to a $N$-particle Fock state $\ket{n_{\nu_1}, n_{\nu_2} \dots}$ with $\sum_j n_{\nu_j} = N$. Therefore we have to go back to description of the $N$-particle Fock state in terms of the single-particle states. For example, take a bosonic state \[ \begin{aligned} \hat F^{(1)} \ket{n_{\nu_1}, n_{\nu_2} \dots} &= \hat F^{(1)} \frac{1}{\mathcal{N}} \sum_P \ket{\nu_{i_1}} \cdots \ket{\nu_{i_N}} \\ &= \sum_{j=1}^N \frac{1}{\mathcal{N}} \sum_P \ket{\nu_{i_1}} \cdots \underbrace{\hat f\ket{\nu_{i_j}}}_{f_{\nu_{i_j}}} \dots \ket{\nu_{i_N}} \\ &= \frac{1}{\mathcal{N}} \sum_P (\sum_{j=1}^N f_{\nu_{i_j}}) \ket{\nu_{i_1}} \dots \ket{\nu_{i_N}} \\ &= (\sum_i f_{\nu_i} \underbrace{n_{\nu_i}}_{\hat n_{\nu_i}}) \ket{n_{\nu_1}, n_{\nu_2} \dots} \end{aligned} \] In the second to last step, we exchanged the order of summation, since the sum over $j$ is the same for any permutation (as we sum over all indices). Then, we can use the definition of the number operator, and since the above relation holds for any state (derivation is the same for a fermionic state), we have \[ F^{(1)} = \sum_\nu f_\nu \hat n_\nu \] More generally, we can work in a non-diagonal basis by using $b_\nu = \sum_\mu \langle{\nu}|{\mu}\rangle b_\mu$ and $\hat f = \sum_\nu f_\nu \ket{\nu}\bra{\nu}$ which gives 1-particle operator in the Fock space \[ \hat F^{(1)} = \sum_{\mu\mu'} \bra{\mu}\hat f\ket{\mu'} b_\mu^\dagger b_{\mu'} \] For fermions, simply change bosonic $b_\nu$ to fermionic $a_\nu$.

In terms of the field operators $\Psi(\vec x, \sigma)$ and $\Psi^\dagger(\vec x, \sigma)$ , we can write the 1-particle operator (assuming local operators $\hat f_i$ ) for fermionic system as $\begin{aligned} \hat F^{(1)} &= \sum_{\mu \nu} \bra{\mu}\hat f\ket{\nu} a_\mu^\dagger a_\nu \\ &= \sum_{\sigma \sigma'}\int d^3 \vec x d^3 \vec x' \sum_{\mu\nu} \braket{\mu}{\vec x, \sigma} \underbrace{\mel{\vec x, \sigma}{\hat f}{\vec x', \sigma'}}_{f_{\sigma \sigma'}(\vec x)\delta^{(3)}(\vec x - \vec x')} \underbrace{\braket{\vec x', \sigma'}{\nu}}_{\varphi_\nu(\vec x', \sigma')} a_\mu^\dagger a_\nu \\ &= \sum_{\sigma \sigma'}\int d^3 \vec x \qty(\sum_\mu \varphi^*_\mu(\vec x, \sigma) a_\mu^\dagger) f_{\sigma \sigma'}(\vec x) \qty(\sum_\nu \varphi_\nu(\vec x, \sigma') a_\nu) \\ &= \sum_{\sigma \sigma'}\int d^3 \vec x \Psi^\dagger(\vec x, \sigma) f_{\sigma \sigma'}(\vec x) \Psi(\vec x, \sigma'). \end{aligned}$ An example of a local operator is the kinetic energy operator for which $f_{\sigma \sigma'}(\vec x) = \delta_{\sigma \sigma'}\qty(-\frac{\hbar^2 \grad^2}{2m})$ .

Action of one-particle operator on a Fock space state \[ \begin{aligned} \hat F^{(1)}\ket{1_{\nu_1}, \dots 1_{\nu_N}} &= \hat F^{(1)} a_{\nu_1}^\dagger \dots a_{\nu_N}^\dagger\ket{0} = \sum_{\mu,\nu} \bra{\mu}\hat f\ket{\nu} a_\mu^\dagger a_\nu a_{\nu_1}^\dagger \dots a_{\nu_N}^\dagger\ket{0} \\ \overset{\text{Wick's th.}}{=}& \sum_\mu \sum_{i=1}^N (-1)^{i-1} \bra{\mu}\hat f\ket{\nu} a_\mu^\dagger a_{\nu_1}^\dagger \dots \underbrace{a_{v_{i-1}}^\dagger a_{v_{i+1}}^\dagger}_{a_{v_{i}}^\dagger \text{missing}} \dots a_{\nu_N}^\dagger\ket{0} \\ &= \sum_\mu \sum_{i=1}^N \bra{\mu}\hat f\ket{\nu_i} a_{\nu_1}^\dagger \dots a_{v_{i-1}}^\dagger a_\mu^\dagger a_{v_{i+1}}^\dagger \dots a_{\nu_N}^\dagger\ket{0}\\ &= \sum_\mu \sum_{i=1}^N \bra{\mu}\hat f\ket{\nu_i} \ket{1_{\nu_1}, \dots 1_{\nu_{i-1}}, 1_{\mu}, 1_{\nu_{i+1}} \dots 1_{\nu_N}}. \end{aligned} \]

Alternative derivation for $\hat F^{(1)}$

Alternatively, we can write $\hat f \ket{\nu}$ in the position basis as $\bra{\vec x, \sigma}\hat f\ket{\nu} = \sum_\mu \braket{\vec x,\sigma}{\mu} \bra{\mu}{\hat f}\ket{\nu} = \sum_\mu \varphi_\mu(\vec x,\sigma)\bra{\mu}{\hat f}\ket{\nu}.$ Next, we assume, for brevity, that $\hat f$ is both local and diagonal in spin basis. That is, $\bra{\vec x,\sigma}\hat f \ket{\vec x,\sigma} = f(\vec x)\delta^{(3)}(\vec x - \vec x')\delta_{\sigma \sigma'}$ . Then, from the above equation by inserting the completeness relation $\hat 1 = \sum_{\sigma'} \int d^3 x' \dyad{\vec x',\sigma'}$ we obtain $f(\vec x) \varphi_\nu(\vec x, \sigma) = \sum_\mu \varphi_\mu(\vec x,\sigma)\bra{\mu}{\hat f}\ket{\nu}.$ Recalling that the fermionic WF is the Slater determinant $\psi_{\nu_1 \dots \nu_N}(1, \dots N) = \frac{1}{\sqrt{N!}} \sum_{P \in S_{N}} \epsilon_{P} \varphi_{\nu_{1}}(i_{1}) \varphi_{\nu_{2}}(i_{2}) \cdots \varphi_{\nu_{N}}(i_N) = \frac{1}{\sqrt{N!}}\sum_P \epsilon_P \prod_{j=1}^N \varphi_{\nu_j}(\underbrace{P(j)}_{i_j(P)})$ we can determine the effect of $\hat F^{(1)}$ on $\psi_{\nu_1 \dots \nu_N}(1, \dots N)$ : $F^{(1)} \psi_{\nu_1 \dots \nu_N}(1, \dots N) =\sum_{i=1}^N f_i \frac{1}{\sqrt{N!}}\sum_P \epsilon_P \prod_{j=1}^N \varphi_{\nu_j}(P(j)).$ Next, we rewrite the sums so that we sum over for each , and use the result above for $f(\vec x) \varphi_\nu(\vec x, \sigma)$ : $\begin{aligned} F^{(1)} \psi_{\nu_1 \dots \nu_N}(1, \dots N) &= \frac{1}{\sqrt{N!}}\sum_P \epsilon_P \sum_{i=1}^N f(P(i)) \prod_{j=1}^N \varphi_{\nu_j}(P(j)) \\ &= \sum_\mu \sum_{i=1}^N \bra{\mu}{\hat f}\ket{\nu_i} \underbrace{\frac{1}{\sqrt{N!}}\sum_P \epsilon_P \varphi_\mu(P(i))\prod_{j \neq i} \varphi_{\nu_j}(P(j))}_{\text{Slater det. from } i \rightarrow \mu}. \end{aligned}$ This matches the earlier result: a particle is first removed from $\ket{\nu_i}$ and then a new one is created in $\ket{\mu}$ .

Control question on calculating the one-particle energy

# ConsidertheoneparticleHamiltonianbrbrhatH0sumnEnhatandaggerhatanbrbrCalculatetheenergyofthestatepsiranglea1daggerbra2dagger0rangleieElanglepsihatH0psibrrangleTheresultis

Two-particle operators

2-particle operators $\hat F^{(2)}$ describe pairwise interactions between particles. In many cases, these pairwise interactions are described by a pair-potential, e.g. Coulomb interaction $V(\vec x, \vec x') \propto |\vec x - \vec x'|^{-1}$ between electrons. Thus, it is convenient to work in the position-basis, i.e. with the field operators $\Psi(\vec x)$. For brevity, we neglect spin in this section.

Consider a symmetric 2-particle interaction $V(\vec x, \vec x') = V(\vec x', \vec x)$. Now, we want to find $\hat V$ in the second quantized form so that it matches the position-basis result \[ \hat V \ket{\vec x_1, \vec x_2 \dots \vec x_N} = \sum_{i > j} V(\vec x_i,\vec x_j) \ket{\vec x_1, \vec x_2 \dots \vec x_N} = \frac{1}{2}\sum_{i\neq j} V(\vec x_i,\vec x_j) \ket{\vec x_1, \vec x_2 \dots \vec x_N}. \] Note how the sums go over all the possible interaction pairs once. Since we know the form of 1-particle operator, we can guess a possible generalization \[ \hat V = K \int d^3 \vec x \int d^3 \vec x' \Psi^\dagger(\vec x) \Psi^\dagger(\vec x') V(\vec x,\vec x') \Psi(\vec x) \Psi(\vec x'). \] For example, an equally plausible integrand $\Psi^\dagger(\vec x)\Psi(\vec x) V(\vec x,\vec x')\Psi^\dagger(\vec x') \Psi(\vec x')$ does not work, so the guess must be fairly good. In any case, the prefactor $K$ is still to be determined.

We can now apply the ansatz to the above equation. For a fermionic system, we have \[ \begin{aligned} \Psi^\dagger(\vec x) \Psi^\dagger(\vec x') \Psi(\vec x) \Psi(\vec x') \ket{\vec x_1, \vec x_2 \dots \vec x_N} &= \Psi^\dagger(\vec x) \Psi^\dagger(\vec x') \Psi(\vec x) \Psi(\vec x') \Psi^\dagger(\vec x_1) \dots \Psi^\dagger(\vec x_N)\ket{0}\\ &= -\sum_{i\neq j} \delta^{(3)}(\vec x - \vec x_i) \delta^{(3)}(\vec x' - \vec x_j) \ket{\vec x_1, \vec x_2 \dots \vec x_N} \end{aligned} \] by the use of the anticommutation relations. By integrating over $\vec x$ and $\vec x'$ we find that $K = - \frac{1}{2}$. The minus sign we may absorb to the ordering of the field operators giving in the end (including spin) \[ \hat F^{(2)} = \frac{1}{2}\sum_{\sigma_1\sigma_2\sigma_3\sigma_4} \int d^3\vec x d^3\vec x' \Psi^\dagger(\vec x,\sigma_1) \Psi^\dagger(\vec x',\sigma_2) V_{\sigma_1\sigma_2\sigma_3\sigma_4}(\vec x,\vec x') \Psi(\vec x',\sigma_3) \Psi(\vec x,\sigma_4) \]

This interaction can be illustrated with the Feynman diagram

$\begin{tikzpicture}[node distance=1cm and 1.5cm] \coordinate[label=left:${(\vec x,\sigma_4)}$] (e1); \coordinate[below right=of e1] (aux1); \coordinate[above right=of aux1,label=right:${(\vec x, \sigma_1)}$] (e2); \coordinate[below=1.25cm of aux1] (aux2); \coordinate[below left=of aux2,label=left:${(\vec x',\sigma_3)}$] (e3); \coordinate[below right=of aux2,label=right:${(\vec x',\sigma_2)}$] (e4); \draw[postaction={decorate}, decoration={markings,mark=at position .5 with {\arrow{triangle 45}}}] (e1) -- (aux1); \draw[postaction={decorate}, decoration={markings,mark=at position .5 with {\arrow{triangle 45}}}] (aux1) -- (e2); \draw[postaction={decorate}, decoration={markings,mark=at position .5 with {\arrow{triangle 45}}}] (e3) -- (aux2); \draw[postaction={decorate}, decoration={markings,mark=at position .5 with {\arrow{triangle 45}}}] (aux2) -- (e4); \draw[decorate, draw=black, decoration={coil,aspect=0}] (aux1) -- node[label=right:${V_{\sigma_1\sigma_2\sigma_3\sigma_4}(\vec x, \vec x')}$] {} (aux2); \end{tikzpicture}$

in which the annihilation and creation operators are represented by the incoming and outgoing lines, respectively. Solid lines are fermions, and the interaction is represented by the wiggly line connecting the two electron lines. Feynman diagrams will become more familiar (and they will be given a more rigorous meaning) in the quantum field theory course, but for now, we just use them to give an intuitive description for the interaction operators.

As before, we can then move to a different basis to obtain a more general form 2-particle operator in the Fock space \[ \hat F^{(2)} = \sum_{\mu\mu'\nu\nu'} v_{\mu \mu' \nu \nu'} a^\dagger_\mu a^\dagger_{\mu'} a_{\nu'}a_{\nu} \] where $v_{\mu \mu' \nu \nu'} = \frac{1}{2}\langle{\mu \mu'}|\hat V |\nu \nu'\rangle$. Generally, $v_{\mu \mu' \nu \nu'} = v_{\mu' \mu \nu' \nu}$. This result holds for bosonic system similarly.

# two-body-interaction-in-a-translational-invariant-system

Two-body interaction in a translational invariant system

We concentrate especially on two-body interactions that depend only on the distance between the fermions. In this case \[V_{\sigma_1\sigma_2\sigma_3\sigma_4}({\vec{r}},{\vec{r}}') \equiv V_2({\vec{r}}-{\vec{r}}')\delta_{\sigma_1 \sigma_4}\delta_{\sigma_2\sigma_3}.\] Then \[\begin{aligned} &H^{(2)}=\frac{1}{2} \sum_{\sigma \sigma'} \int d^{3} {\vec{r}}_1 d^{3} {\vec{r}}_2 \Psi^\dagger({\vec{r}}_1,\sigma) \Psi^\dagger({\vec{r}}_2,\sigma') V_2({\vec{r}}_1-{\vec{r}}_2) \Psi({\vec{r}}_2,\sigma') \Psi({\vec{r}}_1,\sigma)\\ &=\frac{1}{2} \sum_{\sigma\sigma'} \sum_{nmkl} \int d^{3} {\vec{r}}_1 d^{3} {\vec{r}}_2 \varphi_n^*({\vec{r}}_1) \varphi_m^*({\vec{r}}_2) V_2({\vec{r}}_1-{\vec{r}}_2) \varphi_l({\vec{r}}_2)\varphi_k({\vec{r}}_1) a_{n\sigma}^\dagger a_{m\sigma'}^\dagger a_{l\sigma'} a_{k\sigma},\end{aligned}\] where the second line is a transformation from the field operators to those specified by the quantum numbers $n,m,l,k$. In particular, for definite-momentum states the wave functions are plane waves, \[\varphi_n({\vec{r}})=\frac{1}{\sqrt{V}} e^{i{\vec{k}}_n \cdot {\vec{r}}}.\] Now let us define $\Delta {\vec{r}}= {\vec{r}}_1-{\vec{r}}_2$ and ${\vec{R}}=({\vec{r}}_1+{\vec{r}}_2)/2$. We can write the above Hamiltonian as \[\begin{aligned} H^{(2)} = \frac{1}{2V^2} \sum_{\sigma \sigma'} &\sum_{nmlk} \int d^{3} \Delta {\vec{r}}e^{-i({\vec{k}}_n-{\vec{k}}_m-{\vec{k}}_k+{\vec{k}}_l) \cdot \Delta {\vec{r}}} V_2(\Delta {\vec{r}})\\& \times \underbrace{\int d^{3} {\vec{R}} e^{-i({\vec{k}}_n+{\vec{k}}_m-{\vec{k}}_k-{\vec{k}}_l) \cdot {\vec{R}}}}_{V\delta({\vec{k}}_n+{\vec{k}}_m-{\vec{k}}_k-{\vec{k}}_l)} a_{n\sigma}^\dagger a_{m\sigma'}^\dagger a_{l\sigma'} a_{k\sigma}.\end{aligned}\] Now define ${\vec{k}}_k \equiv {\vec{k}}$, ${\vec{k}}_l \equiv {\vec{k}}'$, and $\vec q=({\vec{k}}_n-{\vec{k}}_m-{\vec{k}}_k+{\vec{k}}_l)/2=k_n-k_k=-k_m+k_l$ (for those momenta coupled by the $\delta$-function in the above equation). Moreover, the Fourier transformation of the two-body potential is \[V_2(\vec q)= \int d^{3} {\vec{r}}e^{-i \vec q \cdot {\vec{r}}} V_2({\vec{r}}).\] With these, we can write the most general form of the two-body Hamiltonian in the translationally invariant case as \[H^{(2)} = \frac{1}{2V} \sum_{\sigma \sigma'} \sum_{{\vec{k}},{\vec{k}}',\vec q} V_2(\vec q) a_{{\vec{k}}+\vec{q},\sigma}^\dagger a_{{\vec{k}}'-\vec q,\sigma'}^\dagger a_{{\vec{k}}' \sigma'} a_{{\vec{k}}\sigma}.\] This interaction can be represented with the following Feynman diagram, in which the spin of the individual particles is conserved. Momentum conservation in the diagram shows up so that at every vertex the sum of incoming momenta equals the sum of outgoing momenta.

$\begin{tikzpicture}[node distance=1cm and 1.5cm] \coordinate[label=left:${(\vec k, \sigma)}$] (e1); \coordinate[below right=of e1] (aux1); \coordinate[above right=of aux1,label=right:${(\vec k{+}\vec q,\sigma)}$] (e2); \coordinate[below=1.25cm of aux1] (aux2); \coordinate[below left=of aux2,label=left:${(\vec k',\sigma')}$] (e3); \coordinate[below right=of aux2,label=right:${(\vec k'{-}\vec q,\sigma')}$] (e4); \draw[postaction={decorate}, decoration={markings,mark=at position .5 with {\arrow{triangle 45}}}] (e1) -- (aux1); \draw[postaction={decorate}, decoration={markings,mark=at position .5 with {\arrow{triangle 45}}}] (aux1) -- (e2); \draw[postaction={decorate}, decoration={markings,mark=at position .5 with {\arrow{triangle 45}}}] (e3) -- (aux2); \draw[postaction={decorate}, decoration={markings,mark=at position .5 with {\arrow{triangle 45}}}] (aux2) -- (e4); \draw[decorate, draw=black, decoration={coil,aspect=0}] (aux1) -- node[label=right:$\mkern-65muV_2(\vec q)\mkern30mu\vec q$] {} (aux2); \draw[->] (1.7,-2.0) -- (1.7,-1.2); \end{tikzpicture}$

Non-interacting Fermi gas

Let us consider a gas of nonrelativistic, noninteracting spin- $\frac 1 2$ fermions, a Fermi gas, for which $\hat H = \sum_{i_1}^N \hat H_i^{(1)},\quad\text{with}\quad \hat H_i^{(1)} = \frac{(\hat{\vec p}_i)^2}{2m}.$ I.e. for 1-particle states and wavefunctions $\hat H^{(1)} \psi_{\vec p\sigma}(\vec x) = \varepsilon_{p} \psi_{\vec p\sigma}(\vec x),\quad\text{with}\quad \varepsilon_{p} = \frac{p^2}{2m}.$

In order to have wavefunctions normalized to 1, we put the particles in a box of volume V=L^3 , and require periodic boundary conditions: $\psi(x,y,z,\sigma) = \psi(x+L,y,z,\sigma) = \psi(x,y+L,z,\sigma) = \psi(x,y,z+L,\sigma)$ In a box , the allowed values of momentum are $\vec p = \frac{2\pi\hbar}{L}(n_x,n_y,n_z)$ , for integer n_i .

The Hamiltonian, written in terms of creation and annihilation operators, is $\hat H = \sum_{\vec p,\sigma} \varepsilon_{p} a^\dagger_{\vec p\sigma} a_{\vec p\sigma} = \sum_{\vec p,\sigma} \varepsilon_{p} \hat n_{\vec p\sigma},$ where $\vec p$ is summed over all allowed values, as determined from the boundary conditions and spin is summed over the set $\{\uparrow,\downarrow\}$ .

The total number of particles is given by the expectation value of $\hat N = \sum_{\vec p,\sigma} \hat n_{\vec p\sigma}.$

We observe that $[\hat H, \hat N]=0$ , so the total number of particles is conserved in the noninteracting Fermi gas. This is understandable as there are no interactions in the system which could change the number of particles.

The next task is find the ground state of the Fermi gas. If the particles were bosons, they would all occupy the lowest energy state in the ground state. Now, for the case of identical fermions, there can be at most one particle occupying a state $(\vec p, \sigma)$, i.e. two fermions for each $\vec p$.

Thus, for fermions, the ground state is such that all states with energies below some highest energy $\varepsilon_F=\frac{p_F^2}{2m}$ are occupied by the $N$ particles in the system. $\varepsilon_F$ is the Fermi energy. In the momentum space the occupied states form a sphere of radius $p_F$, where $p_F$ is known as the Fermi momentum.

The $N$-particle ground state can be expressed as \[ |{\rm g}\rangle = \prod_{|\vec p|\leq p_F,\sigma} a^\dagger_{\vec p\sigma} |0\rangle. \]

$p_F$ and $\varepsilon_F$ are fixed by the number of particles. To determine them, we apply the total number operator to the ground state. \[ \hat N |{\rm g}\rangle = \sum_{\vec p,\sigma} \hat n_{\vec p\sigma}|{\rm g}\rangle = 2 \sum_{|\vec p|\leq p_F} |{\rm g}\rangle. \] At the thermodynamic limit $V\to\infty$ we find a particle density \[ n = \frac{N}{V}=2\lim_{V\to\infty}\sum_{|\vec p|\leq p_F} 1 = 2 \int_{p\leq p_F} d\vec p \frac{1}{(2\pi\hbar)^3} = 2\cdot \frac{1}{(2\pi \hbar)^3} \cdot \frac{4}{3}\pi p_F^3, \] where the integral was evaluated by noticing that it describes a three dimensional sphere in momentum space. The surface of the sphere (the states $p=p_F$) is the Fermi surface.

$\begin{tikzpicture} \begin{axis}[ axis x line=middle, axis y line=middle, xlabel = $p_x/p_F$, ylabel = {$\varepsilon_p/\varepsilon_F$}, ymin=-0.2, ymax=2.5 ] %Here the blue parabloa is defined \addplot [ domain=-1:1, samples=100, color=blue, style={ultra thick} ] {x^2}; \addplot [ domain=-1.5:1.5, samples=100, color=blue, ] {x^2}; \addplot [ domain=-1.8:1.8, samples=100, color=red, style=dashed ] {1}; \end{axis} \end{tikzpicture}$

Fig: Dispersion relation of free fermions (blue line). Fermi energy is marked with the red dashed line. The occupied states for -particle ground state are marked with the bold blue line.

$\begin{tikzpicture} \shade[ball color = blue!80, opacity = 0.3] (0,0) circle (2cm); \draw (0,0) circle (2cm); \draw (-2,0) arc (180:360:2 and 0.6); \draw[dashed] (2,0) arc (0:180:2 and 0.6); \fill[fill=black] (0,0) circle (1pt); \draw[<->] (0.1,0.1) -- node[above]{$p_F$} (1.9,0.1); \draw[->] (0,0) -- (0,2.5) node[anchor=north west] {$p_z$}; \draw[->] (0,0) -- (2.5,0) node[anchor=west] {$p_y$}; \draw[->] (0,0) -- (-1.0,-1.0) node[anchor=north west] {$p_x$}; \end{tikzpicture}$

Fig: Occupied states (blue) of the Fermi gas form a sphere in in the momentum space.

Let us also calculate the ground state energy. It is \[ H|g\rangle = \sum_{\vec p, \sigma} \epsilon_{\vec p} n_{\vec p \sigma} |g\rangle = 2\sum_{|\vec p|<p_F} \epsilon_{\vec p} |g\rangle = E_0 |g\rangle. \] In other words, the ground state energy is \[ E_0 = 2\sum_{|\vec p|<p_F} \epsilon_{\vec p}=\frac{2}{2m} \sum_{|\vec p|<p_F} \vec p^2 \overset{V \rightarrow \infty}{=} \frac{1}{m} \frac{V}{(2\pi \hbar)^3} \int_{p\le p_F} d^3 p p^2 = \frac{1}{m} \frac{V}{(2\pi \hbar)^3} 4 \pi \frac{p_F^5}{5}. \] This can also be expressed in terms of the number of particles $N=n V$ by using the relation $p_F = \hbar (3 \pi^2 n)^{1/3}$ between the Fermi momentum and the particle density $n$: \[ E_0 = \frac{V p_F^5}{10\pi^2 \hbar^3 m} = N \frac{3 p_F^2}{10 m} = \frac{3}{5} N \epsilon_F. \] This result is used below.

Control question: average particle energy in a 1D Fermi gas

# ConsideramodelofaonedimensionalsystemwithnoninteractingfermionsCalculatetheaverageparticleenergyTheresultis

Some properties of the Fermi gas

Single-particle correlation function

The single-particle correlation function of the Fermi gas is defined by \[ G_\sigma(\vec r - \vec r') \equiv \frac{2}{n} \langle g | \Psi^\dagger(\vec r,\sigma) \Psi(\vec r',\sigma) |g\rangle. \] This correlation function describes the probability amplitude that annihilating a particle at $\vec r$ and creating another one at $\vec r'$ (in the same spin state) yields the Fermi gas. Alternatively, it can be regarded as the probability amplitude of the transition from state $\Psi(\vec r',\sigma)|g\rangle$ to the state $\Psi(\vec r,\sigma)|g\rangle$. The prefactor $2/n$ (where $n$ is the particle density) makes sure that $G_\sigma(0)=1$. This can be seen by a direct calculation (exercise) or identifying $\Psi^\dagger(\vec r,\sigma) \Psi(\vec r,\sigma)$ as the operating yielding the local density for particles with spin $\sigma$.

Calculating the single-particle correlation function

Let us calculate it: \[ \begin{aligned} G_{\sigma}\left(\vec r-\vec r^{\prime}\right) &=\frac{2}{n} \sum_{\vec{p}, s;\vec p',s';|\vec p|,|\vec p'|<p_F} \underbrace{\varphi^*_{\vec p,s}(\vec r,\sigma)}_{\frac{1}{\sqrt{V}} e^{-i\hbar \vec p \cdot \vec r} \delta_{s\sigma}} \varphi_{\vec p',s'}(\vec r',\sigma) \underbrace{\bra{g} a_{\vec p,s}^\dagger a_{\vec p',s'} \ket{g}}_{\delta_{\vec p,\vec p'} \delta_{ss'}}\\ &=\frac{2}{n} \frac{1}{V} \sum_{|\vec p| < p_F} e^{-i\vec p \cdot (\vec r-\vec r')/\hbar} \overset{V \rightarrow \infty} = \frac{2}{n} \int_{p < p_F} \frac{d^3 p}{(2\pi\hbar)^3} e^{-i \vec p \cdot (\vec r-\vec r')/\hbar} \\ &=\frac{2}{n} \frac{1}{(2\pi\hbar)^3}\int_{p < p_F} dp p^2 \int_0^{2\pi} d\phi \int_{-1}^1 d(\cos\theta) e^{-i p (|\vec r-\vec r'|) \cos \theta/\hbar}. \end{aligned} \] Let us denote $p=\hbar k$ and $r=|\vec r-\vec r'|$. The integral is then \[ \begin{aligned} G_{\sigma}\left(r\right) &= \frac{2}{n} \frac{1}{(2\pi)^2} \int_0^{k_F} dk k^2 \int_{-1}^{1} d \cos \theta e^{-ikr \cos \theta} = \frac{3}{(k_F r)^3}[\sin(k_F r) - (k_F r) \cos(k_F r)] \\&= \frac{3j_1(k_F r)}{k_F r}. \end{aligned} \] Here we used again the result $k_F^3 = 3\pi^2 n$.

We can express the single-particle correlation function of the Fermi gas in terms of the spherical Bessel function $j_1(x)$: \[ G_\sigma(\vec r-\vec r') = G_\sigma(r) = 3 \frac{j_1 (k_F r)}{k_F r}. \] Note that the size of the correlation function depends only on the distance between the two points, not on their position or the orientation of the vector connecting them. This reflects the translation and rotation symmetry of the Fermi gas.

We can plot it:

$\begin{tikzpicture} \begin{axis}[ scale only axis=true, axis x line=middle, axis y line=left, xlabel = $k_F r$, ylabel = {$G_{\sigma}(r)$}, ymin=-0.2, ymax=1.1, style=thick, width=0.6\textwidth, height=0.35\textwidth, ] \addplot [color=blue] coordinates { (0,1) (0.1,0.999) (0.4,0.9841) (0.7,0.9518) (1,0.9035) (1.3,0.8409) (1.6,0.7663) (1.9,0.6826) (2.2,0.5926) (2.5,0.4995) (2.8,0.4063) (3.1,0.3161) (3.4,0.2314) (3.7,0.1545) (4,0.0871) (4.3,0.0305) (4.6,-0.0147) (4.9,-0.0484) (5.2,-0.0708) (5.5,-0.0830) (5.8,-0.0861) (6.1,-0.0817) (6.4,-0.0714) (6.7,-0.0571) (7,-0.0404) (7.3,-0.0231) (7.6,-0.0064) (7.9,0.0083) (8.2,0.0203) (8.5,0.0289) (8.8,0.0340) (9.1,0.0356) (9.4,0.0340) (9.7,0.0298) (10,0.0235) (10.3,0.0160) (10.6,0.0080) (10.9,0.0001) (11.2,-0.0069) (11.5,-0.0127) (11.8,-0.0168) (12.1,-0.0191) (12.4,-0.0195) (12.7,-0.0182) (13,-0.0155) (13.3,-0.0117) (13.6,-0.0073) (13.9,-0.0026) (14.2, 0.0020) (14.5,0.0060) (14.8,0.0092) }; \end{axis} \end{tikzpicture}$ Single-particle correlation function of the Fermi gas oscillates within a scale of the order of the Fermi wavelength $2\pi/k_F = \lambda_F$ and tends to zero for $r \gg \lambda_F$ .

Pair distribution function

We can also define a pair distribution function of a Fermi gas via \[ g_{\sigma\sigma'}(\vec{r}-\vec{r}') \equiv (\frac{2}{n})^2 \bra{g} \Psi^\dagger(\vec r,\sigma) \Psi^\dagger (\vec r',\sigma') \Psi(\vec r',\sigma') \Psi(\vec r,\sigma) \ket{g}. \] The pair distribution function describes the conditional probability amplitude of finding a particle of spin $\sigma$ at $\vec r$, if there is another particle of spin $\sigma'$ at $\vec r'$.

In the exercises, you will show that

If $\sigma \neq \sigma'$, $g_{\sigma\sigma'}(\vec{r}-\vec{r}') = 1.$
If $\sigma=\sigma'$, $g_{\sigma\sigma}(\vec{r}-\vec{r}') = 1-G_\sigma(r)^2$.

Plotting $g_{\sigma\sigma}(r)$:

$\begin{tikzpicture} \begin{axis}[ scale only axis=true, axis x line=middle, axis y line=left, xlabel = $k_F r$, ylabel = {$g_{\sigma\sigma}(r)$}, ymin=-0.2, ymax=1.1, style=thick, width=0.6\textwidth, height=0.35\textwidth, ] \addplot [color=blue] coordinates { (0,0) (0.1,1-0.999) (0.4,1-0.9841) (0.7,1-0.9518) (1,1-0.9035) (1.3,1-0.8409) (1.6,1-0.7663) (1.9,1-0.6826) (2.2,1-0.5926) (2.5,1-0.4995) (2.8,1-0.4063) (3.1,1-0.3161) (3.4,1-0.2314) (3.7,1-0.1545) (4,1-0.0871) (4.3,1-0.0305) (4.6,1+0.0147) (4.9,1+0.0484) (5.2,1+0.0708) (5.5,1+0.0830) (5.8,1+0.0861) (6.1,1+0.0817) (6.4,1+0.0714) (6.7,1+0.0571) (7,1+0.0404) (7.3,1+0.0231) (7.6,1+0.0064) (7.9,1-0.0083) (8.2,1-0.0203) (8.5,1-0.0289) (8.8,1-0.0340) (9.1,1-0.0356) (9.4,1-0.0340) (9.7,1-0.0298) (10,1-0.0235) (10.3,1-0.0160) (10.6,1-0.0080) (10.9,1-0.0001) (11.2,1+0.0069) (11.5,1+0.0127) (11.8,1+0.0168) (12.1,1+0.0191) (12.4,1+0.0195) (12.7,1+0.0182) (13,1+0.0155) (13.3,1+0.0117) (13.6,1+0.0073) (13.9,1+0.0026) (14.2,1-0.0020) (14.5,1-0.0060) (14.8,1-0.0092) }; \end{axis} \end{tikzpicture}$ shows that $g_{\sigma\sigma}(r\rightarrow 0) \rightarrow 0$ , whereas $g_{\sigma\sigma}(r \gg \lambda_F) \rightarrow 1$ . This is a signature of an effective repulsion between pairs of equal-spin fermions even in a non-interacting gas! This repulsion originates from the Pauli exclusion principle.

Excited particle states from the Fermi gas

So far we have only discussed the ground state properties of the Fermi gas. To describe dynamics, we must consider the other states of the system, the excited states.

In the Fermi gas, the simplest possible excited state with $N$ particles is \[ |{\rm e}\rangle = a_{\vec q'\sigma'}^\dagger a_{\vec q\sigma} |{\rm g}\rangle, \] with $|\vec q| \leq p_F$ and $|\vec q'| >p_F$. This state is obtained by taking a single particle ($\vec q,\sigma$) from a ground state and moving it to a state ($\vec q',\sigma'$) above the Fermi surface. We say this state has a particle excitation with ($\vec q',\sigma'$) and a hole excitation with ($\vec q,\sigma$).

The energy of such an excited state is $E=E_0 -\varepsilon_q +\varepsilon_{q'}$ (show this as an exercise!), which is always higher than the ground state energy $E_0$. We can now choose to measure the single particle energies not from the band bottom, but relative to the Fermi energy. We define a new dispersion relation \[ \xi_q = \begin{cases} \varepsilon_q-\varepsilon_F,\quad |\vec q| > p_F \qquad\text{particles}\\ \varepsilon_F-\varepsilon_q,\quad |\vec q| \leq p_F\qquad\text{holes} \end{cases} \] which describes the energy of a particle and hole excitations. Both excitations have a positive energy. In terms of $\xi_q$, the energy of the excited state relative to GS is then a sum of positive quantities, \[ E-E_0= (\varepsilon_F-\varepsilon_q) + (\varepsilon_q'-\varepsilon_F) = \xi_q +\xi_{q'}. \] i.e. the excited state energy is a sum of particle and hole excitation energies. This also applies to more complicated excited states of the Fermi gas.

$\begin{tikzpicture} \begin{axis}[ axis x line=middle, axis y line=middle, xlabel = $p_x/p_F$, ylabel = {$\xi_p/\varepsilon_F$}, ymin=-0.5, ymax=2.5 ] %Here the blue parabloa is defined \addplot [ domain=-1.8:1.8, samples=100, color=blue, style=thick ] {(x^2-1)*sign(x^2-1)}; \end{axis} \end{tikzpicture}$

Fig:Particle-hole excitation spectrum in the p_x direction.

We now define new creation and annihilation operators to match the dispersion $\xi_q$: \[ \begin{aligned} c_{\vec q,\sigma} &\equiv a_{-\vec q,-\sigma}^\dagger,\quad |\vec q|\leq p_F, \qquad\text{hole}\\ c_{\vec q,\sigma}^\dagger &\equiv a_{-\vec q,-\sigma},\quad |\vec q|\leq p_F, \qquad\text{hole}\\ c_{\vec q,\sigma} &\equiv a_{\vec q,\sigma},\quad |\vec q|>p_F, \qquad\text{particle}\\ c_{\vec q,\sigma}^\dagger &\equiv a_{\vec q,\sigma}^\dagger,\quad |\vec q|>p_F. \qquad\text{particle} \end{aligned} \]

The creation and annihilation operators for particles are the same as the original operators. For holes, the annihilation and creation operators switch their roles and we also invert the momentum and spin, so that $b_{\vec q,\sigma}^\dagger$ always increases the total momentum and spin of the system by $\vec p$ and $\sigma/2$, respectively.

The new operators have the same anticommutation relations as the original operators: \[ \{ c_{\vec q,\sigma}, c_{\vec q',\sigma'}^\dagger\}=\delta_{\vec q,\vec q'} \delta_{\sigma\sigma'},\\ \{ c_{\vec q,\sigma}^\dagger, c_{\vec q',\sigma'}^\dagger\} = \{ c_{\vec q,\sigma}, c_{\vec q',\sigma'}\} = 0. \]

The vacuum has the property $a_{\vec q,\sigma} |0\rangle=0$, i.e. it is quenched by all original annihilation operators. Similarly $c_{\vec q,\sigma}|{\rm g}\rangle = 0$ for all particle-hole annihilation operators. Below the Fermi surface this follows from the Pauli exclusion principle ($(a^\dagger_{\vec q,\sigma})^2=0$). The ground state of the Fermi gas is an effective vacuum for the particle-hole operators. This has practical consequences, as it makes the calculation of expectation values of particle-hole operator products easy. For example, instead of an expectation value $(2N+4)$-operator product, we might only need to evaluate an expectation value of $4$-operator product.

The Hamiltonian of the Fermi gas also needs to be transformed to the particle-hole basis. To do this, we first divide it into two parts, of which one will describe particles and the other holes: \[ H= \sum_{|\vec p|>p_F,\sigma} \varepsilon_{p} a^\dagger_{\vec p,\sigma} a_{\vec p,\sigma} + \sum_{|\vec p|\leq p_F,\sigma} \varepsilon_{p} a^\dagger_{\vec p,\sigma} a_{\vec p,\sigma}. \] The two parts can then be separately expressed in terms of the new dispersion $\xi_p$ and particle-hole operators. We get \[ \begin{aligned} H &= \sum_{|\vec p|>p_F,\sigma} (\varepsilon_{p}-\varepsilon_F) c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma} + \varepsilon_F \sum_{|\vec p|>p_F,\sigma} c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma}\\ &+\sum_{|\vec p|\leq p_F,\sigma} (\varepsilon_F-\varepsilon_{p}) c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma} - \varepsilon_F \sum_{|\vec p|\leq p_F,\sigma} c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma} + \sum_{|\vec p|\leq p_F,\sigma} \varepsilon_p\\ &= \sum_{\vec p,\sigma} \xi_p c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma}\;+\; \varepsilon_F(\hat N_p- \hat N_h)\; +\; E_0, \end{aligned} \] where $\hat N_p = \sum_{|p|>p_F,\sigma} c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma}$ and $\hat N_h = \sum_{|p|\leq p_F,\sigma} c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma}$ are the particle and hole number operators. The particle number operator of the original particles is \[ \hat N = \sum_{\vec p,\sigma} a^\dagger_{\vec p\sigma} a_{\vec p\sigma} = N + \hat N_p - \hat N_h \] Thus, the Hamiltonian is \[ H=\sum_{\vec p,\sigma} \xi_p c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma}\;+\; \varepsilon_F(\hat N-N)\; +\; E_0. \] For a system with constant number of particles, the second term is 0, and the last term is a constant energy shift which can be scaled away. We end up with a simple Hamiltonian \[ H=\sum_{\vec p,\sigma} \xi_p c^\dagger_{\vec p,\sigma} c_{\vec p,\sigma}. \] After dropping the $E_0$ term, the expectation value of the Hamiltonian corresponds to the energy of the system relative to the effective vacuum $|{\rm g}\rangle$, and not relative to the original vacuum $|0\rangle$.

Above, we calculated the energy of the simple excited state $|{\rm e}\rangle \equiv a^\dagger_{\vec q',\sigma'} a_{\vec q,\sigma}|{\rm g}\rangle = c^\dagger_{\vec q',\sigma'}c_{-\vec q,-\sigma}^\dagger|{\rm g}\rangle,$ where $|\vec q'|>p_F$ and $|\vec q|\leq p_F$ . By using the transformed Hamiltonian and Wick's theorem, we find $\begin{align*} E-E_0 &= \langle {\rm e} | H |{\rm e}\rangle\\ &= \sum_{\vec p,\sigma''} \xi_p \bigg(\langle {\rm g}| \contraction[1.05ex]{}{c}{\vphantom{c}_{-\vec q,-\sigma} c_{\vec q',\sigma'}}{c} \contraction[1.05ex]{c_{-\vec q,-\sigma} c_{\vec q',\sigma'} c^\dagger_{\vec p\sigma''}}{c}{\vphantom{c}_{\vec p\sigma''} c^\dagger_{\vec q',\sigma'}}{c} \contraction[2.05ex]{c_{-\vec q,-\sigma}}{c}{\vphantom{c}_{\vec q',\sigma'} c^\dagger_{\vec p\sigma''}c_{\vec p\sigma''}}{c} c_{-\vec q,-\sigma} c_{\vec q',\sigma'} c^\dagger_{\vec p\sigma''} c_{\vec p\sigma''} c^\dagger_{\vec q',\sigma'}c_{-\vec q,-\sigma}^\dagger |{\rm g}\rangle\\ &\mkern63mu+ \langle {\rm g}| \contraction[2.05ex]{}{c}{\vphantom{c}_{-\vec q,-\sigma} c_{\vec q',\sigma'} c^\dagger_{\vec p\sigma''}c_{\vec p\sigma''} c^\dagger_{\vec q',\sigma'}}{c} \contraction[1.05ex]{c_{-\vec q,-\sigma}}{c}{\vphantom{c}_{\vec q',\sigma'}}{c} \contraction[1.05ex]{c_{-\vec q,-\sigma} c_{\vec q',\sigma'} c^\dagger_{\vec p\sigma''}}{c}{\vphantom{c}_{\vec p\sigma''}}{c} c_{-\vec q,-\sigma} c_{\vec q',\sigma'} c^\dagger_{\vec p\sigma''} c_{\vec p\sigma''} c^\dagger_{\vec q',\sigma'}c_{-\vec q,-\sigma}^\dagger |{\rm g}\rangle \bigg)\\ &= \xi_q + \xi_{q'}, \end{align*}$ which agrees with the earlier calculation.

Interacting Fermi system

Let us now add the two-body interactions into the Hamiltonian so that the total Hamiltonian is of the form \[ H = \sum_{k,l} \bra{k} \hat H_0 \ket{l} a_k^\dagger a_l + \frac{1}{2} \sum_{ijkl} \langle ij |\hat V_2 |kl\rangle a_i^\dagger a_j^\dagger a_l a_k. \]

In the case where the single-particle eigenstates are spin degenerate and described by some momentum $\vec k$, and when the two-body interaction is independent of translations, we may simplify this to the form

\[ H=\sum_{\vec{k},\sigma} \epsilon_{\vec{k}} a_{\vec{k} \sigma}^\dagger a_{\vec{k}\sigma} + \frac{1}{2V} \sum_{\sigma \sigma'} \sum_{{\vec{k}},{\vec{k}}',\vec q} V_2(\vec q) a_{{\vec{k}}+\vec{q},\sigma}^\dagger a_{{\vec{k}}'-\vec q,\sigma'}^\dagger a_{{\vec{k}}' \sigma'} a_{{\vec{k}}\sigma}. \]

Finding the ground state of this Hamiltonian for a generic $N$-particle system (with a large $N$) turns out to be very hard if not impossible. This generic problem is discussed for example in Ch. 1 of this book by Xiao-Gang Wen. When there are only a few particles, such exact ground states may be found computationally. However, in the case of a large number of particles majority of approaches at least start from an approximative scheme introduced below: that of a Hartree-Fock or mean field method, whose idea is to try to pose the interacting Hamiltonian in terms of the closest analogue non-interacting Hamiltonian, i.e., in terms of an effective Fermi gas, where the "other" particles act as an effective potential term in the single-particle problem. In some cases this emerging effective potential spontaneously breaks some symmetry of the system. Such symmetry breaking can be associated with phase transitions. Below we describe two such models, related with spontaneous symmetry breaking in the spin sector (and hence describing magnetism) or one where we have to redefine the fermionic states as in the BCS theory of superconductivity.

Hartree-Fock (mean field) method

The Hartree-Fock method is perhaps the most often used approach to finding an approximate ground state for a system of a large number $N$ of pair-wise interacting fermions. The main idea is to assume that the ground state can still be written in terms of some Slater determinant, \[ |\Phi\rangle = \prod_{k=1}^N a_k^\dagger |0\rangle = |1_1 1_2 \dots 1_N 0_{N+1} 0_{N+2} \dots\rangle, \] where the quantum numbers $k$ that are included in this product are chosen such that the energy \[ E_{\Phi}=\bra{\Phi} \hat H \ket{\Phi} \] is minimized. The approach is thus similar to the one in the variational principle discussed in earlier quantum mechanics courses.

In the Ansatz, we assume that the quantum numbers are organized such that the states $k \in \{1,\dots,N\}$ correspond to occupied states, and the others ($k>N$) to unoccupied states.

Now suppose our Hamiltonian is \[ H = \sum_{k,l} \bra{k} \hat H_0 \ket{l} a_k^\dagger a_l + \frac{1}{2} \sum_{ijkl} \langle ij |\hat V_2 |kl\rangle a_i^\dagger a_j^\dagger a_l a_k. \]

Let us compute the expectation value of the energy $E_\Phi$: \[ \bra{\Phi} \hat H \ket{\Phi} = \sum_{kl} \langle k | \hat H_0 | l \rangle \bra{\Phi} a_k^\dagger a_l \ket{\Phi} + \frac{1}{2}\sum_{ijkl}\langle ij | \hat V_2 |kl \rangle \bra{\Phi} a_i^\dagger a_j^\dagger a_l a_k \ket{\Phi}. \] These can be simplified under the assumption that $\ket{\Phi}$ is a Slater determinant. \[ \bra{\Phi} a_k^\dagger a_l \ket{\Phi} =\begin{cases} 0, \text{ if } k \text{ or } l > N\\\delta_{kl}, \text{ for } k,l \in \{1,\dots,N\} \end{cases} \] and \[ \bra{\Phi} a_i^\dagger a_j^\dagger a_l a_k \ket{\Phi} = 0 \text{ if one of the indices } > N. \] We also notice that the latter term is non-zero only if (i) $i=l$ and $j=k$ or (ii) $i=k$ and $j=l$ (and all are $\in \{1,\dots,N\}$). In other words \[ \begin{aligned} \bra{\Phi} a_i^\dagger a_j^\dagger a_l a_k \ket{\Phi} &=\delta_{ik} \delta_{jl} \bra{\Phi} a_i^\dagger a_j^\dagger a_j a_k \ket{\Phi} + \delta_{il} \delta_{jk}\bra{\Phi} a_i^\dagger a_j^\dagger a_k a_l \ket{\Phi}\\&= \delta_{ik} \delta_{jl} - \delta_{il} \delta_{jk}. \end{aligned} \] When showing this, it is enough to assume $i\neq j$ as for $i=j$ we would be removing two fermions from the same state. Then we simply need to shift the annihilation operator next to the corresponding creation operator and notice that we are calculating the expectation value of $\hat n_i \hat n_j$. At this point you can probably infer the reason for the sign change?

Plugging this result to the expression for the energy and using the Kronecker $\delta$-functions to get rid of some of the sums, the result is \[ E_\Phi = \bra{\Phi} H \ket{\Phi} = \sum_{k=1}^N \langle k | \hat H_0 |k\rangle + \frac{1}{2} \sum_{kl} (\langle kl | \hat V_2 |kl\rangle - \langle kl |\hat V_2 | lk\rangle). \] The two pairwise interaction terms are frequently called the "direct" and "exchange" terms or the "Hartree" and the "Fock" terms.

Minimizing $E_\Phi$ in general leads to a set of non-linear equations for the Ansatz states. We do not write them down here but instead concentrate below on a couple of example cases. However, the generic Hartree-Fock method is widely used for example for the

Atomic shell model
To describe an interacting electron gas (or an "electron liquid")
Nuclear shell model
and many other systems

Self-consistent mean field method

An alternative formulation of the Hartree-Fock theory can be made so that we replace the two-body interaction term by an effective one-body potential that depends on the state of the system. In other words, we approximate in the two-body interaction \[ \begin{aligned} a_i^\dagger a_j^\dagger a_l a_k \approx & \langle \Phi| a_i^\dagger a_k | \Phi \rangle a_j^\dagger a_l + \langle \Phi| a_j^\dagger a_l |\Phi \rangle a_i^\dagger a_k \\&- \langle \Phi| a_i^\dagger a_l |\Phi\rangle a_j^\dagger a_k-\langle \Phi| a_j^\dagger a_k |\Phi \rangle a_i^\dagger a_l. \end{aligned} \] The two-body interaction term is hence of the form \[ \frac{1}{2}\sum_{ijkl} \langle ij | \hat V_2 |kl\rangle \langle \Phi| a_i^\dagger a_k | \Phi \rangle a_j^\dagger a_l + \langle \Phi| a_j^\dagger a_l |\Phi \rangle a_i^\dagger a_k - \langle \Phi| a_i^\dagger a_l |\Phi\rangle a_j^\dagger a_k-\langle \Phi| a_j^\dagger a_k |\Phi \rangle a_i^\dagger a_l. \] Since for indistinguishable particles the two-body interaction satisfies \[ \langle ij |\hat V_2 |kl\rangle = \langle ji |\hat V_2 |lk\rangle, \] we can interchange the dummy indices $i \leftrightarrow j$ and $k \leftrightarrow l$ in the second and the fourth term. As a result, we get an effective single-particle Hamiltonian \[ H = \sum_{jl} \left[\bra{j} \hat H_0 \ket{l} + \underbrace{\sum_{ik} \bra{ij} \hat V_2 \ket{kl} \bra{\Phi} a_i^\dagger a_k \ket{\Phi}}_{\langle j | \hat V_H | l \rangle}- \underbrace{\sum_{ik} \bra{ij} \hat V_2 \ket{lk} \bra{\Phi} a_i^\dagger a_k \ket{\Phi}}_{\langle j|\hat V_F |k\rangle}\right] a_j^\dagger a_l. \] This procedure defines two "potentials" denoted by operators $\hat V_H$ and $\hat V_F$ corresponding to the Hartree and Fock interaction terms. In principle the problem thus reduces to finding the single-particle spectrum under such potentials. The ground state $|\Phi\rangle$ then corresponds to filling $N$ lowest-energy eigenstates of this new problem. However, now the potentials $\hat V_H$ and $V_F$ themselves depend on $|\Phi\rangle$. The usual procedure for solving this problem is that of a self-consistency iteration:

Start from an Ansatz $|\Phi\rangle$
Compute $\hat V_H$ and $\hat V_F$ (or their representation in a chosen basis)
Find the spectrum of the resulting single-particle Hamiltonian $\hat H_0 + \hat V_H + \hat V_F$
Choose the $N$ lowest-energy states as a new Ansatz $|\Phi\rangle$
Iterate until convergence

Examples of such self-consistent field method are discussed below.

# Selfconsmeanfield

System with translational invariance

Let us consider the case where the system contains translational invariance and the two-body interaction is spin-conserving, as discussed above. In that case the effective one-body Hamiltonian is (note that strictly speaking the following assumes that spin $\sigma$ is a good quantum number - this assumption can also be lifted)

\[ \begin{aligned} H&=\sum_{\vec k \sigma} \epsilon_k a_{\vec k \sigma}^\dagger a_{\vec k \sigma} + \frac{1}{V}\sum_{\sigma \sigma'} \sum_{\vec k\vec k'\vec q} V_2(\vec q) (\underbrace{\langle \Phi| a_{\vec k+\vec q \sigma}^\dagger a_{\vec k \sigma} |\Phi\rangle}_{\delta_{\vec q \vec 0} n_{\vec k \sigma}} a_{\vec k'-\vec q \sigma'}^\dagger a_{\vec k' \sigma'} - \underbrace{\bra{\Phi} a_{\vec k+\vec q \sigma}^\dagger a_{\vec k'\sigma'} \ket{\Phi}}_{\delta_{\vec k+\vec q,\vec k'}\delta_{\sigma \sigma'} n_{\vec k' \sigma}} a_{\vec k'-\vec q \sigma'}^\dagger a_{\vec k \sigma})\\ &=\sum_{\vec k \sigma} [\epsilon_k+V_H-V_F(\vec k,\sigma)] a_{\vec k \sigma}^\dagger a_{\vec k \sigma}, \end{aligned} \] where we defined two self-consistent fields, containing the constant Hartree term \[ V_H=\frac{1}{V} \sum_{\vec k' \sigma} V_2(0) n_{\vec k' \sigma} \] and the momentum and spin dependent Fock term \[ V_F(\vec k,\sigma) = \frac{1}{V}\sum_{\vec q} V_2(\vec q) n_{\vec k+\vec q \sigma}. \] In what follows, we show that the spin dependence in the Fock term may lead to a spontaneous symmetry breaking to a state with non-vanishing spin polarization.