13: The Born-Oppenheimer approximation and the hydrogen molecule ion


Using quantum mechanics to predict the chemical bonding patterns, optimal geometries, and physical and chemical properties of molecules is a large and active field of research known as molecular quantum mechanics or more commonly as quantum chemistry. The density functional theory referred to in the previous lecture, for which the chemistry Nobel prize was given in 1998, has had a tremendous impact in quantum chemistry, with some of the papers in this subject having acquired some 10,000 citations each since their publication. In fact, the 1998 chemistry Nobel prize was shared between Walter Kohn, one of the inventors of density functional theory and John Pople, the developer of a commonly used quantum chemistry software package.

Quantum chemistry calculations allow the geometries of molecules to be computed as well as a wide range of properties. Quantum chemistry can also be used in a novel way, in which the electrons are treated using quantum mechanics but the nuclei are treated as classical particles. We use quantum mechanics to calculate the internuclear forces but then use these forces in Newton's Second Law to study the motion of the nuclei during chemical reactions. This gives us a microscopic window into the specific motions, the complex dance, executed by the nuclei during a simple or complex chemical process.

The methods of quantum chemistry have become very sophisticated, and there are various software packages that can be downloaded for carrying out the calculations of quantum chemistry. It should be noted that these packages use a series of approximations to solve the Schrödinger equation because for all but the simplest of molecules, exact solutions are not available. We will discuss some of these methods, but first we need to introduce some of the underlying theory.

The Born-Oppenheimer Approximation

In this section, we will discuss one of the most important and fundamental approximations in molecular quantum mechanics. This approximation was developed by Max Born and J. Robert Oppenheimer in 1927. We will consider a very general molecule with $$N$$ nuclei and $$M$$ electrons. The coordinates of the nuclei are $$R_1 ,...,R_N$$. The coordinates of the electrons are $$r_1 ,...,r_M$$, and their spin variables are $$S_{z,1},...,S_{z,M}$$. For shorthand, we will denote the complete set of nuclear coordinates as $$R$$ and the set of electrons coordinates as $$r$$ and $$x$$ the complete set of electron coordinates and spin variables. The total molecular wave function $$\Psi (x,R)$$ depends on $$3N+4M$$ variables, which makes it a very cumbersome object to deal with. The Born-Oppenheimer approximation leads to a very important simplification of the wave function.

If we could neglect the electron-nuclear interaction, then the wave function would be a simple product $$\Psi (x,R)=\psi_{elec}(r)\psi_{nucl}(R)$$. However, we cannot neglect this term, but it might still be possible to write the wave function as a product. We note, first, that most nuclei are 3-4 orders of magnitude heaver than an electron. For the lightest nucleus, the proton,

$m_p \approx 2000 m_e$

This mass difference is large enough to have important physical consequences. Let us think classically about this mass difference first. If two particles interact in some way, and one is much heavier than the other, the light particle will move essentially as a slave'' of the heavy particle. That is, it will simply follow the heavy particle wherever it goes, and, it will move rapidly in response to the heavy particle motion. As an illustration of this phenomenon, consider the simple mechanical system pictured below:

Considering this as a classical system, we expect that the motion will be dominated by the large heavy particle, which is attached to a fixed wall by a spring. The small, light particle, which is attached to the heavy particle by a spring will simply follow the heavy particle and execute rapid oscillations around it. The figure below (bottom panel) illustrates this:

In this illustration, when the heavy particle moves even a tiny amount, the light particle executes many oscillations around the heavy particle. Thus, we see that the light particle moves quite a bit, and we can think of the light particle as generating on the fly'' an effective potential $$V_{eff}$$ for the heavy particle, at least at the present position of the heavy particle. What we can conclude, therefore, is that we can approximately fix the position of the heavy particle and just allow the light particle to move around before we advance the heavy particle to its next position.

How do the different mass scales of electrons and nuclei manifest themselves in the time-independent Schrödinger equation? Light particles, such as electrons, tend to have very diffuse, delocalized wave functions, while heavy particles, such as the nuclei, tend to have wave functions that are very localized about the classical positions. This is illustrated in the figure below:

Thus, taking the idea of solving the electronic problem for fixed nuclear positions seriously, we can write the molecular wave function as

$\Psi^{(mol)}(x,R) = \psi^{(elec)}(x,R)\psi^{(nuc)}(R)$

which suggests that the electronic wave function $$\psi^(elec)(x,R)$$ is solved using the nuclear positions $$R$$ simply as parameters that characterize the wave function in the same way that the electron mass and electron charge do. This is not an exact wave function for the molecule but an approximate one. Hence, we call this the Born-Oppenheimer approximation to the wave function. The total classical energy of the molecule is

$K_e + K_n + V_{ee}(r) + V_{en}(r,R) + V_{nn}(R) = E$

where $$K_e$$ and $$K_n$$ are the electronic and nuclear kinetic energies. Therefore, the total Hamiltonian of the molecule is

$\hat{H} = \hat{K}_e + \hat{K}_n + V_{ee}(r) + V_{en}(r,R) + V_{nn}(R)$

where $$\hat{K}_e$$ and $$\hat{K}_n$$ are the kinetic energy operators that result from substituting momenta for derivatives. Thus, $$\hat{K}_e$$ and $$\hat{K}_n$$ contain second derivatives with respect to electronic and nuclear coordinates, respectively. The details of these operators are not that important for this discussion.

We now write the total molecular Hamiltonian as

$\hat{H} = \hat{H}_{elec} + \hat{H}_{nucl}$

where we assign

\begin{align*}\hat{H}_{elec} &= \hat{K}_e +V_{ee}(r)+V_{en}(r,R)\\ \hat{H}_{nucl} &= \hat{K}_n +V_{nn}(R)\end{align*}

If we now substitute this into the Schrödinger equation, we obtain

$[\hat{K}_e +\hat{K}_n +V_{ee}(r)+V_{nn}(R)+V_{en}(r,R)]\psi^{(elec)}(x,R)\psi^{(nuc)}(R)=E\psi^{(elec)}(x,R)\psi^{(nuc)}(R)$

Let us consider how the electronic and nuclear kinetic energy operators act on the wave function. $$\hat{K}_e$$ contains derivatives with respect only to electronic coordinates. Thus, it has no effect on the nuclear wave function, and we can write

$\hat{K}_e\psi^{(elec)}(x,R)\psi^{(nuc)}(R) =\psi^{(nuc)}(R)\hat{K}_e\psi^{(elec)}(x,R)$

The operator $$\hat{K}_n$$ contains derivatives with respect to nuclear coordinates, thus it affects both parts of the wave function:

$\hat{K}_n \psi^{(elec)}(x,R)\psi^{(nuc)}(R)=\psi^{(nuc)}(R)\hat{K}_n \psi^{(elec)}(x,R)+\psi^{(elec)}(x,R)\hat{K}_n \psi^{(nuc)}(R)+mixed \ derivate \ term$

The last term arises from an application of the product rule using the fact that $$\hat{K}_n \sim d^2/dR^2$$. What is important here is that if the nuclear wave function is highly localized, then its curvature is high in the vicnity of the nucleus, and hence, $$\hat{K}_n\psi^{(nuc)}$$ is very large compared to $$\hat{K}_n \psi^{(elec)}$$ since $$\psi^{(elec)}$$ changes much less dramatically spatially due to the delocalized nature of the electronic wave function. Thus, we need only keep the term in the above equation that contains $$\hat{K}_n \psi^{(nuc)}$$:

$\hat{K}_n \psi^{(elec)}(x,R)\psi^{(nuc)}(R) \approx \psi^{(elec)}(x,R)\hat{K}_n \psi^{(nuc)}(R)$

Using this fact in the Schrödinger equation, we can write the equation as

$\psi^{(nuc)}(R)[\hat{K}_e +V_{ee}(r)+V_{en}(r,R)]\psi^{(elec)}(x,R)+\psi^{(elec)}(x,R)[\hat{K}_n +V_{nn}(R)]\psi^{(nuc)}(R)=E\psi^{(elec)}(x,R)\psi^{(nuc)}(R)$

or

$\psi^{(nuc)}(R)\hat{H}_{elec}(r,R)\psi^{(elec)}(x,R)+\psi^{(elec)}(x,R)\hat{H}_{nuc}(R)\psi^{(nuc)}(R)=E\psi^{(elec)}(x,R)\psi^{(nuc)}(R)$

We now divide both sides by $$\psi^{(elec)}(x,R)\psi^{(nuc)}(R)$$, which gives

$\frac{\hat{H}_{elec}(r,R)\psi^{(elec)}(x,R)}{\psi^{(elec)}(x,R)}=E-\dfrac{\hat{H}_{nuc}\psi^{(nuc)}(R)}{\psi^{(nuc)}(R)}$

Note that the right side depends only on the nuclear coordinates $$R$$, which means that we can write compactly as some function $$\varepsilon (R)$$. Hence, we obtain

$\frac{\hat{H}_{elec}(r,R)\psi^{(elec)}(x,R)}{\psi^{(elec)}(x,R)}=\varepsilon (R)$

or

$\hat{H}_{elec}(r,R)\psi^{(elec)}(x,R)=\varepsilon (R)\psi^{(elec)}(x,R)$

which is called the electronic Schrödinger equation. it yields a set of wave functions $$\psi_{\alpha}^{(elec)}(x,R)$$ and energies $$\varepsilon_\alpha (R)$$ characterized by a set of quantum numbers , which is the full set of quantum numbers needed to characterize the wave function, e.g. $$n,l,m$$ for hydrogen. Thus, the electronic Schrödinger equation can be written as

$\hat{H}_{elec}(r,R)\psi_{\alpha}^{(elec)}(x,R)=\varepsilon_\alpha (R)\psi_{\alpha}^{(elec)}(x,R)$

Note that the potential energy $$V_{en}(r,R)$$ depends on the electronic coordinates and the fixed nuclear coordinates $$R$$. Thus, the nuclear coordinates can be thought of as additional parameters in the potential energy. Alternatively, we can think of the nuclear as providing a fixed background potential energy in addition to the internal electron-electron repulsion $$V_{ee}(r)$$ that governs the shape of the electronic wave function.

Now the nuclear part comes from setting the right side of the full Schrödinger equation equal to $$\varepsilon_\alpha (R)$$, which yields:

$E-\frac{\hat{H}_{nuc}\psi^{(nuc)}(R)}{\psi^{(nuc)}(R)}=\varepsilon_\alpha (R)$

and if we multiply by $$\psi^{(nuc)}(R)$$, we obtain the nuclear Schrödinger equation

\begin{align*}[\hat{H}_{nucl}+\varepsilon_\alpha (R)]\psi^{(nuc)}(R) &= E\psi^{(nuc)}(R)\\ [\hat{K}_n +V_{nn}(R)+\varepsilon_\alpha (R)]\psi^{(nuc)}(R) &= E\psi^{(nuc)}(R)\end{align*}

The Schrödinger equations for the electronic and nuclear part of the wave function are not exact but approximate based on the separation of electronic and nuclear degrees of freedom we assume can be made based on the mass differential.

Two interesting facts about the approximate Schrödinger equations should be noted. First, the electronic Schrödinger equation yields a set of energy levels $$\varepsilon_\alpha (R)$$ that depend on the nuclear configuration. Thus, at each nuclear configuration we choose, we get a different set of energy levels. These energy levels change continuously with the nuclear configuration, in fact. At the same time, the nuclear Schrödinger equation involves a potential energy with two terms $$V_{nn}(R)+\varepsilon_\alpha (R)$$, i.e. the nuclear-nuclear repulsion plus the electronic energy level $$\varepsilon_\alpha$$ as a function of $$R$$. Therefore, there is a different nuclear Schrödinger equation for each electronic energy level. The potential energies $$V_{nn}(R)+\varepsilon_\alpha (R)$$ are called Born-Oppenheimer potential energies or Born-Oppenheimer surfaces. They are surfaces (or hypersurfaces) because they must be viewed as continuous functions of $$R$$ (Figure $$\PageIndex{1}$$):

Figure $$\PageIndex{1}$$​​​​​​: Potentials energies surfaces

The potential energy $$V_{nn}(R)+\varepsilon_\alpha (R)$$ determines the geometry of the molecule by minimization, and the nuclear probability distribution by solution of the nuclear Schrödinger equation. It also determines the nuclear dynamics if we solve the nuclear time-dependent Schrödinger equation. Thus, $$\varepsilon_\alpha (R)$$, the electronic energy levels, are of utmost importance in molecular quantum mechanics. For this reason, we need to understand the distributions of electrons in molecules.

Incidentally, we could now introduce a further approximation in which we treat the nuclei as classical particles with a potential energy $$V_{nn}(R)+\varepsilon_\alpha (R)$$. The nuclei would then move as classical point particles via Newton's Second Law, which would now take the form

$M\frac{d^2 R}{dt^2}=-\dfrac{d}{dR}[V_{nn}(R)+\varepsilon_\alpha (R)]$

Solving The Hydrogen Molecule Cation Exactly

The hydrogen molecule ion $$H_{2}^{+}$$ is the only molecule for which we can solve the electronic Schrödinger equation exactly. Note that it has just one electron! In fact, there are no multi-electron molecules we can solve exactly. Thus from $$H_2$$ on up to more complicated molecules, we only have approximate solutions for the allowed electronic energies and wave functions. Before we discuss these, however, let us examine the exact solutions for $$H_{2}^{+}$$ starting with a brief outline of how the exact solution is carried out.

The figure below shows the geometry of the $$H_{2}^{+}$$ molecule ion and the coordinate system we will use.

The two protons are labeled A and B, and the distances from each proton to the one electron are $$r_A$$ and $$r_B$$, respectively. Let $$R$$ be the distance between the two protons (this is the only nuclear degree of freedom that is important, and the electronic wave function will depend parametrically only on $$R$$). The coordinate system is chosen so that the protons lie on the z-axis one at a distance $$R/2$$ above the $$xy$$ plane and one a distance $$-R/2$$ below the $$xy$$ plane (see figure). The classical energy of the electron is

$\frac{p^2}{2m_e}-\frac{e^2}{4\pi \epsilon_0}\left ( \frac{1}{r_A}+\frac{1}{r_B} \right )=E$

The nuclear-nuclear term $$V_{nn}(R)=e^2 /4\pi \epsilon_0 R$$ is a constant, and we can define the potential energy relative to this quantity.

The energy is not a simple function of energies for the $$x$$, $$y$$, and $$z$$ directions, so we try another coordinate system to see if we can simplify the problem. In fact, this problem has a natural cylindrical symmetry (analogous to the spherical symmetry of the hydrogen atom) about the z-axis. Thus, we try cylindrical coordinates. In cylindrical coordinates the distance of the electron from the z-axis is denoted $$\rho$$, the angle $$\phi$$ is the azimuthal angle, as in spherical coordinates, and the last coordinate is just the Cartesian $$z$$ coordinate. Thus,

$x=\rho cos\phi \;\;\;\; y=\rho sin\phi \;\;\;\; z=z$

Using right triangles, the distance $$r_A$$ and $$r_B$$ can be shown to be

$r_A =\sqrt{\rho^2 +(R/2-z)^2}\;\;\;\; r_B =\sqrt{\rho^2 +(R/2 +z)^2}$

The classical energy becomes

$\frac{p_{\rho}^{2}}{2m_e}+\frac{p_{\phi}^{2}}{2m_e \rho^2}+\frac{p_{z}^{2}}{2m_e}-\frac{e^2}{4\pi \epsilon_0}\left ( \frac{1}{\sqrt{\rho^2 +(R/2-z)^2}}+\frac{1}{\sqrt{\rho^2 +(R/2+z)^2}} \right ) =E$

First, we note that the potential energy does not depend on $$\phi$$, and the classical energy can be written as a sum $$E=\varepsilon_{\rho ,z}+\varepsilon_\phi$$ of a $$\rho$$ and $$z$$ dependent term and an angular term. Moreover, angular momentum is conserved as it is in the hydrogen atom. However, in this case, only one component of the angular momentum is actually conserved, and this is the z-component, since the motion is completely symmetric about the z-axis. Thus, we have only certain allowed values of $$L_z$$, which are $$m\bar{h}$$, as in the hydrogen atom, where $$m=0,\pm 1,\pm 2,...$$. The electronic wave function (now dropping the "(elec)" label, since it is understood that we are discussing only the electronic wave function), can be written as a product

$\psi (\rho ,\phi ,z)=G(\rho ,z)y(\phi)$

and $$y(\phi )$$ is given by

$y_m (\phi )=\frac{1}{\sqrt{2\pi}}e^{im\phi}$

which satisfies the required boundary condition $$y_m (0)=y_m (2\pi )$$. Unfortunately, what is left in $$\rho$$ and $$z$$ is still not that simple. But if we make one more change of coordinates, the problem simplifies. We introduce two new coordinates $$\mu$$ and $$\nu$$ defined by

$\mu =\frac{r_A +r_B}{R}\;\;\;\; \nu =\frac{r_A -r_B}{R}$

Note that when $$\nu =0$$, the electron is in the $$xy$$ plane. Thus, $$\nu$$ is analogous to $$z$$ in that it varies most as the electron moves along the z-axis. The presence or absence of a node in the $$xy$$ plane will be an important indicator of wave functions that lead to a chemical bond in the molecule or not. The coordinate $$\mu$$, on the other hand,is minimum when the electron is on the z-axis and grows as the distance of the electron from the z-axis increases. Thus, $$mu$$ is analogous to $$\rho$$. The advantage of these coordinates is that the wave function turns out to be a simple product

$\psi (\mu ,\nu ,\phi ,R)=M(\mu )N(\nu )y(\phi )$

which greatly simplifies the problem and allows the exact solution.

The mathematical structure of the exact solutions is complex and nontransparent, so we will only look at these graphically, where we can gain considerably insight. First, we note that the quantum number $$m$$ largely determines how the solutions appear. First, let us introduce the nomenclature for designating the orbitals (solutions of the Schrödinger equation, wave functions) of the system

1. If $$m=0$$, the orbitals are called $$\sigma$$ orbitals, analogous to the $$s$$ orbitals in hydrogen.
2. If $$m=1$$, the orbitals are called $$\pi$$ orbitals, analogous to the $$p$$ orbitals in hydrogen.
3. If $$m=2$$, the orbitals are called $$\delta$$ orbitals, analogous to the $$d$$ orbitals in hydrogen.
4. If $$m=3$$, the orbitals are called $$\phi$$ orbitals, analogous to the $$f$$ orbitals in hydrogen.

These orbitals are known as molecular orbitals because the describe the electronic wave functions of an entire molecule. There are four designators that we use to express each molecular orbital:

1. A greek letter, $$\sigma$$, $$\pi$$, $$\delta$$, $$\phi$$, .... depending on the quantum number $$m$$.
2. A subscript qualifier g or u depending on how an orbital $$\psi$$ behaves with respect to a spatial reflection or parity operation $$r\rightarrow -r$$. If $\psi (-r)=\psi (r)$ then $$\psi (r)$$ is an even function of $$r$$, so we use the $$"g"$$ designator, where $$g$$ stands for the German word gerade, meaning "even". If $\psi (-r)=-\psi (r)$ then $$\psi (r)$$ is an odd function of $$r$$, and we use the $$"u"$$ designator, where $$u$$ stands for the German word ungerade, meaning "odd".
3. An integer $$n$$ in front of the Greek letter to designate the energy level. This is analogous to the integer we use in atomic orbitals $$(1s, 2s, 2p,...)$$.
4. An asterisk or no asterisk depending on the presence or absence of nodes between the nuclei. If there is significant amplitude between the nuclei, then the orbital favors a chemical bond, and the orbital is called a bonding orbital. If there is a node between the nuclei, the orbital does not favor bonding, and the orbital is called an antibonding orbital.

So, the first few orbitals, in order of increasing energy are:

$1\sigma_g ,1\sigma_{u}^{*},2\sigma_g ,2\sigma_{u}^{*}, 1\pi_u ,3\sigma_g ,1\pi_{g}^{*},3\sigma_u$

These orbitals are depicted in the figure below:

The $$1\sigma_g$$ and $$1\sigma_{u}^{*}$$ orbitals are the lowest in energy, however, note that the $$1\sigma_{u}^{*}$$ contains one more node than the $$1\sigma_g$$ orbital, hence it has a higher energy. Similarly for the $$2\sigma_g$$ and $$2\sigma_{u}^{*}$$ orbitals. The former has two nodes while the latter has three and, therefore, it is of higher energy. The next set of orbitals are displayed in the right panel. In this set of orbitals, the $$1\pi_u$$ is the lowest energy with a single node. The number of nodes increases as we go up in energy in this subset of orbitals. In addition, in all of the orbitals, the nodal structure reveals the bonding/antibonding character of each orbital. When there is significant amplitude between the nuclei, the orbital is a bonding orbital, otherwise, when there is a node there, it is an antibonding orbital.

What do bonding and antibonding orbitals mean in terms of the corresponding energy levels. Consider just the first two energy levels $$\varepsilon_0 (R)$$ and $$\varepsilon_1 (R)$$ corresponding to the $$1\sigma_g$$ and $$1\sigma_{u}^{*}$$ orbitals. The ground-state orbital $$1\sigma_g$$ is bonding and the first excited state $$1\sigma_{u}^{*}$$ is antibonding. In the figure below, we plot the energy levels $$\varepsilon_0 (R)$$ and $$\varepsilon_1 (R)$$ as functions of the internuclear separation $$R$$:

The picture reveals that the curve $$\varepsilon_0 (R)$$ (the red curve) has a well-defined minimum, corresponding to the equilibrium bond length. This is due to the bonding character of the orbital. On the other hand, the first excited state $$\varepsilon_1 (R)$$ (the blue curve) has no such minimum. There is no stable bond length in this orbital, which reflects the antibonding character. Thus, exciting the molecule into this electronic state causes it to dissociate. The energy needed to do this depends on the internuclear separation $$R$$. When $$R=R_e$$, the equilibrium bond length, an energy of several $$Ry$$ would be needed. At larger separations, the energy decreases as the two curves approach each other.