# 8.3: Hartree-Fock Equations are Solved by the Self-Consistent Field Method

The Hartree method is used to determined the wavefunction and the energy of a quantum multi-electron system in a stationary state. The Hartree method often assumes that the exact, N-body wave function of the system can be approximated by a product of single-electron wavefucntions. By invoking the variational method, one can derive a set of N-coupled equations for the N spin orbitals. A solution of these equations yields the Hartree wavefunction and energy of the system. It is one step better than the "Ignorance is Bliss approach", but still far from a good approximation.

### The Unsolvable System

The generic multielectron atom including terms for the additional electrons with a general charge \(Z\); e.g.

\[V_{nuclear-electron}(r_1) = -\dfrac {Z}{\left\vert\mathbf{r} - \mathbf{R}\right\vert} \label{8.3.1}\]

in atomic units with \(\left\vert\mathbf{r} - \mathbf{R}\right\vert\) is the distance between the electron and the nucleus, The Hamiltonian must also have terms for electron-electron repulsion (also in atomic units)

\[V_{electron-electron}(r_{12}) = \dfrac {1}{\left\vert \mathbf{r} - \mathbf{r}^{\prime} \right\vert} \label{8.3.2}\]

with \(\left\vert \mathbf{r} - \mathbf{r}^{\prime} \right\vert\) is the distance between electron 1 and electron 2. So the proper mulit-electron Hamiltonian can be constructed

\[\hat {H} (r_1, r_2, ... r_n) = -\dfrac {\hbar ^2}{2m_e} \sum _i \nabla ^2_i + \sum _i V_{nuclear-electron} (r_i) + \sum _{i \ne j} V_{electron-electron} (r_{ij}) \label{8.3.3}\]

Given what we have learned from the previous quantum mechanical systems we’ve studied, we predict that exact solutions to the multi-electron Schrödinger equation would consist of a family of multi-electron wavefunctions, each with an associated energy eigenvalue. These wavefunctions and energies would describe the ground and excited states of the multi-electron atom, just as the hydrogen wavefunctions and their associated energies describe the ground and excited states of the hydrogen atom. We would predict quantum numbers to be involved, as well.

The fact that electrons interact through their electron-electron repulsion (final term in Equation \(\ref{8.3.3}\)) means that an exact wavefunction for a multi-electron system would be a single function that depends **simultaneously **upon the coordinates of all the electrons; i.e., a multi-electron wavefunction:

\[|\Psi (r_1, r_2, \cdots r_i) \rangle \label{8.3.4}\]

Unfortunately, the electron-electron repulsion terms make it impossible to find an exact solution to the Schrödinger equation for many-electron atoms.

### The Hartree Approximation

The method for finding best possible one-electron wavefunctions that was published by Hartree in 1948 and improved two years later by Fock. For the Schrödinger equation to be solvable, the variables must **be separable**. The variables are the coordinates of the electrons. To separate the variables in a way that retains information about electron-electron interactions, the electron-electron term (Equation \(\ref{8.3.1}\)) must be approximated so it depends only on the coordinates of one electron. Such an approximate Hamiltonian can account for the interaction of the electrons in **an average way**. The exact one-electron eigenfunctions of this approximate Hamiltonian then can be found by solving the Schrödinger equation. These functions are the best possible **one-electron functions**.

The Hartree approximation starts by invoking an initial *ansatz* that the multi-electron wavefucntion in Equation \(\ref{8.3.4}\) can be expanded as a product of single-electron wavefunctions

\[\Psi(\mathbf{r}_1,\mathbf{r}_2, \ldots, \mathbf{r}_N) \approx \psi_{1}(\mathbf{r}_1)\psi_{2}(\mathbf{r}_2) \ldots \psi_{N}(\mathbf{r}_N) \label{2.3}\]

from which it follows that the electrons are *independent*, and interact only via the mean-field Coulomb potential. This yields one-electron Schrödinger equations of the form

\[-\dfrac{\hbar^{2}}{2m} \nabla^{2}\psi_{i}(\mathbf{r}) + V(\mathbf{r})\psi_{i}(\mathbf{r}) = \epsilon_{i}\psi_{i}(\mathbf{r}) \label{2.4}\]

or

\[H_e(r) \psi_{i}(\mathbf{r}) = \epsilon_{i}\psi_{i}(\mathbf{r}) \label{2.4B}\]

where \(V(r)\) is the potential in which the electron moves; this includes both the nuclear-electron interaction

\[V_{nucleus}(\mathbf{r}) = -Ze^{2}\sum_{R} \dfrac{1}{\left\vert\mathbf{r} - \mathbf{R}\right\vert} \label{2.5}\]

and the mean field arising from the \(N-1\) other electrons. We smear the other electrons out into a smooth negative charge density \(\rho(\mathbf{r}')\) leading to a potential of the form

\[V_{electron}(\mathbf{r}) = -e\int d\mathbf{r}^{\prime} \rho(\mathbf{r}^{\prime}) \dfrac{1}{\left\vert \mathbf{r} - \mathbf{r}^{\prime} \right\vert} \label{2.6}\]

where

\[\rho(\mathbf{r}) = \sum_{i}^{\text{occupied}}\vert\psi(\mathbf{r})\vert^{2}.\]

The sum over runs over all occupied states.

The wavefunctions that from this approach with the Hamiltonian \(H_e(r)\) involve possess three kinds of energies:

#### Kinetic Energy

The Kinetic energy of the electron has an average value is computed by taking the expectation value of the kinetic energy operator \[-\dfrac{\hbar^2}{2m} \nabla^2\] with respect to any particular solution \(\phi_j(r)\) to the Schrödinger equation: \[KE = \langle\phi_j| -\dfrac{\hbar^2}{2m} \nabla^2 |\phi_j\rangle \]

#### Nuclear-Electron Coulombic Attraction Energy

Coulombic attraction energy with the nucleus of charge \(Z\): \[\langle\phi_j| \dfrac{-Z e^2}{\left\vert\mathbf{r} - \mathbf{R}\right\vert} |\phi_j\rangle\]

#### Electron-Electron Coulombic Repulsion Energy

Coulomb repulsion energies with all of the \(N-1\) other electrons, which are assumed to occupy other atomic orbitals denoted \(\phi_K\), with this energy computed as

\[\sum_{j\neq k} \langle\phi_j(r) \phi_k(r’) |\dfrac{e^2}{|r-r’|} | \phi_j(r) \phi_k(r’)\rangle.\label{8.3.8}\]

The Dirac notation \(\langle\phi_j(r) \phi_k(r’) |\dfrac{e^2}{|r-r’|} | \phi_j(r) \phi_k(r’)\rangle\) is used to represent the two-electron (six-dimensional) Coulomb integral

\[J_{j,k} = \int |\phi_j(r)|^2 |\phi_k(r’)|^2 \dfrac{e^2}{r-r'} dr dr’ \label{8.3.9}\]

that describes the Coulomb repulsion between the charge density \(|\phi_j(r)|^2\) for the electron in \(\phi_j\) and the charge density \(|\phi_k(r’)|^2\) for the electron in \(\phi_k\). Of course, the sum over \(k\) must be limited to exclude \(k=j\) to avoid counting a “self-interaction” of the electron in orbital \(\phi_j\) with itself.

The total energy \(\epsilon_j\) of the orbital \(\phi_j\), is the sum of the above three contributions:

\[\epsilon_J = \langle\phi_j| - \dfrac{\hbar^2}{2m} \nabla^2 |\phi_j\rangle + \langle\phi_j| \dfrac{-Z e^2}{\left\vert\mathbf{r} - \mathbf{R}\right\vert} |\phi_j\rangle + \sum_{j\neq k} \langle\phi_j(r) \phi_k(r’) |\dfrac{e^2}{|r-r’|} | \phi_j(r) \phi_k(r’)\rangle.\label{8.3.10}\]

This treatment of the electrons and their orbitals is referred to as the Hartree-level of theory.

When screened hydrogenic atomic orbitals are used to approximate the \(\phi_j\) and \(\phi_K\) orbitals, the resultant \(\epsilon_J\) values do **not **produce accurate predictions. For example, the negative of \(\epsilon_J\) should approximate the ionization energy for removal of an electron from the orbitals \(\phi_j\). Such ionization potentials (IP s) can be measured, and the measured values do not agree well with the theoretical values when a crude screening approximation is made for the atomic orbitals.

### The Self-Consistant Field (SCF) Approach to the Variational Method

The Hartree Equations are nonlinear and must be solved iteratively. This is because if particles interact, that interaction must be in the Hamiltonian. So until we know where the particles are, we cannot write down the Hamiltonian, but until we know the Hamiltonian, we can’t tell where the particles are.

The idea is to solve the Schrödinger equation for an electron moving in the potential of the nucleus and all the other electrons. We start with a guess for the trial electron charge density, solve N/2 one-particle Schrödinger equations (initially identical) to obtain N electron wavefunctions. Then we construct the potential for each wavefunction from that of the nucleus and that of all the other electrons, symmetrize it, and solve the N/2 Schrödinger equations again. This method is ideal for a computer, because it is easily written as an algorithm Figure \(\PageIndex{1}\).

**Figure \(\PageIndex{1}\)**: Algorithm for Self-consistent field theory.

Although we are concerned here with atoms, the same methodology is used for molecules or even solids (with appropriate potential symmetries and boundary conditions). This is a variational method, so wherever we refer to wavefunctions, we assume that they are expanded in some appropriate basis set.

Fock improved on Hartree’s method by using proper "antisymmetrized wavefunctions" (called the Hartree-Fock method) instead of simple one-electron wavefunctions.

### Shielding and Effective Charges is better than Nothing, but not Good Enough

For an atom, one can approximate the orbitals by using the solutions of the hydrogenic Schrödinger equation. These orbitals are labeled by \(n\), \(l\), and \(m\) quantum numbers for the bound states and by \(l\) and \(m\) quantum numbers and the energy \(E\) for the continuum states. These hydrogen-like orbitals provide qualitative descriptions of orbitals of atoms with more than a single electron. By introducing the concept of screening as a way to represent the repulsive interactions among the electrons of an atom, an effective nuclear charge \(Z_{\rm eff}\) can be used in place of \(Z\) in the hydrogenic \(\psi_{n,l,m}\) and \(E_{n,l}\) formulas to generate approximate atomic orbitals to be filled by electrons in a many-electron atom. For example, in the crudest approximation of a carbon atom, the two \(1s\) electrons experience the full nuclear attraction so \(Z_{\rm eff} =6\) for them, whereas the \(2s\) and \(2p\) electrons are screened by the two \(1s\) electrons, so \(Z_{\rm eff}= 4\) for them. Within this approximation, one then occupies two \(1s\) orbitals with \(Z=6\), two \(2s\) orbitals with \(Z=4\) and two \(2p\) orbitals with \(Z=4\) in forming the full six-electron product wave function of the lowest-energy state of carbon

\[\psi(1, 2, …, 6) = \psi_{1s} \psi_{1s} \psi_{2s} … \psi_{1p}(0). \label{8.3.6}\]

However, such approximate orbitals are not sufficiently accurate to be of use in quantitative simulations of atomic and molecular structure. In particular, their energies do not properly follow the trends in atomic orbital (AO) energies that are taught in introductory chemistry classes (Figure \(\PageIndex{2}\)).

**Figure \(\PageIndex{2}\):** Energies of Atomic Orbitals as Functions of the number of electrons in Neutral Atoms.

For example, the relative energies of the \(3d\) and \(4s\) orbitals are not adequately described in a model that treats electron repulsion effects in terms of a simple screening factor. So, now it is time to examine how we can move beyond the screening model and take the electron repulsion effects, which cause the inter-electronic couplings that render the Schrödinger equation insoluble, into account in a more reliable manner.

### Contriubtors

Adapted from "Quantum States of Atoms and Molecules" by David M. Hanson, Erica Harvey, Robert Sweeney, Theresa Julia Zielinski

Jack Simons (Henry Eyring Scientist and Professor of Chemistry, U. Utah) Telluride Schools on Theoretical Chemistry

Graeme Ackland (University of Edinburgh)