9.1: Interacting Electrons in Atomic Shells

Previously, you have encountered the Schrödinger equation for a particle in a central potential, which can be interpreted as an electron bound to a proton in a hydrogen atom. We found that the energy levels are quantised with quantum number \(n\). In addition, the spin and orbital angular momentum is quantised with quantum numbers \(m_{s}\), \(l\), and \(m\), respectively (we assume that \(s=\frac{1}{2}\) since we are considering electrons). We can denote the set of quantum numbers \(n\), \(m_{s}\), \(l\), and \(m\) by greek indices \(\alpha\), \(\beta\), . . . , and the wave-functions \(u_{\alpha}(\mathbf{r})\) then form a complete orthonormal basis for the bound electron.

It is tempting to keep this complete orthonormal basis for other atoms as well, and assume that the ground state of an \(N\)-electron atom is the tensor product of the \(N\) lowest energy eigenstates, appropriately anti-symmetrized via the Slater determinant. Indeed, the periodic table is based on this assumption. However, this ignores the fact that the electrons interact with each other, and the ground state of a many-electron atom is different. The Hartree-Fock method is designed to take this into account. It is a constrained variational approach, in which the trial state that is optimised over is forced to be a Slater determinant in order to keep the correct particle statistics. In this section we present the Hartree-Fock method, and arrive at the Hartree-Fock equations, which can be solved iteratively. We follow the derivation given by Bransden and Joachain (Physics of Atoms and Molecules, 1983 pp. 320-339).

First, we specify the Hamiltonian. Using the notation \(r_{i}=\left|\mathbf{r}_{i}\right|\) for the distance of the \(i^{\text {th }}\) electron from the nucleus, and \(r_{i j}=\left|\mathbf{r}_{i}-\mathbf{r}_{j}\right|\) for the distance between electrons \(i\) and \(j\), we find that

\[H=H_{1}+H_{2}=-\sum_{i=1}^{N}\left(\frac{\hbar^{2}}{2 m} \nabla_{i}^{2}+\frac{Z}{r_{i}}\right)+\sum_{i<j=1}^{N} \frac{1}{r_{i j}},\tag{9.1}\]

where we used units in which \(e / 4 \pi \epsilon_{0}=1\) and the Hamiltonian is divided into the single electron Hamiltonian \(\left(H_{1}\right)\) and the inter-electron Hamiltonian \(\left(H_{2}\right)\). We choose as a normalised trial quantum state \(\Psi\left(\mathbf{r}_{1}, \ldots, \mathbf{r}_{N}\right)\) the Slater determinant

\[\Psi\left(\mathbf{r}_{1}, \ldots, \mathbf{r}_{N}\right)=\frac{1}{\sqrt{N !}}\left|\begin{array}{cccc}
u_{\alpha}\left(\mathbf{r}_{1}\right) & u_{\beta}\left(\mathbf{r}_{1}\right) & \ldots & u_{v}\left(\mathbf{r}_{1}\right) \\
u_{\alpha}\left(\mathbf{r}_{2}\right) & u_{\beta}\left(\mathbf{r}_{2}\right) & \ldots & u_{v}\left(\mathbf{r}_{2}\right) \\
\vdots & \vdots & \ddots & \vdots \\
u_{\alpha}\left(\mathbf{r}_{N}\right) & u_{\beta}\left(\mathbf{r}_{N}\right) & \ldots & u_{v}\left(\mathbf{r}_{N}\right)
\end{array}\right|,\tag{9.2}\]

and calculate the expectation value of the Hamiltonian \(H\). This must be larger or equal to the ground state \(\boldsymbol{E}_{0}\)

\[\langle\Psi|H| \Psi\rangle \geq E_{0},\tag{9.3}\]

and varying the trial state then allows us to minimise the expectation value. This will get us close to the ground state energy.

Since the Slater determinant is a rather large expression, it will save us quite a bit of writing if we introduce the anti-symmetrisation operator \(\mathscr{A}\), such that

\[\Psi\left(\mathbf{r}_{1}, \ldots, \mathbf{r}_{N}\right)=\sqrt{N !} \mathscr{A} u_{\alpha}\left(\mathbf{r}_{1}\right) u_{\beta}\left(\mathbf{r}_{2}\right) \ldots u_{v}\left(\mathbf{r}_{N}\right) \equiv \sqrt{N !} \mathscr{A} \Phi_{\mathrm{H}},\tag{9.4}\]

where \(\Phi_{\mathrm{H}}\left(\mathbf{r}_{1}, \ldots, \mathbf{r}_{N}\right)\) is the Hartree wave function. The operator \(\mathscr{A}\) can then be written as a sum over permutations \(P\) of the labels \(\alpha, \beta, \ldots v\):

\[\mathscr{A}=\frac{1}{N !} \sum_{P}(-1)^{P} P,\tag{9.5}\]

and \(\mathscr{A}\) is a projection operator: \(\mathscr{A}^{2}=\mathscr{A}=\mathscr{A}^{\dagger}\). Both \(H_{1}\) and \(H_{2}\) commute with \(\mathscr{A}\).

Next, we calculate the expectation values \(\left\langle\Psi\left|H_{1}\right| \Psi\right\rangle\) and \(\left\langle\Psi\left|H_{2}\right| \Psi\right\rangle\). Since \(\left[H_{1}, \mathscr{A}\right]=0\), we can write

\[\left\langle\Psi\left|H_{1}\right| \Psi\right\rangle=N !\left\langle\Phi_{\mathrm{H}}\left|\mathscr{A} H_{1} \mathscr{A}\right| \Phi_{\mathrm{H}}\right\rangle=N !\left\langle\Phi_{\mathrm{H}}\left|H_{1} \mathscr{A}^{2}\right| \Phi_{\mathrm{H}}\right\rangle=N !\left\langle\Phi_{\mathrm{H}}\left|H_{1} \mathscr{A}\right| \Phi_{\mathrm{H}}\right\rangle.\tag{9.6}\]

A permutation of the labels \(\alpha, \beta, \ldots v\) leads to an orthonormal state and \(H_{1}\) is a sum over oneelectron Hamiltonians. We can therefore write this as

\[\begin{aligned}
\left\langle\Psi\left|H_{1}\right| \Psi\right\rangle &=N !\left\langle\Phi_{\mathrm{H}}\left|H_{1} \mathscr{A}\right| \Phi_{\mathrm{H}}\right\rangle=\sum_{P}(-1)^{P}\left\langle\Phi_{\mathrm{H}}\left|H_{1} P\right| \Phi_{\mathrm{H}}\right\rangle=\sum_{i=1}^{N}\left\langle\Phi_{\mathrm{H}}\left|\hat{h}_{i}\right| \Phi_{\mathrm{H}}\right\rangle \\
&=\sum_{\alpha}\left\langle u_{\alpha}\left(\mathbf{r}_{i}\right)\left|\hat{h}_{i}\right| u_{\alpha}\left(\mathbf{r}_{i}\right)\right\rangle \equiv \sum_{\alpha} I_{\alpha}
\end{aligned}\tag{9.7}\]

where \(\hat{h}_{i}\) is the single electron Hamiltonian

\[\hat{h}_{i}=-\frac{\hbar^{2}}{2 m} \nabla_{i}^{2}-\frac{Z}{r_{i}}.\tag{9.8}\]

Next, we calculate the expectation value \(\left\langle\Psi\left|H_{2}\right| \Psi\right\rangle\) of the two-electron interaction Hamiltonians. Using the same reasoning as in Eq. (9.6), we find that

\[\left\langle\Psi\left|H_{2}\right| \Psi\right\rangle=N !\left\langle\Phi_{\mathrm{H}}\left|H_{2} \mathscr{A}\right| \Phi_{\mathrm{H}}\right\rangle.\tag{9.9}\]

Substituting the explicit form of \(\mathscr{A}\), we find

\[\left\langle\Psi\left|H_{2}\right| \Psi\right\rangle=\sum_{i<j} \sum_{P}(-1)^{P}\left\langle\Phi_{\mathrm{H}}\left|\frac{P}{r_{i j}}\right| \Phi_{\mathrm{H}}\right\rangle=\sum_{i<j}\left\langle\Phi_{\mathrm{H}}\left|\frac{1-P_{i j}}{r_{i j}}\right| \Phi_{\mathrm{H}}\right\rangle,\tag{9.10}\]

where \(P_{i j}\) is the exchange operator of electrons \(i\) and \(j\). This expression allows us to write

\[\begin{array}{r}
\left\langle\Psi\left|H_{2}\right| \Psi\right\rangle=\frac{1}{2} \sum_{\alpha, \beta}\left[\left\langle u_{\alpha}\left(\mathbf{r}_{i}\right) u_{\beta}\left(\mathbf{r}_{j}\right)\left|\frac{1}{r_{i j}}\right| u_{\alpha}\left(\mathbf{r}_{i}\right) u_{\beta}\left(\mathbf{r}_{j}\right)\right\rangle\right. \\
\left.-\left\langle u_{\alpha}\left(\mathbf{r}_{i}\right) u_{\beta}\left(\mathbf{r}_{j}\right)\left|\frac{1}{r_{i j}}\right| u_{\beta}\left(\mathbf{r}_{i}\right) u_{\alpha}\left(\mathbf{r}_{j}\right)\right\rangle\right]
\end{array}.\tag{9.11}\]

Note the swap of \(\alpha\) and \(\beta\) in the last ket. This expectation value consists of two terms, namely the direct term

\[J_{\alpha \beta} \equiv\left\langle u_{\alpha}\left(\mathbf{r}_{i}\right) u_{\beta}\left(\mathbf{r}_{j}\right)\left|\frac{1}{r_{i j}}\right| u_{\alpha}\left(\mathbf{r}_{i}\right) u_{\beta}\left(\mathbf{r}_{j}\right)\right\rangle\tag{9.12}\]

and the exchange term

\[K_{\alpha \beta} \equiv\left\langle u_{\alpha}\left(\mathbf{r}_{i}\right) u_{\beta}\left(\mathbf{r}_{j}\right)\left|\frac{1}{r_{i j}}\right| u_{\beta}\left(\mathbf{r}_{i}\right) u_{\alpha}\left(\mathbf{r}_{j}\right)\right\rangle\tag{9.13}\]

The total expectation value therefore becomes

\[\langle\Psi|H| \Psi\rangle=\sum_{\alpha} I_{\alpha}+\frac{1}{2} \sum_{\alpha, \beta}\left(J_{\alpha \beta}-K_{\alpha \beta}\right).\tag{9.14}\]

The matrix elements \(J_{\alpha \beta}\) and \(K_{\alpha \beta}\) are real and symmetric in \(\alpha\) and \(\beta\).

The second step towards the Hartree-Fock method is to find the minimum of \(E \equiv\langle\Psi|H| \Psi\rangle\) by varying the \(u_{\alpha}\left(\mathbf{r}_{i}\right)\). This means finding \(\delta E=0\). However, we must keep the functions \(u_{\alpha}\left(\mathbf{r}_{i}\right)\) orthonormal to each other, and this imposes \(N^{2}\) constraints. We can incorporate these constraints in the variational procedure by including Lagrange multipliers \(\epsilon_{\alpha \beta}\), and the variational equation becomes

\[\delta E-\sum_{\alpha, \beta} \epsilon_{\alpha \beta} \delta\left\langle u_{\alpha}(\mathbf{r}) \mid u_{\beta}(\mathbf{r})\right\rangle=0.\tag{9.15}\]

There is no explicit reference to electron positions \(\mathbf{r}_{i}\) in \(\left\langle u_{\alpha}(\mathbf{r}) \mid u_{\beta}(\mathbf{r})\right\rangle\) since we are only interested in its orthonormality properties. The Lagrange multipliers \(\epsilon_{\alpha \beta}\) form the elements of a Hermitian matrix.

The variational approach ultimately leads to a set of \(N\) coupled equations:

\[\begin{gathered}
E_{\alpha} u_{\alpha}\left(\mathbf{r}_{i}\right)=\hat{h}_{i} u_{\alpha}\left(\mathbf{r}_{i}\right)+\sum_{\beta}\left\langle u_{\beta}\left(\mathbf{r}_{j}\right)\left|r_{i j}^{-1}\right| u_{\beta}\left(\mathbf{r}_{j}\right)\right\rangle u_{\alpha}\left(\mathbf{r}_{i}\right) \\
-\sum_{\beta}\left\langle u_{\beta}\left(\mathbf{r}_{j}\right)\left|r_{i j}^{-1}\right| u_{\alpha}\left(\mathbf{r}_{j}\right)\right\rangle u_{\beta}\left(\mathbf{r}_{i}\right),
\end{gathered}\tag{9.16}\]

which are known as the Hartree-Fock equations. These can be solved by iteration up to any desired precision.