For the hydrogen atom, the exact wave function is known. For helium and lithium, very accurate wave functions have been calculated by including interelectronic distances in the variation functions. For atoms of higher atomic number, one way to find an accurate wave function is to first find an approximate wave function using the Hartree-Fock procedure, which we shall outline in this section. The Hartree-Fock method is the basis for the use of atomic and molecular orbitals in many-electron systems.
The Hamiltonian operator for an $n$-electron atom is
\( \begin{equation} \hat{H}=-\frac{\hbar^{2}}{2 m{e}} \sum{i=1}^{n} \nabla{i}^{2}-\sum{i=1}^{n} \frac{Z e^{2}}{4 \pi \varepsilon{0} r{i}}+\sum{i=1}^{n-1} \sum{j=i+1}^{n} \frac{e^{2}}{4 \pi \varepsilon{0} r{i j}} \tag{11.1} \end{equation} \)
where an infinitely heavy point nucleus was assumed (Section 6.6). The first sum in (11.1) contains the kinetic-energy operators for the $n$ electrons. The second sum is the potential energy (6.58) for the attractions between the electrons and the nucleus of charge $Z e$. For a neutral atom, $Z=n$. The last sum is the potential energy of the interelectronic repulsions. The restriction $j>i$ avoids counting each interelectronic repulsion twice and avoids terms like $e^{2} / 4 \pi \varepsilon{0} r{i i}$. The Hamiltonian (11.1) is incomplete, because it omits spin-orbit and other interactions. The omitted terms are small (except for atoms with high $Z$ ) and will be considered in Sections 11.6 and 11.7.
The Hartree SCF Method
Because of the interelectronic repulsion terms $e^{2} / 4 \pi \varepsilon{0} r{i j}$, the Schrödinger equation for an atom is not separable. Recalling the perturbation treatment of helium (Section 9.3), we can obtain a zeroth-order wave function by neglecting these repulsions. The Schrödinger equation would then separate into $n$ one-electron hydrogenlike equations. The zeroth-order wave function would be a product of $n$ hydrogenlike (one-electron) orbitals:
\( \begin{equation} \psi^{(0)}=f{1}\left(r{1}, \theta{1}, \phi{1}\right) f{2}\left(r{2}, \theta{2}, \phi{2}\right) \cdots f{n}\left(r{n}, \theta{n}, \phi{n}\right) \tag{11.2} \end{equation} \)
where the hydrogenlike orbitals are
\( \begin{equation} f=R{n l}(r) Y{l}^{m}(\theta, \phi) \tag{11.3} \end{equation} \)
For the ground state of the atom, we would feed two electrons with opposite spin into each of the lowest orbitals, in accord with the Pauli exclusion principle, giving the ground-state configuration. Although the approximate wave function (11.2) is qualitatively useful, it is gravely lacking in quantitative accuracy. For one thing, all the orbitals use the full nuclear charge $Z$. Recalling our variational treatments of helium and lithium, we know we can get a better approximation by using different effective atomic numbers for the different orbitals to account for screening of electrons. The use of effective atomic numbers gives considerable improvement, but we are still far from having an accurate wave function. The next step is to use a variation function that has the same form as (11.2) but is not restricted to hydrogenlike or any other particular form of orbitals. Thus we take
\( \begin{equation} \phi=g{1}\left(r{1}, \theta{1}, \phi{1}\right) g{2}\left(r{2}, \theta{2}, \phi{2}\right) \cdots g{n}\left(r{n}, \theta{n}, \phi{n}\right) \tag{11.4} \end{equation} \)
and we look for the functions $g{1}, g{2}, \ldots, g{n}$ that minimize the variational integral $\langle\phi| \hat{H}|\phi\rangle /\langle\phi \mid \phi\rangle$. Our task is harder than in previous variational calculations, where we guessed a trial function that included some parameters and then varied the parameters. In (11.4) we must vary the functions $g{i}$. [After we have found the best possible functions $g_{i}$, Eq. (11.4) will still be only an approximate wave function. The many-electron Schrödinger equation is not separable, so the true wave function cannot be written as the product of $n$ one-electron functions.]
To simplify matters somewhat, we approximate the best possible atomic orbitals with orbitals that are the product of a radial factor and a spherical harmonic:
\( \begin{equation} g{i}=h{i}\left(r{i}\right) Y{l{i}}^{m{i}}\left(\theta{i}, \phi{i}\right) \tag{11.5} \end{equation} \)
This approximation is generally made in atomic calculations. The procedure for finding the functions $g_{i}$ was introduced by Hartree in 1928 and is called the Hartree self-consistent-field (SCF) method. Hartree arrived at the SCF procedure by intuitive physical arguments. The proof that Hartree's procedure gives the best possible variation function of the form (11.4) was given by Slater and by Fock in 1930. [For the proof and a review of the SCF method, see S. M. Blinder, Am. J. Phys., 33, 431 (1965).]
Hartree's procedure is as follows. We first guess a product wave function
\( \begin{equation} \phi{0}=s{1}\left(r{1}, \theta{1}, \phi{1}\right) s{2}\left(r{2}, \theta{2}, \phi{2}\right) \cdots s{n}\left(r{n}, \theta{n}, \phi_{n}\right) \tag{11.6} \end{equation} \)
where each $s{i}$ is a normalized function of $r$ multiplied by a spherical harmonic. A reasonable guess for $\phi{0}$ would be a product of hydrogenlike orbitals with effective atomic numbers. For the function (11.6), the probability density of electron $i$ is $\left|s{i}\right|^{2}$. We now focus attention on electron 1 and regard electrons $2,3, \ldots, n$ as being smeared out to form a fixed distribution of electric charge through which electron 1 moves. We are thus averaging out the instantaneous interactions between electron 1 and the other electrons. The potential energy of interaction between point charges $Q{1}$ and $Q{2}$ is $V{12}=Q{1} Q{2} / 4 \pi \varepsilon{0} r{12}$ [Eq. (6.58)]. We now take $Q{2}$ and smear it out into a continuous charge distribution such that $\rho{2}$ is the charge density, the charge per unit volume. The infinitesimal charge in the infinitesimal volume $d v{2}$ is $\rho{2} d v{2}$, and summing up the interactions between $Q{1}$ and the infinitesimal elements of charge, we have
\( V{12}=\frac{Q{1}}{4 \pi \varepsilon{0}} \int \frac{\rho{2}}{r{12}} d v{2} \)
For electron 2 (with charge $-e$ ), the charge density of the hypothetical charge cloud is $\rho{2}=-e\left|s{2}\right|^{2}$, and for electron $1, Q_{1}=-e$. Hence
\( V{12}=\frac{e^{2}}{4 \pi \varepsilon{0}} \int \frac{\left|s{2}\right|^{2}}{r{12}} d v_{2} \)
Adding in the interactions with the other electrons, we have
\( V{12}+V{13}+\cdots+V{1 n}=\sum{j=2}^{n} \frac{e^{2}}{4 \pi \varepsilon{0}} \int \frac{\left|s{j}\right|^{2}}{r{1 j}} d v{j} \)
The potential energy of interaction between electron 1 and the other electrons and the nucleus is then
\( \begin{equation} V{1}\left(r{1}, \theta{1}, \phi{1}\right)=\sum{j=2}^{n} \frac{e^{2}}{4 \pi \varepsilon{0}} \int \frac{\left|s{j}\right|^{2}}{r{1 j}} d v{j}-\frac{Z e^{2}}{4 \pi \varepsilon{0} r_{1}} \tag{11.7} \end{equation} \)
We now make a further approximation beyond assuming the wave function to be a product of one-electron orbitals. We assume that the effective potential acting on an electron in an atom can be adequately approximated by a function of $r$ only. This central-field approximation can be shown to be generally accurate. We therefore average $V{1}\left(r{1}, \theta{1}, \phi{1}\right)$ over the angles to arrive at a potential energy that depends only on $r_{1}$ :
\( \begin{equation} V{1}\left(r{1}\right)=\frac{\int{0}^{2 \pi} \int{0}^{\pi} V{1}\left(r{1}, \theta{1}, \phi{1}\right) \sin \theta{1} d \theta{1} d \phi{1}}{\int{0}^{2 \pi} \int_{0}^{\pi} \sin \theta d \theta d \phi} \tag{11.8} \end{equation} \)
We now use $V{1}\left(r{1}\right)$ as the potential energy in a one-electron Schrödinger equation,
\( \begin{equation} \left[-\frac{\hbar^{2}}{2 m{e}} \nabla{1}^{2}+V{1}\left(r{1}\right)\right] t{1}(1)=\varepsilon{1} t_{1}(1) \tag{11.9} \end{equation} \)
and solve for $t{1}(1)$, which will be an improved orbital for electron 1 . In (11.9), $\varepsilon{1}$ is the energy of the orbital of electron 1 at this stage of the approximation. Since the potential energy in (11.9) is spherically symmetric, the angular factor in $t{1}(1)$ is a spherical harmonic involving quantum numbers $l{1}$ and $m{1}$ (Section 6.1). The radial factor $R{1}(1)$ in $t{1}$ is the solution of a one-dimensional Schrödinger equation of the form (6.17). We get a set of solutions $R{1}(1)$, where the number of nodes $k$ interior to the boundary points ( $r=0$ and $\infty$ ) starts at zero for the lowest energy and increases by 1 for each higher energy (Section 4.2). We now define the quantum number $n$ as $n \equiv l+1+k$, where $k=0,1,2, \ldots$. We thus have $1 s$, $2 s, 2 p$, and so on, orbitals (with orbital energy $\varepsilon$ increasing with $n$ ) just as in hydrogenlike atoms, and the number of interior radial nodes $(n-l-1)$ is the same as in hydrogenlike atoms (Section 6.6). However, since $V{1}\left(r{1}\right)$ is not a simple Coulomb potential, the radial factor $R{1}\left(r{1}\right)$ is not a hydrogenlike function. Of the set of solutions $R{1}\left(r{1}\right)$, we take the one that corresponds to the orbital we are improving. For example, if electron 1 is a $1 s$ electron in the beryllium $1 s^{2} 2 s^{2}$ configuration, then $V{1}\left(r{1}\right)$ is calculated from the guessed orbitals of one $1 s$ electron and two $2 s$ electrons, and we use the radial solution of (11.9) with $k=0$ to find an improved $1 s$ orbital.
We now go to electron 2 and regard it as moving in a charge cloud of density
\( -e\left[\left|t{1}(1)\right|^{2}+\left|s{3}(3)\right|^{2}+\left|s{4}(4)\right|^{2}+\cdots+\left|s{n}(n)\right|^{2}\right] \)
due to the other electrons. We calculate an effective potential energy $V{2}\left(r{2}\right)$ and solve a one-electron Schrödinger equation for electron 2 to get an improved orbital $t_{2}(2)$. We continue this process until we have a set of improved orbitals for all $n$ electrons. Then we go back to electron 1 and repeat the process. We continue to calculate improved orbitals until there is no further change from one iteration to the next. The final set of orbitals gives the Hartree self-consistent-field wave function.
How do we get the energy of the atom in the SCF approximation? It seems natural to take the sum of the orbital energies of the electrons, $\varepsilon{1}+\varepsilon{2}+\cdots+\varepsilon{n}$, but this is wrong. In calculating the orbital energy $\varepsilon{1}$, we iteratively solved the one-electron Schrödinger equation (11.9). The potential energy in (11.9) includes, in an average way, the energy of the repulsions between electrons 1 and 2,1 and $3, \ldots, 1$ and $n$. When we solve for $\varepsilon{2}$, we solve a one-electron Schrödinger equation whose potential energy includes repulsions between electrons 2 and 1,2 and $3, \ldots, 2$ and $n$. If we take $\sum{i} \varepsilon_{i}$, we will count each interelectronic repulsion twice. To correctly obtain the total energy $E$ of the atom, we must take
\( \begin{align} E & =\sum{i=1}^{n} \varepsilon{i}-\sum{i=1}^{n-1} \sum{j=i+1}^{n} \iint \frac{e^{2}\left|g{i}(i)\right|^{2}\left|g{j}(j)\right|^{2}}{4 \pi \varepsilon{0} r{i j}} d v{i} d v{j} \tag{11.10}\ E & =\sum{i} \varepsilon{i}-\sum{i} \sum{j>i} J_{i j} \end{align} \)
where the average repulsions of the electrons in the Hartree orbitals of (11.4) were subtracted from the sum of the orbital energies, and where the notation $J_{i j}$ was used for Coulomb integrals [Eq. (9.99)].
The set of orbitals belonging to a given principal quantum number $n$ constitutes a shell. The $n=1,2,3, \ldots$ shells are the $K, L, M, \ldots$ shells, respectively. The orbitals belonging to a given $n$ and a given $l$ constitute a subshell. Consider the sum of the Hartree probability densities for the electrons in a filled subshell. Using (11.5), we have
\( \begin{equation} 2 \sum{m=-l}^{l}\left|h{n, l}(r)\right|^{2}\left|Y{l}^{m}(\theta, \phi)\right|^{2}=2\left|h{n, l}(r)\right|^{2} \sum{m=-l}^{l}\left|Y{l}^{m}(\theta, \phi)\right|^{2} \tag{11.11} \end{equation} \)
where the factor 2 comes from the pair of electrons in each orbital. The spherical-harmonic addition theorem (Merzbacher, Section 9.7) shows that the sum on the right side of (11.11) equals $(2 l+1) / 4 \pi$. Hence the sum of the probability densities is $[(2 l+1) / 2 \pi]\left|h_{n, l}(r)\right|^{2}$, which is independent of the angles. A closed subshell gives a spherically symmetric probability density, a result called Unsöld's theorem. For a half-filled subshell, the factor 2 is omitted from (11.11), and here also we get a spherically symmetric probability density.
The Hartree-Fock SCF Method
The alert reader may have realized that there is something fundamentally wrong with the Hartree product wave function (11.4). Although we have paid some attention to spin and the Pauli exclusion principle by putting no more than two electrons in each spatial orbital, any approximation to the true wave function should include spin explicitly and should be antisymmetric to interchange of electrons (Chapter 10). Hence, instead of the spatial orbitals, we must use spin-orbitals and must take an antisymmetric linear combination of products of spin-orbitals. This was pointed out by Fock (and by Slater) in 1930, and an SCF calculation that uses antisymmetrized spin-orbitals is called a Hartree-Fock calculation. We have seen that a Slater determinant of spin-orbitals provides the proper antisymmetry. For example, to carry out a Hartree-Fock calculation for the lithium ground state, we start with the function (10.54), where $f$ and $g$ are guesses for the $1 s$ and $2 s$ orbitals. We then carry out the SCF iterative process until we get no further improvement in $f$ and $g$. This gives the lithium ground-state Hartree-Fock wave function.
The differential equations for finding the Hartree-Fock orbitals have the same general form as (11.9):
\( \begin{equation} \hat{F} u{i}=\varepsilon{i} u_{i}, \quad i=1,2, \ldots, n \tag{11.12} \end{equation} \)
where $u{i}$ is the $i$ th spin-orbital, the operator $\hat{F}$, called the Fock (or Hartree-Fock) operator, is the effective Hartree-Fock Hamiltonian, and the eigenvalue $\varepsilon{i}$ is the orbital energy of spin-orbital $i$. However, the Hartree-Fock operator $\hat{F}$ has extra terms as compared with the effective Hartree Hamiltonian given by the bracketed terms in (11.9). The
Hartree-Fock expression for the total energy of the atom involves exchange integrals $K_{i j}$ in addition to the Coulomb integrals that occur in the Hartree expression (11.10). See Section 14.3. [Actually, Eq. (11.12) applies only when the Hartree-Fock wave function can be written as a single Slater determinant, as it can for closed-subshell atoms and atoms with only one electron outside closed subshells. When the Hartree-Fock wave function contains more than one Slater determinant, the Hartree-Fock equations are more complicated than (11.12).]
The orbital energy $\varepsilon_{i}$ in the Hartree-Fock equations (11.12) can be shown to be a good approximation to the negative of the energy needed to ionize a closed-subshell atom by removing an electron from spin-orbital $i$ (Koopmans' theorem; Section 15.5).
Originally, Hartree-Fock atomic calculations were done by using numerical methods to solve the Hartree-Fock differential equations (11.12), and the resulting orbitals were given as tables of the radial functions for various values of $r$. [The Numerov method (Sections 4.4 and 6.9) can be used to solve the radial Hartree-Fock equations for the radial factors in the Hartree-Fock orbitals; the angular factors are spherical harmonics. See D. R. Hartree, The Calculation of Atomic Structures, Wiley, 1957; C. Froese Fischer, The Hartree-Fock Method for Atoms, Wiley, 1977.]
In 1951, Roothaan proposed representing the Hartree-Fock orbitals as linear combinations of a complete set of known functions, called basis functions. Thus for lithium we would write the Hartree-Fock $1 s$ and $2 s$ spatial orbitals as
\( \begin{equation} f=\sum{i} b{i} \chi{i}, \quad g=\sum{i} c{i} \chi{i} \tag{11.13} \end{equation} \)
where the $\chi{i}$ functions are some complete set of functions, and where the $b{i}$ 's and $c{i}$ 's are expansion coefficients that are found by the SCF iterative procedure. Since the $\chi{i}$ (chi $i$ ) functions form a complete set, these expansions are valid. The Roothaan expansion procedure allows one to find the Hartree-Fock wave function using matrix algebra (see Section 14.3 for details). The Roothaan procedure is readily implemented on a computer and is often used to find atomic Hartree-Fock wave functions and nearly always used to find molecular Hartree-Fock wave functions.
A commonly used set of basis functions for atomic Hartree-Fock calculations is the set of Slater-type orbitals (STOs) whose normalized form is
\( \begin{equation} \frac{\left(2 \zeta / a{0}\right)^{n+1 / 2}}{[(2 n)!]^{1 / 2}} r^{n-1} e^{-\zeta r / a{0}} Y_{l}^{m}(\theta, \phi) \tag{11.14} \end{equation} \)
The set of all such functions with $n, l$, and $m$ being integers obeying (6.96)-(6.98) but with $\zeta$ having all possible positive values forms a complete set. The parameter $\zeta$ is called the orbital exponent. To get a truly accurate representation of the Hartree-Fock orbitals, we would have to include an infinite number of Slater orbitals in the expansions. In practice, one can get very accurate results by using only a few judiciously chosen Slater orbitals. (Another possibility is to use Gaussian-type basis functions; see Section 15.4.)
Clementi and Roetti did Hartree-Fock calculations for the ground state and some excited states of the first 54 elements of the periodic table [E. Clementi and C. Roetti, At. Data Nucl. Data Tables, 14, 177 (1974); Bunge and co-workers have recalculated these wave functions; C. F. Bunge et al., At. Data Nucl. Data Tables, 53, 113 (1993); Phys. Rev. A, 46, 3691 (1992); these atomic wave functions can be found at www.ccl.net/cca/data/ atomic-RHF-wavefunctions/index.shtml]. For example, consider the Hartree-Fock groundstate wave function of helium, which has the form [see Eq. (10.41)]
\( f(1) f(2) \cdot 2^{-1 / 2}[\alpha(1) \beta(2)-\alpha(2) \beta(1)] \)
Clementi and Roetti expressed the $1 s$ orbital function $f$ as the following combination of five $1 s$ Slater-type orbitals:
\( f=\pi^{-1 / 2} \sum{i=1}^{5} c{i}\left(\frac{\zeta{i}}{a{0}}\right)^{3 / 2} e^{-\zeta{i} r / a{0}} \)
where the expansion coefficients $c{i}$ are $c{1}=0.76838, c{2}=0.22346, c{3}=0.04082$, $c{4}=-0.00994, c{5}=0.00230$ and where the orbital exponents $\zeta{i}$ are $\zeta{1}=1.41714$, $\zeta{2}=2.37682, \zeta{3}=4.39628, \zeta{4}=6.52699, \zeta{5}=7.94252$. [The largest term in the expansion has an orbital exponent that is similar to the orbital exponent (9.62) for the simple trial function (9.56).] The Hartree-Fock energy is -77.9 eV , as compared with the true nonrelativistic energy, -79.0 eV . The $1 s$ orbital energy corresponding to $f$ was found to be -25.0 eV , as compared with the experimental helium ionization energy of 24.6 eV .
For the lithium ground state, Clementi and Roetti used a basis set consisting of two 1 s STOs (with different orbital exponents) and four $2 s$ STOs (with different orbital exponents). The lithium $1 s$ and $2 s$ Hartree-Fock orbitals were each expressed as a linear combination of all six of these basis functions. The Hartree-Fock energy is -202.3 eV , as compared with the true energy -203.5 eV .
Electron densities calculated from Hartree-Fock wave functions are quite accurate. Figure 11.1 compares the radial distribution function of argon (found by integrating the electron density over the angles $\theta$ and $\phi$ and multiplying the result by $r^{2}$ ) calculated by the Hartree-Fock method with the experimental radial distribution function found by electron diffraction. (Recall from Section 6.6 that the radial distribution function is proportional to the probability of finding an electron in a thin spherical shell at a distance $r$ from the nucleus.) Note the electronic shell structure in Fig. 11.1. The high nuclear charge in ${ }_{18} \mathrm{Ar}$ makes the average distance of the $1 s$ electrons from the nucleus far less than in H or He . Thus there is only a moderate increase in atomic size as we go down a given group in the periodic table. Calculations show that the radius of a sphere containing $98 \%$ of the Hartree-Fock electron probability density gives an atomic radius in good agreement with the empirically determined van der Waals radius. [See C. W. Kammeyer and D. R. Whitman, J. Chem. Phys., 56, 4419 (1972).]
Although the radial distribution function of an atom shows the shell structure, the electron probability density integrated over the angles and plotted versus $r$ does not oscillate. Rather, for ground-state atoms this probability density is a maximum at the nucleus (because of the $s$ electrons) and continually decreases as $r$ increases. Similarly, in molecules the maxima in electron probability density usually occur at the nuclei; see, for example, Fig. 13.7. [For further discussion, see H. Weinstein, P. Politzer, and S. Srebnik, Theor. Chim. Acta, 38, 159 (1975).]
FIGURE 11.1 Radial distribution function in Ar as a function of $r$. The broken line is the result of a Hartree-Fock calculation. The solid line is the result of electron-diffraction data. [Reprinted figure with permission from L.S. Bartell and L. O. Brockway, Physical Review Series II, Vol 90, 833, 1953. Copyright 1953 by the American Physical Society.]
Accurate representation of a many-electron atomic orbital (AO) requires a linear combination of several Slater-type orbitals. For rough calculations, it is convenient to have simple approximations for AOs. We might use hydrogenlike orbitals with effective nuclear charges, but Slater suggested an even simpler method: to approximate an AO by a single function of the form (11.14) with the orbital exponent $\zeta$ taken as
\( \begin{equation} \zeta=(Z-s) / n \tag{11.15} \end{equation} \)
where $Z$ is the atomic number, $n$ is the orbital's principal quantum number, and $s$ is a screening constant calculated by a set of rules (see Prob. 15.62). A Slater orbital replaces the polynomial in $r$ in a hydrogenlike orbital with a single power of $r$. Hence a single Slater orbital does not have the proper number of radial nodes and does not represent well the inner part of an orbital.
A lot of computation is required to perform a Hartree-Fock SCF calculation for a many-electron atom. Hartree did several SCF calculations in the 1930s, when electronic computers were not in existence. Fortunately, Hartree's father, a retired engineer, enjoyed numerical computation as a hobby and helped his son. Nowadays computers have replaced Hartree's father.