The many-electron problem & Born–Oppenheimer

Everything qc-rs computes is, at bottom, an approximate solution of one equation — the Schrödinger equation for a molecule. This chapter writes that equation down, explains why it cannot be solved exactly, and takes the first essential step toward making it tractable: separating the motion of the nuclei from that of the electrons.

A molecule, quantum-mechanically

A molecule is a collection of positively charged nuclei and negatively charged electrons, held together by electrostatic (Coulomb) forces. In quantum mechanics its stationary states are described by a wavefunction \(\Psi\) and obey the time-independent Schrödinger equation

\[ \hat H\,\Psi = E\,\Psi , \]

where \(\hat H\) is the Hamiltonian operator (the total energy) and \(E\) is the energy of the state. Solve this and you know, in principle, everything: the energy, the structure, and every property.

Note

Atomic units Quantum chemistry works in atomic units, chosen so that the electron mass, elementary charge, \(\hbar\), and the Coulomb constant are all \(1\). Two you will see constantly: energies are in hartree (\(1\ E_h \approx 27.211\ \text{eV} \approx 627.5\ \text{kcal/mol}\)) and lengths in bohr (\(1\ a_0 \approx 0.529\ \text{Å}\)). This is why the quickstart energy was a plain number like -76.026772 — it is in hartree. When you pass unit="angstrom", qc-rs converts your coordinates to bohr internally.

The molecular Hamiltonian

For \(M\) nuclei (labelled \(A,B,\dots\), with charges \(Z_A\) and masses \(M_A\)) and \(N\) electrons (labelled \(i,j,\dots\)), the Hamiltonian in atomic units is a sum of five physically distinct terms:

\[ \hat H = \underbrace{-\sum_{A}\frac{1}{2M_A}\nabla_A^2}_{\text{nuclear kinetic}} \;\underbrace{-\sum_{i}\frac{1}{2}\nabla_i^2}_{\text{electron kinetic}} \;\underbrace{-\sum_{A,i}\frac{Z_A}{r_{Ai}}}_{\text{electron–nucleus attraction}} \;\underbrace{+\sum_{i<j}\frac{1}{r_{ij}}}_{\text{electron–electron repulsion}} \;\underbrace{+\sum_{A<B}\frac{Z_A Z_B}{R_{AB}}}_{\text{nucleus–nucleus repulsion}} , \]

where \(r_{Ai}=|\mathbf r_i-\mathbf R_A|\) is an electron–nucleus distance, \(r_{ij}=|\mathbf r_i-\mathbf r_j|\) an electron–electron distance, and \(R_{AB}\) a nucleus–nucleus distance. Every term is either kinetic energy or a Coulomb interaction — nothing exotic. The difficulty is entirely in how they are coupled.

Why it cannot be solved exactly

The villain is the electron–electron repulsion term \(\sum_{i<j} 1/r_{ij}\). It ties the coordinates of every electron to those of every other, so the \(N\)-electron wavefunction \(\Psi(\mathbf r_1,\mathbf r_2,\dots,\mathbf r_N)\) does not separate into a product of one-electron functions. This coupling — that each electron’s best move depends on where all the others are — is called electron correlation, and it is the reason there is no closed-form solution for anything past the one-electron hydrogen atom. Worse, the wavefunction lives in \(3N\) spatial dimensions: even storing it on a grid is hopeless for more than a few electrons. The whole enterprise of quantum chemistry is the search for good, controlled approximations.

Electrons carry one more essential rule. They are fermions, so the wavefunction must be antisymmetric: swapping the full coordinates (space and spin) of any two electrons flips its sign,

\[ \Psi(\dots,\mathbf x_i,\dots,\mathbf x_j,\dots) = -\,\Psi(\dots,\mathbf x_j,\dots,\mathbf x_i,\dots) . \]

This is the Pauli principle; the next chapters build it in from the start (via Slater determinants).

The Born–Oppenheimer approximation

The first great simplification exploits a huge disparity: a nucleus is at least ~1800 times heavier than an electron, so electrons move far faster and effectively “see” the nuclei as stationary. We may therefore clamp the nuclei in place and solve for the electrons alone.

Definition 1 (Born–Oppenheimer approximation)

Treat the nuclear positions \(\{\mathbf R_A\}\) as fixed parameters. The electronic Hamiltonian drops the nuclear kinetic term and treats the constant nuclear repulsion separately,

\[ \hat H_{\text{elec}} = -\sum_{i}\frac{1}{2}\nabla_i^2 - \sum_{A,i}\frac{Z_A}{r_{Ai}} + \sum_{i<j}\frac{1}{r_{ij}} , \]

and one solves the electronic Schrödinger equation for a fixed geometry,

\[ \hat H_{\text{elec}}\,\Psi_{\text{elec}} = E_{\text{elec}}\,\Psi_{\text{elec}} . \]

The total energy at that geometry is the electronic energy plus the (now constant) nuclear repulsion,

\[ E(\mathbf R) = E_{\text{elec}}(\mathbf R) + \sum_{A<B}\frac{Z_A Z_B}{R_{AB}} . \]

Two ideas fall out of this that you already met in Part I:

• The potential energy surface (PES). Because \(E(\mathbf R)\) depends on the fixed geometry, sweeping the nuclei traces out an energy landscape. A minimum of this surface is an equilibrium structure — exactly what the geometry-optimization tutorial walked downhill to find. Its slope is the force on the nuclei (the gradient).

• Nuclear repulsion is just an additive constant. For a single geometry the \(\sum_{A<B} Z_A Z_B/R_{AB}\) term is a number added at the end. This is why the total energy the quickstart printed already includes it — the electronic problem is where all the difficulty lives.

Where this leaves us

Born–Oppenheimer reduces “solve a molecule” to “solve the electrons for a fixed set of nuclei,” but the electronic Schrödinger equation is still an intractable many-body problem, thanks to correlation and the \(3N\) dimensions. The next chapter shows the standard escape: stop trying to find \(\Psi\) exactly, and instead expand it in a finite basis and minimize the energy — turning calculus into linear algebra.