Comprehensive Study Notes

The hydrogen molecule is the simplest molecule containing an electron-pair bond. The purely electronic Hamiltonian (13.5) for $\mathrm{H}_{2}$ is in atomic units

\(
\begin{equation}
\hat{H}=-\frac{1}{2} \nabla{1}^{2}-\frac{1}{2} \nabla{2}^{2}-\frac{1}{r{a 1}}-\frac{1}{r{a 2}}-\frac{1}{r{b 1}}-\frac{1}{r{b 2}}+\frac{1}{r_{12}} \tag{13.91}
\end{equation}
\)

where 1 and 2 are the electrons and $a$ and $b$ are the nuclei (Fig. 13.18). Just as in the helium atom, the $1 / r_{12}$ interelectronic-repulsion term prevents the Schrödinger equation from being separable. We therefore use approximation methods.

We start with the molecular-orbital approach. The ground-state electron configuration of $\mathrm{H}{2}$ is $\left(\sigma{g} 1 s\right)^{2}$, and we can write an approximate wave function as the Slater determinant

\(
\begin{align}
\frac{1}{\sqrt{2}}\left|\begin{array}{ll}
\sigma{g} 1 s(1) \alpha(1) & \sigma{g} 1 s(1) \beta(1) \
\sigma{g} 1 s(2) \alpha(2) & \sigma{g} 1 s(2) \beta(2)
\end{array}\right| & =\sigma{g} 1 s(1) \sigma{g} 1 s(2) \cdot 2^{-1 / 2}[\alpha(1) \beta(2)-\beta(1) \alpha(2)] \
& =f(1) f(2) \cdot 2^{-1 / 2}[\alpha(1) \beta(2)-\beta(1) \alpha(2)] \tag{13.92}
\end{align}
\)

FIGURE 13.18 Interparticle distances in $\mathrm{H}_{2}$.

which is similar to $(10.26)$ for the helium atom. To save time, we write $f$ instead of $\sigma_{g} 1 s$. As we saw in Section 10.4, omission of the spin factor does not affect the variational integral for a two-electron problem. Hence we want to choose $f$ so as to minimize

\(
\frac{\iint f^{}(1) f^{}(2) \hat{H} f(1) f(2) d v{1} d v{2}}{\iint|f(1)|^{2}|f(2)|^{2} d v{1} d v{2}}
\)

where the integration is over the spatial coordinates of the two electrons. Ideally, $f$ should be found by an SCF calculation. For simplicity we can use an $\mathrm{H}{2}^{+}$-like MO. (The $\mathrm{H}{2}$ Hamiltonian becomes the sum of two $\mathrm{H}{2}^{+}$Hamiltonians if we omit the $1 / r{12}$ term.) We saw in Section 13.5 that the function [Eq. (13.57)]

\(
\frac{k^{3 / 2}}{(2 \pi)^{1 / 2}\left(1+S{a b}\right)^{1 / 2}}\left(e^{-k r{a}}+e^{-k r_{b}}\right)
\)

gives a good approximation to the ground-state $\mathrm{H}{2}^{+}$wave function. Hence we try as a variation function $\phi$ for $\mathrm{H}{2}$ the product of two such LCAO functions, one for each electron:

\(
\begin{gather}
\phi=\frac{\zeta^{3}}{2 \pi\left(1+S{a b}\right)}\left(e^{-\zeta r{a 1}}+e^{-\zeta r{b 1}}\right)\left(e^{-\zeta r{a 2}}+e^{-\zeta r{b 2}}\right) \tag{13.93}\
\phi=\frac{1}{2\left(1+S{a b}\right)}\left[1 s{a}(1)+1 s{b}(1)\right]\left[1 s{a}(2)+1 s{b}(2)\right] \tag{13.94}
\end{gather}
\)

where the effective nuclear charge $\zeta$ will differ from $k$ for $\mathrm{H}_{2}^{+}$. Since

\(
\hat{H}=\hat{H}{1}^{0}+\hat{H}{2}^{0}+1 / r_{12}
\)

where $\hat{H}{1}^{0}$ and $\hat{H}{2}^{0}$ are $\mathrm{H}_{2}^{+}$Hamiltonians for each electron, we have

\(
\iint \phi^{*} \hat{H} \phi d v{1} d v{2}=2 W{1}+\iint \frac{\phi^{2}}{r{12}} d v{1} d v{2}
\)

where $W{1}$ is given by (13.63) with $k$ replaced by $\zeta$. The evaluation of the $1 / r{12}$ integral is complicated and is omitted [see Slater, Quantum Theory of Molecules and Solids, Volume 1, page 65, and Appendix 6]. Coulson performed the variational calculation in 1937, using (13.93). [For the literature references of the $\mathrm{H}{2}$ calculations mentioned in this and later sections, see the bibliography in A. D. McLean et al., Rev. Mod. Phys., 32, 211 (1960).] Coulson found $R{e}=0.732 \AA$, which is close to the true value $0.741 \AA$; the minimum in the calculated $U(R)$ curve gave $D{e}=3.49 \mathrm{eV}$, as compared with the true value
4.75 eV (Table 13.2). (Of course, the percent error in the total electronic energy is much less than the percent error in $D{e}$, but $D{e}$ is the quantity of chemical interest.) The value of $\zeta$ at $0.732 \AA$ is 1.197 , which is less than $k$ for $\mathrm{H}{2}^{+}$. We attribute this to the screening of the nuclei from each electron by the other electron.

How can we improve on the above simple MO result? We can look for the best possible MO function $f$ in (13.92) to get the Hartree-Fock wave function for $\mathrm{H}_{2}$. This was done by Kolos and Roothaan [W. Kolos and C. C. J. Roothaan, Rev. Mod. Phys., 32, 219 (1960)]. They expanded $f$ in elliptic coordinates [Eq. (13.34)]. Since $m=0$ for the ground state, the $e^{i m \phi}$ factor in the SCF MO is equal to 1 and $f$ is a function of $\xi$ and $\eta$ only. The expansion used is

\(
f=e^{-\alpha \xi} \sum{p, q} a{p q} \xi^{p} \eta^{q}
\)

where $p$ and $q$ are integers and $\alpha$ and $a{p q}$ are variational parameters. The Hartree-Fock results are $R{e}=0.732 \AA$ and $D{e}=3.64 \mathrm{eV}$, which is not much improvement over the value 3.49 eV given by the simple LCAO molecular orbital. The correlation energy for $\mathrm{H}{2}$ is thus -1.11 eV , close to the value -1.14 eV for the two-electron helium atom (Section 11.3). To get a truly accurate binding energy, we must go beyond the SCF approximation of writing the wave function in the form $f(1) f(2)$. We can use the same methods we used for atoms: configuration interaction and introduction of $r_{12}$ into the trial function.

First, consider configuration interaction (CI). To reach the exact ground-state wave function, we include contributions from SCF (or other) functions for all the excited states with the same symmetry as the ground state. In the first approximation, only contributions from the lowest-lying excited states are included. The first excited configuration of $\mathrm{H}{2}$ is $\left(\sigma{g} 1 s\right)\left(\sigma{u}^{*} 1 s\right)$, which gives the terms ${ }^{1} \Sigma{u}^{+}$and ${ }^{3} \Sigma{u}^{+}$. (We have one $g$ and one $u$ electron, so the terms are of odd parity.) The ground-state configuration $\left(\sigma{g} 1 s\right)^{2}$ is a ${ }^{1} \Sigma{g}^{+}$state. Hence we do not get any contribution from the $\left(\sigma{g} 1 s\right)\left(\sigma{u}^{*} 1 s\right)$ states, since they have different parity from the ground state. Next consider the configuration $\left(\sigma{u}^{} 1 s\right)^{2}$. This is a closed-shell configuration having the single state ${ }^{1} \Sigma{g}^{+}$. This is of the right symmetry to contribute to the ground-state wave function. As a simple CI trial function, we can take a linear combination of the MO wave functions for the $\left(\sigma{g} 1 s\right)^{2}$ and $\left(\sigma_{u}^{} 1 s\right)^{2}$ configurations. To simplify things, we will use the LCAO-MOs as approximations to the MOs. Thus we take

\(
\begin{equation}
\phi=\sigma{g} 1 s(1) \sigma{g} 1 s(2)+c \sigma_{u}^{} 1 s(1) \sigma_{u}^{} 1 s(2) \tag{13.95}
\end{equation}
\)

where $\sigma{g} 1 s$ and $\sigma{u}^{*} 1 s$ are given by (13.57) and (13.58) with a variable orbital exponent and $c$ is a variational parameter. This calculation was performed by Weinbaum in 1933. The result is a bond length of $0.757 \AA$ and a dissociation energy of 4.03 eV , which is a considerable improvement over the Hartree-Fock result $D{e}=3.64 \mathrm{eV}$. The orbital exponent has the optimum value 1.19 . We can improve on this result by using a better form for the MOs of each configuration and by including more configuration functions. Hagstrom did a CI calculation in which the MOs were represented by expansions in elliptic coordinates. With 33 configuration functions, he found $D{e}=4.71 \mathrm{eV}$, close to the true value 4.75 eV [S. Hagstrom and H. Shull, Rev. Mod. Phys., 35, 624 (1963)].

Now consider the use of $r{12}$ in $\mathrm{H}{2}$ trial functions. The first really accurate calculation of the hydrogen-molecule ground state was done by James and Coolidge in 1933. They used the trial function

\(
\exp \left[-\delta\left(\xi{1}+\xi{2}\right)\right] \sum c{m n j k p}\left[\xi{1}^{m} \xi{2}^{n} \eta{1}^{j} \eta{2}^{k}+\xi{1}^{n} \xi{2}^{m} \eta{1}^{k} \eta{2}^{j}\right] r{12}^{p}
\)

where the summation is over integral values of $m, n, j, k$, and $p$. The variational parameters are $\delta$ and the $c{m n j k p}$ coefficients. The James and Coolidge function is symmetric with
respect to interchange of electrons 1 and 2, as it should be, since we have an antisymmetric ground-state spin function. With 13 terms in the sum, James and Coolidge found $D{e}=4.72 \mathrm{eV}$, only 0.03 eV in error. Their work has been extended by Kolos, Wolniewicz, and co-workers, who used as many as 279 terms in the sum. Since it is $D{0}$ that is determined from the observed electronic spectrum, they used the Cooley-Numerov method (Section 13.2) to calculate the vibrational levels from their theoretical $U(R)$ curve and then calculated $D{0}$. Including relativistic corrections and corrections to the Born-Oppenheimer approximation, they found $D{0} / h c=36118.1 \mathrm{~cm}^{-1}$, in agreement with the spectroscopically determined value $36118.1 \mathrm{~cm}^{-1}$ [W. Kolos et al., J. Chem. Phys., 84, 3278 (1986); L. Wolniewicz, J. Chem. Phys., 99, 1851 (1993)]. An even more precise $D{0}$ was calculated by workers who did high-precision calculations of relativistic corrections, corrections to the Born-Oppenheimer approximation, and quantumelectrodynamics corrections using a Born-Oppenheimer $U(R)$ found from wave functions having as many as 7000 terms to get $D_{0} / h c=36118.0695(10) \mathrm{cm}^{-1}$, where the number in parentheses is the estimate of the uncertainty of the last digits [K. Piszczatowski et al., J. Chem. Theory Comput., 5, 3039 (2009); www.fuw.edu.pl/~krp/papers/D0.pdf]. The experimental value is $36118.0696(4) \mathrm{cm}^{-1}$. (Quantum electrodynamics is the relativistic quantum-mechanical theory of the interaction of radiation and matter formulated by Feynman, Schwinger, and Tomonaga in the late 1940s and is an improvement on the quantum field theory of Dirac.)