The Hydrogen Molecular Ion and the Chemical Bond

Derivation

The Hydrogen Molecular Ion and the Chemical Bond

D-311 Home PU-306 Threads energy · matter · symmetry Depends on The Born-Oppenheimer Approximation, variational-principle

Statement

For the one-electron molecular ion \(\mathrm{H_2^+}\) in the Born–Oppenheimer approximation, we build the electronic wavefunction as a linear combination of the two atomic \(1s\) orbitals, \(\psi = c_A\phi_A + c_B\phi_B\). The variational principle and the homonuclear symmetry force \(c_A = \pm c_B\), yielding a symmetric (bonding) state \(\psi_+\) and an antisymmetric (antibonding) state \(\psi_-\) with energies \(E_\pm(R) = E_{1s} + \dfrac{1}{R} - \dfrac{J \pm K}{1 \pm S}\) (atomic units). The bonding curve \(E_+(R)\) develops a minimum below the separated-atom energy: electron sharing lowers the energy and binds the nuclei.

Why it matters

\(\mathrm{H_2^+}\) is the simplest molecule — two nuclei, one electron — and the only molecule whose electronic Schrödinger equation is exactly separable (in confocal elliptic coordinates). It is therefore the benchmark against which every approximate theory of chemical bonding is calibrated. The LCAO treatment here is the first rung of molecular-orbital theory: the same bonding/antibonding splitting, overlap integral, and resonance integral reappear in every diatomic and, generalised, in band structure of solids.

Physically, the derivation answers a foundational question: why do atoms stick together? The naive electrostatic picture (electron "between" the nuclei) is quantitatively wrong; the true mechanism is quantum — constructive interference of the two atomic amplitudes builds electron density in the internuclear region and, more subtly, relaxes the electron's kinetic-energy confinement. LCAO makes both effects computable.

Assumptions

Born–Oppenheimer separation holds.If dropped, the nuclei cannot be clamped at fixed \(R\); the electronic energy \(E_\pm(R)\) is no longer a well-defined potential-energy surface and vibronic coupling mixes electronic states. The basis is exactly two atomic \(1s\) orbitals.If dropped (adding \(2p_z\), polarisation, or a floating basis), the variational energy drops further; the two-term result is a rigorous upper bound, not the exact curve. The atomic orbitals carry unit effective nuclear charge, \(\zeta = 1\).If dropped and \(\zeta\) is variationally optimised, the well deepens from \(1.76\,\mathrm{eV}\) toward \(2.35\,\mathrm{eV}\) and \(R_e\) contracts toward the exact value. The trial function is real and the Hamiltonian Hermitian, so \(H_{AB}=H_{BA}\) and \(S_{AB}=S_{BA}\).If dropped (complex orbitals or a non-Hermitian model potential), the secular problem need not have real eigenvalues and the simple \((\alpha\pm\beta)/(1\pm S)\) factorisation fails. The electronic ground state is nondegenerate and the two centres are identical.If dropped (heteronuclear, \(H_{AA}\neq H_{BB}\)), the eigenvectors are no longer the symmetric/antisymmetric combinations and \(c_A\neq\pm c_B\).

Derivation

\[ \hat{H}_{\mathrm{el}} = -\tfrac{1}{2}\nabla^2 - \frac{1}{r_A} - \frac{1}{r_B}, \qquad E(R) = \frac{\langle \psi|\hat{H}_{\mathrm{el}}|\psi\rangle}{\langle \psi|\psi\rangle} + \frac{1}{R} \]

In Born–Oppenheimer, clamp the protons at separation \(R\); the electron sees two fixed Coulomb wells (atomic units: \(\hbar=m_e=e=4\pi\varepsilon_0=1\)). The constant nuclear repulsion \(1/R\) adds directly to the electronic expectation value. A

\[ \psi(\mathbf{r}) = c_A\,\phi_A(\mathbf{r}) + c_B\,\phi_B(\mathbf{r}), \qquad \phi_{A,B} = \frac{1}{\sqrt{\pi}}\,e^{-r_{A,B}} \]

LCAO ansatz: near nucleus \(A\) the electron should look like a hydrogen \(1s\) on \(A\), and likewise for \(B\). The \(c_{A,B}\) are variational parameters. A

\[ \frac{\partial E}{\partial c_A}=0,\quad \frac{\partial E}{\partial c_B}=0 \;\Longrightarrow\; \sum_{n} \big(H_{mn} - E\,S_{mn}\big)\,c_n = 0 \]

By the variational principle the best coefficients minimise the Rayleigh quotient. Differentiating \(E\,\langle\psi|\psi\rangle=\langle\psi|\hat H|\psi\rangle\) and setting derivatives to zero gives the linear (secular) equations, with \(H_{mn}=\langle\phi_m|\hat H|\phi_n\rangle\), \(S_{mn}=\langle\phi_m|\phi_n\rangle\). B

\[ \det\!\begin{pmatrix} H_{AA}-E & H_{AB}-E S \\ H_{AB}-E S & H_{AA}-E \end{pmatrix} = 0 \]

Non-trivial \((c_A,c_B)\) exist only if the secular determinant vanishes. Homonuclear symmetry sets \(H_{AA}=H_{BB}\), \(S_{AA}=S_{BB}=1\), and Hermiticity sets \(H_{AB}=H_{BA}\), \(S_{AB}=S_{BA}=S\). B

\[ (H_{AA}-E)^2 = (H_{AB}-ES)^2 \;\Longrightarrow\; H_{AA}-E = \pm\,(H_{AB}-ES) \]

The determinant is a difference of squares. Taking the square root gives two branches. This is pure algebra on the \(2\times2\) determinant. A

\[ E_\pm = \frac{H_{AA} \pm H_{AB}}{1 \pm S} \equiv \frac{\alpha \pm \beta}{1 \pm S} \]

Solve each branch for \(E\). Define the Coulomb (site) integral \(\alpha \equiv H_{AA}\) and the resonance (hopping) integral \(\beta \equiv H_{AB}\). Since \(\beta<0\), the \(+\) branch lies lowest. A

\[ \psi_+ = \frac{\phi_A+\phi_B}{\sqrt{2(1+S)}}, \qquad \psi_- = \frac{\phi_A-\phi_B}{\sqrt{2(1-S)}} \]

Back-substituting each eigenvalue gives \(c_A=+c_B\) for \(E_+\) and \(c_A=-c_B\) for \(E_-\); the prefactors follow from \(\langle\psi_\pm|\psi_\pm\rangle=1\) using \(\langle\phi_A|\phi_B\rangle=S\). B

\[ \alpha = E_{1s} - J,\qquad \beta = E_{1s}\,S - K \]

Because \(\big(-\tfrac12\nabla^2 - 1/r_A\big)\phi_A = E_{1s}\phi_A\) with \(E_{1s}=-\tfrac12\), the extra term is the attraction of the electron on \(A\) to the other nucleus: \(J=\langle\phi_A|1/r_B|\phi_A\rangle\) (Coulomb) and \(K=\langle\phi_A|1/r_B|\phi_B\rangle\) (exchange/resonance). C

\[ S = e^{-R}\!\left(1+R+\tfrac{R^2}{3}\right),\quad J = \frac{1}{R}-e^{-2R}\!\left(1+\frac{1}{R}\right),\quad K = e^{-R}\!\left(1+R\right) \]

Closed forms obtained by evaluating the two-centre integrals in confocal elliptic (prolate spheroidal) coordinates, where the exponentials separate. These are exact for hydrogenic \(1s\) orbitals with \(\zeta=1\). C

\[ E_\pm(R) = \frac{\alpha\pm\beta}{1\pm S} + \frac{1}{R} = E_{1s} - \frac{J \pm K}{1 \pm S} + \frac{1}{R} \]

Insert Step 8 into Step 6 and add the nuclear repulsion from Step 1. The factor \((1\pm S)^{-1}\) cancels the \(E_{1s}(1\pm S)\) piece cleanly, leaving the compact potential-energy curves. B

Result

\[ \boxed{\,E_\pm(R) = E_{1s} + \frac{1}{R} - \frac{J(R) \pm K(R)}{1 \pm S(R)}\,},\qquad E_{1s}=-\tfrac12\ \text{hartree} \]

Reading. The bonding curve \(E_+(R)\) uses \(+\) throughout: the resonance integral \(K>0\) enters as \(-(J+K)/(1+S)\), a large lowering, and produces a minimum at \(R_e\approx 2.5\,a_0\) that lies \(D_e\approx1.76\,\mathrm{eV}\) below the separated \(\mathrm{H}+\mathrm{H^+}\) limit. The antibonding curve \(E_-(R)\) uses \(-\); the resonance term now raises the energy and the curve is purely repulsive, with a nodal plane midway between the nuclei where \(\psi_-=0\). Bonding piles electron amplitude between the nuclei; antibonding evacuates it.

Units check. In atomic units every term is an energy in hartree: \(E_{1s}\) and \(1/R\) (with \(R\) in bohr) are energies; \(J,K\) arise from \(\langle 1/r\rangle\) which has dimension inverse length = energy in atomic units; \(S\) is dimensionless, so \((J\pm K)/(1\pm S)\) is an energy. Restoring SI, multiply by the hartree \(E_h = m_e e^4/(4\pi\varepsilon_0)^2\hbar^2 = 27.2114\,\mathrm{eV}\) and measure \(R\) in bohr \(a_0=0.52918\,\text{Å}\).

Limiting cases

Separated atoms \(R\to\infty\): \(S,J,K\to0\) and \(1/R\to0\), so \(E_\pm\to E_{1s}=-\tfrac12\) hartree — a lone H atom plus a bare proton, as it must.
United atom \(R\to0\): the two protons merge into a \(Z=2\) nucleus; \(E_+\) should correlate to the He\(^+\) \(1s\) level \(-2\) hartree, but the frozen \(\zeta=1\) LCAO fails here — the \(1/R\) repulsion diverges and dominates.
Weak overlap \(S\ll1\): \(E_\pm\approx E_{1s}+1/R-(J\pm K)\), so the splitting is symmetric, \(E_\pm - \bar E \approx \mp K\), set purely by the resonance integral.
Optimised \(\zeta\): replacing \(\zeta=1\) by the variational \(\zeta\approx1.24\) at equilibrium contracts the orbitals, deepening \(D_e\) to \(\approx2.35\,\mathrm{eV}\) and shifting \(R_e\) to \(\approx2.0\,a_0\).

Breaks when

Small \(R\) (united-atom limit). The fixed \(\zeta=1\) minimal basis cannot describe the electron collapsing onto a \(Z=2\) centre; the LCAO curve turns over incorrectly and misses the He\(^+\) correlation. Only a flexible or scaled basis recovers the right small-\(R\) behaviour.
Breakdown of Born–Oppenheimer. Near curve crossings or for the loosely bound highest vibrational levels, nuclear kinetic energy is not negligible; \(E_\pm(R)\) ceases to be a rigid potential and non-adiabatic coupling terms must be restored.
Many-electron molecules. For anything beyond one electron the derivation omits electron–electron repulsion and exchange entirely; a two-orbital LCAO of \(\mathrm{H_2}\) already needs the Coulomb/exchange two-electron integrals and, at large \(R\), dissociates to the wrong (ionic-contaminated) products.
Strong external fields or relativistic regime. A large field mixes in high-\(\ell\) character absent from the \(1s\)-only basis; spin–orbit and relativistic corrections, negligible for hydrogen, are simply not in this Hamiltonian.

Failure modes

Forgetting the \(1/R\) nuclear repulsion. Students plot \(\alpha\pm\beta\) only; without \(+1/R\) the "bonding" curve has no minimum and dives to \(-\infty\).
Confusing overlap \(S\) with the resonance integral \(K\). \(S=\langle\phi_A|\phi_B\rangle\) is dimensionless; \(K=\langle\phi_A|1/r_B|\phi_B\rangle\) is an energy. They are different integrals.
Sign error on \(\beta\). Treating \(\beta>0\) and assigning \(E_-\) as bonding. The bonding (nodeless, symmetric) state is the lower one, \(E_+\).
Normalising as if orbitals were orthogonal. Writing \(\psi_+=(\phi_A+\phi_B)/\sqrt2\) and dropping the \((1+S)\). The cross term \(2S\) is exactly what makes bonding work.
"Electrostatic glue" misreading. Claiming the bond energy is just the electron sitting between two protons; the dominant stabilisation in \(\mathrm{H_2^+}\) is kinetic (delocalisation), not classical electrostatics.
Using \(E_{1s}=-13.6\,\mathrm{eV}\) inside atomic-unit formulas. Mixing eV numbers into hartree expressions; keep one unit system until the final conversion.

Discussion

The bonding state lowers the energy for two cooperating reasons. Constructive interference of \(\phi_A\) and \(\phi_B\) builds electron density in the internuclear region — the term \(2\phi_A\phi_B\) in \(|\psi_+|^2\) — which the electron feels as extra Coulomb attraction to both nuclei simultaneously. This is the origin of the negative resonance integral contribution. The antibonding state has the opposite sign of interference, a node between the nuclei, and depleted density there, so it is destabilised.

The two states track the exact \(\mathrm{H_2^+}\) eigenstates, which are labelled by their parity under inversion: \(\psi_+\) is \(1\sigma_g\) (gerade, even), \(\psi_-\) is \(1\sigma_u^*\) (ungerade, odd). This gerade/ungerade symmetry is exact for any homonuclear diatomic and is the reason \(c_A=\pm c_B\) emerged without solving anything — group theory fixes the eigenvectors, the variational calculation only fixes the energies. The same one-electron reasoning, iterated over a lattice, becomes the tight-binding model of solids: \(\alpha\) is the on-site energy and \(\beta\) the hopping amplitude that broadens levels into bands.

A subtler and initially counter-intuitive point, established by Hellmann, Ruedenberg, and later Kutzelnigg from careful energy decomposition: the driving force of covalent bonding in \(\mathrm{H_2^+}\) is a lowering of the electron's kinetic energy, not a lowering of potential energy. Delocalisation over two centres relaxes the confinement of the wavefunction along the bond axis, reducing \(\langle p_z^2\rangle\). The naive LCAO with frozen \(\zeta=1\) partly obscures this because it violates the virial theorem \(2\langle T\rangle = -\langle V\rangle\) at every \(R\); only after variationally scaling \(\zeta(R)\) — which contracts the orbitals and restores the virial ratio — does the kinetic-energy lowering appear cleanly. This is why the scaled calculation both deepens the well and shortens the bond.

Common misconceptions. (i) The bond is not "the electron glued electrostatically between the protons" — a purely classical charge in the midplane is a saddle, not a stable configuration, and the real stabilisation is quantum-mechanical. (ii) The overlap \(S\) is not a small correction to be dropped; it is precisely what distinguishes the two energies and normalisations. (iii) "Antibonding" does not mean "unoccupied" — it is a genuine excited electronic state, purely repulsive, that governs photodissociation of \(\mathrm{H_2^+}\).

Worked examples

Example 1 — Bonding binding energy at \(R=2.5\,a_0\).

\[ S = e^{-R}\!\left(1+R+\tfrac{R^2}{3}\right) = e^{-2.5}\!\left(1+2.5+\tfrac{6.25}{3}\right) = 0.082085\times5.5833 = 0.4583 \]

Evaluate the overlap integral at \(R=2.5\) bohr. A

\[ J = \frac{1}{R}-e^{-2R}\!\left(1+\frac1R\right) = 0.4 - e^{-5}(1.4) = 0.4-0.009433 = 0.39057 \]

Coulomb integral (hartree). A

\[ K = e^{-R}(1+R) = e^{-2.5}(3.5) = 0.082085\times3.5 = 0.28730 \]

Resonance integral (hartree). A

\[ E_+ = E_{1s} + \frac1R - \frac{J+K}{1+S} = -0.5 + 0.4 - \frac{0.67787}{1.4583} = -0.5648\ \text{hartree} \]

Insert into the boxed result with the \(+\) branch. B

\[ D = E_{1s} - E_+ = -0.5 -(-0.5648) = 0.0648\ \text{hartree} = 0.0648\times27.2114\ \mathrm{eV} \]

Binding energy = separated-atom energy minus molecule energy. A

\[ D_e \approx 1.76\ \mathrm{eV},\qquad R_e \approx 2.5\,a_0 = 1.32\ \text{Å} \]

Reading. The simple LCAO binds \(\mathrm{H_2^+}\) by \(1.76\,\mathrm{eV}\) at a bond length of \(1.32\) Å, versus the exact \(2.79\,\mathrm{eV}\) at \(1.06\) Å — right physics, quantitatively soft because of the frozen minimal basis.

Units check. hartree \(\times\,27.2114\,\mathrm{eV\,hartree^{-1}} = \mathrm{eV}\); bohr \(\times\,0.52918\,\text{Å}\,a_0^{-1} = \text{Å}\).

Example 2 — The antibonding state is repulsive at the same geometry.

\[ 1-S = 1-0.4583 = 0.5417,\qquad J-K = 0.39057-0.28730 = 0.10327 \]

Reuse the \(R=2.5\) integrals from Example 1; the \(-\) branch needs the differences. A

\[ E_- = E_{1s} + \frac1R - \frac{J-K}{1-S} = -0.5 + 0.4 - \frac{0.10327}{0.5417} = -0.2907\ \text{hartree} \]

Boxed result, \(-\) branch. B

\[ E_- - E_{1s} = -0.2907-(-0.5) = +0.2094\ \text{hartree} = +0.2094\times27.2114\ \mathrm{eV} \]

Energy relative to the separated-atom asymptote; a positive value means no binding. A

\[ E_- - E_{1s} \approx +5.70\ \mathrm{eV}\ >\ 0 \]

Reading. At the same \(R\) where the bonding state is \(1.76\,\mathrm{eV}\) below threshold, the antibonding state sits \(5.7\,\mathrm{eV}\) above it — strongly repulsive, with a node in the internuclear plane. Occupying \(\psi_-\) would break the bond, exactly the channel accessed in \(\mathrm{H_2^+}\) photodissociation.

Units check. Difference of two hartree energies converted with \(27.2114\,\mathrm{eV\,hartree^{-1}}\); sign convention: positive = above the \(\mathrm{H}+\mathrm{H^+}\) limit.

Problems

Compute the overlap integral \(S\) for two hydrogen \(1s\) orbitals (\(\zeta=1\)) at \(R=2.0\,a_0\).
Solution
\(S=e^{-2}\left(1+2+\tfrac{4}{3}\right)=0.135335\times4.3333=0.5865.\) Overlap is substantial (nearly 0.6) at this separation, which is why the \((1+S)\) normalisation cannot be neglected.
Show that normalising \(\psi_+ = N(\phi_A+\phi_B)\) gives \(N = [2(1+S)]^{-1/2}\).
Solution
\(\langle\psi_+|\psi_+\rangle = N^2\big(\langle\phi_A|\phi_A\rangle + \langle\phi_B|\phi_B\rangle + 2\langle\phi_A|\phi_B\rangle\big) = N^2(1+1+2S)=2N^2(1+S).\) Setting this equal to 1 gives \(N=1/\sqrt{2(1+S)}\). For \(\psi_-\) the cross term changes sign, giving \(1/\sqrt{2(1-S)}\).
Show analytically that \(E_\pm(R)\to E_{1s}\) as \(R\to\infty\), and identify the slowest-decaying correction to \(E_+\).
Solution
As \(R\to\infty\): \(S,K\sim e^{-R}\to0\), \(J\to1/R\), and \(1/R\to0\). Then \(E_+\to E_{1s}+1/R - J = E_{1s}+1/R-1/R = E_{1s}.\) The \(1/R\) from nuclear repulsion is exactly cancelled by the leading term of \(J\) (screening by the electron). The residual is set by the exponentially small \(K\) and the \(-e^{-2R}(1+1/R)\) tail of \(J\); the correction decays as \(\sim e^{-R}\), not as a power law — a signature of exchange/resonance binding rather than a van der Waals tail.
At \(R=2.0\,a_0\) compute \(E_+\), \(E_-\), and the bonding–antibonding splitting \(\Delta=E_--E_+\) in eV. Use \(S=0.5865\), and evaluate \(J,K\) yourself.
Solution
\(J=0.5-e^{-4}(1.5)=0.5-0.018316\times1.5=0.47253.\) \(K=e^{-2}(3)=0.135335\times3=0.40600.\) Bonding: \(E_+=-0.5+0.5-\tfrac{0.87853}{1.5865}=-0.55374\) hartree. Antibonding: \(E_-=-0.5+0.5-\tfrac{J-K}{1-S}=0-\tfrac{0.06653}{0.4135}=-0.16090\) hartree. Splitting \(\Delta = -0.16090-(-0.55374)=0.39284\) hartree \(=10.69\,\mathrm{eV}.\) The splitting is large and grows as the atoms approach.
The variational optimisation of the orbital exponent (Finkelstein–Horowitz) gives \(\zeta\approx1.238\) at equilibrium. Explain qualitatively why \(\zeta>1\) improves the energy, and state the resulting \(R_e\) and \(D_e\).
Solution
A bonding electron is attracted by two protons, so it experiences an effective nuclear charge greater than one; contracting the \(1s\) orbitals (\(\zeta>1\)) lets them concentrate density in the internuclear region and lowers the energy. By the scaling \(r\to\zeta r\) the integrals become functions of \(\rho=\zeta R\), and minimising \(E(\zeta,R)\) jointly yields \(\zeta\approx1.238\), \(R_e\approx2.00\,a_0\ (1.06\ \text{Å})\), and \(D_e\approx0.0865\) hartree \(\approx2.35\,\mathrm{eV}\) — much closer to the exact \(2.79\,\mathrm{eV}\), and it restores the virial ratio \(2\langle T\rangle=-\langle V\rangle\) at the minimum.