Variational principle, LCAO & basis sets#

The previous chapter left us with the electronic Schrödinger equation — exact, but unsolvable. This chapter introduces the strategy that makes quantum chemistry computable: stop searching for the exact wavefunction, restrict it to a form we can handle, and minimize the energy within that form. Two ingredients make this concrete — the variational principle (why minimizing works) and basis sets (the form we restrict to). Together they turn a differential equation into a matrix problem a computer can solve.

The variational principle#

Suppose we cannot find the true ground-state wavefunction \(\Psi_0\), so we guess a trial wavefunction \(\Psi\) instead. How good is our guess? The variational principle gives a remarkable answer: the energy of any trial wavefunction is an upper bound to the true ground-state energy.

Theorem 1 (Variational principle)

For any well-behaved trial wavefunction \(\Psi\), the Rayleigh quotient

\[ E[\Psi] = \frac{\langle\Psi|\hat H|\Psi\rangle}{\langle\Psi|\Psi\rangle} \]

satisfies \(E[\Psi] \ge E_0\), where \(E_0\) is the exact ground-state energy, with equality iff \(\Psi = \Psi_0\).

This single fact reshapes the whole problem. Because a lower energy always means a better wavefunction, we do not need to solve the Schrödinger equation directly — we can parameterize \(\Psi\) with some adjustable numbers and minimize \(E[\Psi]\) over them. Calculus (solving a differential equation) becomes optimization (finding the lowest energy). Every method in this manual — Hartree–Fock, DFT, and the correlated methods — is a particular choice of trial form plus this minimization.

Molecular orbitals and the LCAO idea#

The parameters have to come from somewhere. The organizing idea is the molecular orbital (MO): a one-electron function \(\psi_i(\mathbf r)\) that describes how a single electron is spread over the whole molecule. We will build the many-electron wavefunction out of MOs (the next chapter makes this precise with a Slater determinant); for now the question is just how to represent one MO.

We cannot store a general 3-D function exactly, so we expand each MO in a fixed set of simple functions — the basis set. The functions \(\{\phi_\mu\}\) are atom-centered “atomic orbitals,” and an MO is a Linear Combination of Atomic Orbitals (LCAO):

\[ \psi_i(\mathbf r) = \sum_{\mu=1}^{K} C_{\mu i}\,\phi_\mu(\mathbf r) . \]

Now the unknowns are just the numbers \(C_{\mu i}\) — the MO coefficients. Minimizing the energy over them is a finite problem: an MO is fixed once you know its \(K\) coefficients.

What the basis functions are#

In practice the atomic orbitals \(\phi_\mu\) are built from Gaussian-type orbitals (GTOs). A primitive Gaussian centered on atom \(A\) has the form

\[ g(\mathbf r) = N\,x_A^{\,a} y_A^{\,b} z_A^{\,c}\; e^{-\alpha\,|\mathbf r - \mathbf A|^2}, \]

where \((x_A,y_A,z_A)=\mathbf r-\mathbf A\), the exponents \(a,b,c\) set the angular shape (their sum is the angular momentum: \(0=s\), \(1=p\), \(2=d\), …), \(\alpha\) controls the width, and \(N\) normalizes it. Real basis functions are contracted — fixed linear combinations of several primitives that mimic an atomic orbital,

\[ \phi_\mu(\mathbf r) = \sum_{k} d_{k\mu}\,g_k(\mathbf r), \]

with the contraction coefficients \(d_{k\mu}\) baked into the basis-set definition.

Note

Why Gaussians? A hydrogen-like orbital actually decays like \(e^{-\zeta r}\) (a Slater function), which Gaussians reproduce only imperfectly (they are too flat at the nucleus and fall off too fast). Gaussians win anyway for one decisive reason: the product of two Gaussians on different centers is a single Gaussian on a point between them. That identity makes the millions of two-electron integrals a molecule needs analytically fast — which is exactly why qc-rs (via the libcint library) uses them.

Basis-set families#

Bigger, more flexible basis sets describe the electrons better — at higher cost. The vocabulary you will meet in ao=:

  • Minimal (e.g. sto-3g) — one function per occupied atomic orbital. Cheap, qualitative.

  • Split-valence / multiple-zeta (e.g. cc-pvdz = double-zeta, cc-pvtz = triple-zeta) — two, three, … functions per valence orbital, letting orbitals shrink or expand in different bonding environments.

  • Polarization functions — higher angular momentum than the atom needs in isolation (a \(d\) on oxygen, a \(p\) on hydrogen), letting orbitals distort along bonds. The “P” in cc-pVDZ is these.

  • Diffuse functions — extra broad functions (a aug- prefix) for anions, lone pairs, and long-range interactions.

The correlation-consistent family cc-pVXZ (\(X=\) D, T, Q, …) is designed to approach the exact answer smoothly as \(X\) grows.

From coefficients to a matrix problem#

Insert the LCAO expansion into the Rayleigh quotient and minimize over the coefficients \(C_{\mu i}\). Because the atomic orbitals are not orthogonal (neighbouring atoms’ functions overlap), the result is a generalized matrix eigenvalue problem,

\[ \mathbf{F}\,\mathbf{C} = \mathbf{S}\,\mathbf{C}\,\boldsymbol{\varepsilon}, \qquad S_{\mu\nu} = \langle \phi_\mu | \phi_\nu \rangle , \]

where \(\mathbf S\) is the overlap matrix, \(\mathbf F\) encodes the energy, \(\mathbf C\) holds the MO coefficients, and \(\boldsymbol\varepsilon\) the orbital energies. The exact form of \(\mathbf F\) is the subject of the next chapter; the point here is that the impossible differential equation has become a \(K\times K\) matrix problem, and \(K\) — the number of basis functions — is what you choose when you pick ao=.

Basis-set convergence, in numbers#

Here is the whole idea in one experiment. The RHF energy of water, in three basis sets of growing size:

ao=

kind

basis functions (\(K\))

RHF energy / \(E_h\)

sto-3g

minimal

7

−74.963023

cc-pvdz

double-zeta

24

−76.026772

cc-pvtz

triple-zeta

58

−76.057127

The energy drops steadily toward the true Hartree–Fock limit as the basis grows — exactly as the variational principle promises (a bigger, more flexible trial form can only lower the energy). But notice \(K\) grows fast (7 → 24 → 58), and the cost of the two-electron integrals grows much faster than \(K\), so in practice you choose the smallest basis that is accurate enough. This basis-set incompleteness is one half of the “basis and method” caveat behind the slightly-too-large water dipole in the tutorial; the other half — the method — is next.

Tip

What qc-rs does with ao= ao="cc-pvdz" tells qc-rs which contracted GTOs to place on each atom. ao_rep="spherical" (the default) uses the \(2\ell+1\) spherical harmonics per shell (5 \(d\)-functions, not 6); ao_rep="cartesian" keeps the Cartesian set. The count \(K\) above is for the spherical default.

With the wavefunction reduced to a finite matrix problem, we can finally specify \(\mathbf F\) — the Hartree–Fock method.