4.2 Theoretical Background

4.2.1 SCF and LCAO Approximations

The fundamental equation of non-relativistic quantum chemistry is the time-independent Schrödinger equation,

H^(𝐑,𝐫)Ψ(𝐑,𝐫)=E(𝐑)Ψ(𝐑,𝐫). (4.1)

In quantum chemistry, this equation is solved as a function of the electronic variables (𝐫), for fixed values of the nuclear coordinates (𝐑). The Hamiltonian operator in Eq. (4.1), in atomic units, is

H^=-12i=1Ni2-12A=1M1MAA2-i=1NA=1MZAriA+i=1Nj>iN1rij+A=1MB>AMZAZBRAB (4.2)

where

2=2x2+2y2+2z2. (4.3)

In Eq. (4.2), Z is the nuclear charge, MA is the ratio of the mass of nucleus A to the mass of an electron, RAB=|𝐑A-𝐑B| is the distance between nuclei A and B, rij=|𝐫i-𝐫j| is the distance between the ith and jth electrons, riA=|𝐫i-𝐑A| is the distance between the ith electron and the Ath nucleus, M is the number of nuclei and N is the number of electrons. The total energy E is an eigenvalue of H^, with a corresponding eigenfunction (wave function) Ψ.

Separating the motions of the electrons from that of the nuclei, an idea originally due to Born and Oppenheimer,109 yields the electronic Hamiltonian operator:

H^elec=-12i=1Ni2-i=1NA=1MZAriA+i=1Nj>iN1rij (4.4)

The solution of the corresponding electronic Schrödinger equation,

H^elecΨelec=EelecΨelec, (4.5)

affords the total electronic energy, Eelec, and electronic wave function, Ψelec, which describes the distribution of the electrons for fixed nuclear positions. The total energy is obtained by simply adding the nuclear–nuclear repulsion energy [the fifth term in Eq. (4.2)] to the total electronic energy:

Etot=Eelec+Enuc. (4.6)

Solving the eigenvalue problem in Eq. (4.5) yields a set of eigenfunctions (Ψ0, Ψ1, Ψ2) with corresponding eigenvalues E0E1E2.

Our interest lies in determining the lowest eigenvalue and associated eigenfunction which correspond to the ground state energy and wave function of the molecule. However, solving Eq. (4.5) for other than the most trivial systems is extremely difficult and the best we can do in practice is to find approximate solutions.

The first approximation used to solve Eq. (4.5) is the independent-electron (mean-field) approximation, in which the wave function is approximated as an antisymmetrized product of one-electron functions, namely, the MOs. Each MO is determined by considering the electron as moving within an average field of all the other electrons. This affords the well-known Slater determinant wave function849, 850

Ψ=1n!|χ1(1)χ2(1)χn(1)χ1(2)χ2(2)χn(2)χ1(n)χ2(n)χn(n)| (4.7)

where χi, a spin orbital, is the product of a molecular orbital ψi and a spin function (α or β).

One obtains the optimum set of MOs by variationally minimizing the energy in what is called a “self-consistent field” or SCF approximation to the many-electron problem. The archetypal SCF method is the Hartree-Fock (HF) approximation, but these SCF methods also include KS-DFT (Chapter 5). All SCF methods lead to equations of the form

f^(i)χ(𝐱i)=εχ(𝐱i), (4.8)

where the Fock operator f^(i) for the ith electron is

f^(i)=-12i2+υeff(i). (4.9)

Here 𝐱i are spin and spatial coordinates of the ith electron, the functions χ are spin orbitals and υeff is the effective potential “seen” by the ith electron, which depends on the spin orbitals of the other electrons. The nature of the effective potential υeff depends on the SCF methodology, i.e., on the choice of density-functional approximation.

The second approximation usually introduced when solving Eq. (4.5) is the introduction of an AO basis {ϕμ} linear combinations of which will then determine the MOs. There are many standardized, atom-centered Gaussian basis sets and details of these are discussed in Chapter 8.

After eliminating the spin components in Eq. (4.8) and introducing a finite basis,

ψi=μcμiϕμ, (4.10)

Eq. (4.8) reduces to the Roothaan-Hall matrix equation

𝐅𝐂=𝜺𝐒𝐂. (4.11)

Here, 𝐅 is the Fock matrix, 𝐂 is a square matrix of molecular orbital coefficients, 𝐒 is the AO overlap matrix with elements

Sμν=ϕμ(𝐫)ϕν(𝐫)𝑑𝐫 (4.12)

and 𝜺 is a diagonal matrix containing the orbital energies. Generalizing to an unrestricted formalism by introducing separate spatial orbitals for α and β spin in Eq. (4.7) yields the Pople-Nesbet equations750

𝐅α𝐂α = 𝜺α𝐒𝐂α (4.13)
𝐅β𝐂β = 𝜺β𝐒𝐂β (4.14)

In SCF methods, an initial guess is for the MOs is first determined, and from this, an average field seen by each electron can be calculated. A new set of MOs can be obtained by solving the Roothaan-Hall or Pople-Nesbet eigenvalue equations, resulting in the restricted or unrestricted finite-basis SCF approximation. This procedure is repeated until the new MOs differ negligibly from those of the previous iteration. The Hartree-Fock approximation for the effective potential in Eq. (4.9) inherently neglects the instantaneous electron-electron correlations that are averaged out by the SCF procedure, and while the chemistry resulting from HF calculations often offers valuable qualitative insight, quantitative energetics are often poor. In principle, the DFT methodologies are able to capture all the correlation energy, i.e., the difference in energy between the HF energy and the true energy. In practice, the best-available density functionals perform well but not perfectly, and conventional post-HF approaches to calculating the correlation energy (see Chapter 6) are often required.

That said, because SCF methods often yield acceptably accurate chemical predictions at low- to moderate computational cost, self-consistent field methods are the cornerstone of most quantum-chemical programs and calculations. The formal costs of many SCF algorithms is 𝒪(N4), that is, they grow with the fourth power of system size, N. This is slower than the growth of the cheapest conventional correlated methods, which scale as 𝒪(N5) or worse, algorithmic advances available in Q-Chem can reduce the SCF cost to 𝒪(N) in favorable cases, an improvement that allows SCF methods to be applied to molecules previously considered beyond the scope of ab initio quantum chemistry.

Types of ground-state energy calculations currently available in Q-Chem are summarized in Table 4.1.

Calculation $rem Variable JOBTYPE
Single point energy (default) SINGLE_POINT or SP
Force (energy + gradient) FORCE
Equilibrium structure search OPTIMIZATION or OPT (Ch. 10)
Transition structure search TS (Ch. 10)
Intrinsic reaction pathway RPATH (Sec. 10.5)
Potential energy scan PES_SCAN (Sec. 10.4)
Vibrational frequency calculation FREQUENCY or FREQ
Polarizability and relaxed dipole POLARIZABILITY, DIPOLE (Sec. 11.14.1)
NMR chemical shift NMR (Sec. 11.13.1)
Indirect nuclear spin-spin coupling ISSC (Sec. 11.13.1)
Ab initio molecular dynamics AIMD (Sec. 10.7)
Ab initio path integrals PIMD, PIMC (Sec. 10.8)
BSSE (counterpoise) correction BSSE (Sec. 13.4.3)
Energy decomposition analysis EDA (Sec. 13.5)
Symmetry-adapted perturbation theory XSAPT (Secs. 13.12 and 13.13)
Table 4.1: The type of calculation to be run by Q-Chem is controlled by the $rem variable JOBTYPE.