5.

Basic Methods of Computational Chemistry

Hartree–Fock method, DFT methods, theoretical foundations of electron correlation methods.
One-electron basis sets.

The Hartree–Fock Method

Expectation Value of the Energy (with a Determinantal Wave Function)

The expectation value of the energy of a (two-electron) system is:

The determinantal form of the wave function can be inserted into the above equation, and the
energy terms can be expressed as the expectation value of the individual terms of the
Hamiltonian. The result is:

1
(The exchange integral is zero for electrons with different spin.)

Hii is , the one-electron energy of the electron in orbital φi (electron–nuclei attraction

plus kinetic energy). Hnuc is a constant term that comes from the repulsion between the nuclei.
The energy is essentially a functional because everything can be expressed with the orbitals.
Formulation of the Hartree–Fock equations

The resulting equations are the so-called Hartree–Fock equations:

where the Fock operator is:

2
where the two operators are the Coulomb and exchange operators, respectively.
This is not yet a real eigenvalue equation, since the Lagrange multipliers εij couple the
equations. There is, however, a particular transformation that makes the matrix of the Lagrange
multipliers diagonal. This transformation transforms the orbitals but leaves the wave function
unchanged. The orbitals corresponding to this particular choice are called canonical orbitals
and can be obtained from the eigenvalue equation form of the canonical Hartree-Fock
equations.

(Notice the following:)

The Hartree–Fock–Roothaan method

Because the Hartree–Fock equations are still very complicated equations, we express functions
as linear combinations of basis functions and operators as their matrix representations in the
same basis. In practice, atomic orbitals are used for the basis set, which leads to the well-known
LCAO-MO method:

This step implies a further approximation. The optimization now means finding the best
coefficients in a finite basis set. If we insert the above form of the orbitals into the HF equations,
we get the Hartree–Fock–Roothaan equations:

where Fba = , Sba = , and Ci being the vector of the coefficients Cai, defining
the ith MO. A compact form of the above equations can be obtained if we collect all Ci vectors
into a matrix and all εi into a diagonal matrix:

3
Interpretation of the Results of the Hartree–Fock Method
Density Matrix
The probability of finding an electron (any electron) at a given point in space (r) is the electron
density:

For a determinantal wave function, this it the sum of the electron densities of the individual
orbitals.

Elements Pab are the elements of the density matrix 𝑃.

Mulliken Population Analysis

Integration of the electron density over the enetire space gives the number of electrons:

4
Orbital Energy, Total Energy
The quantity εi in the Hartree-Fock equations is the exact definition of the orbital energy. This
is the expectation value of the Fock operator for orbital φi.
εi =

(Jii = Kii, they cancel.) The sum of orbital energies is not equal to the total energy, because the
electron-electron interaction would be included twice. Thus, the total energy is:

E=

Ionization Energy, Koopmans’ Theorem

If we assume that the orbital energies are unchanged upon removal of an electron
(approximation) from an orbital of the molecule, then the (first) ionization energy of the
molecule is the negative of the energy of the orbital where the electron is removed from. This
is Koopmans’ theorem.
Similarly, electron affinity could be calculated as the energy of the empty orbital where the
electron is placed. This is, however, a much poorer approximation because the empty orbital is
not well described.

Density Functional Theory (DFT)

The wave function of an N-electron system is a function of 4N variables, which determines all
properties of the system, e.g., the energy. On the other hand, the density is a three-dimensional
functional, which gives the probability of finding an electron (any of them) at a given point.
The density can be calculated from the wave function. Does the density also define all
properties of the system? After all, it defines the distribution of electrons in the molecule.
The Hamiltonian can be partitioned as:

Where 𝑇̂ is the kinetic energy operator and Vee is the electron-electron repulsion. V(r) is the so-
called “external potential”.

5
The Hohenberg–Kohn Theorems

First Hohenberg–Kohn theorem: A V(r) external potential is defined completely by the

(ground state) electron density.
In other words: no two different potentials of the nuclei can give the same electron density (with
N fixed).
Therefore, the external potential can be obtained from the electron density:

6
Second Hohenber–Kohn theorem: E[ρ] assumes its minimum value for the correct ρ(r) if the
admissible functions are restricted by the condition ∫ 𝜌(𝑟) 𝑑𝑟 = 𝑁.
This means that the “best” density can be obtained by minimization, provided the form of the
functional is known.

Unfortunately, the form of the functional is not known, and HK2 does not say anything about
how the ground state electron density can be obtained.
The Kohn-Sham (KS) DFT
We introduce orbitals to approximate the kinetic energy:

TS is not the same as T[ρ]:

with correction Ec the “correlation contribution”, which is not known.

The two unknown contributions to the energy functional, the correlation contribution Ec and
the exchange contribution Ex are often combined into a single term:
Exc[ρ] = Ec[ρ] + Ex[ρ]
This is the exchange-correlation functional, which collects all unknown contributions.
Unfortunately, theory tells us virtually nothing about the form of this functional.
Thus, the final form of the Kohn-Sham energy functional is:

We can apply the variational principle (according to HK2) to this functional to obtain the “best”
density, the one which gives the lowest energy. Since the density in this method is given by the
Kohn-Sham orbitals, in practice we are looking for the “best” orbitals. Thus, again, a functional
of the energy with these orbitals should be set up and minimized. Again, the optimal orbitals
will be solutions of one-particle equations. This is called the Kohn-Sham equation:

7
The equations are solved in an iterative procedure. We get orbitals and can calculate the density
from them. In principle, all properties are determined by the ground state density, but the
functional is not known, not even for the energy.

DFT Functionals
There are several levels of approximation (Jacob’s ladder):
1. Homogeneous electron gas: Vxc can be explicitly calculated from the density
2. Gradient-corrected functionals: Vxc depends on the density and its gradients; its form
includes empirical parameters.
3. Hybrid methods: the exchange part of the potential is that of the Hartree-Fock method;
new empirical parameters.
4. Correction for long distance interactions: empirical form, obtained from wave function
methods.
The following table lists the names of some commonly used functionals according to the leel
of approximation. The name of the complete exchange-correlation functional is the combination
of the names of the exchange and correlation parts (e.g., B3-LYP).

LYP: Lee, Young and Parr correlation functional

B3: Becke’s hybrid functional
(The B3LYP method has really good accuracy.)

8
Wave Function Methods: Inclusion of “Electron Correlation”
General Considerations About Electron Correlation
The Hartree–Fock and DFT methods are passed on the independent electron model. This is a
basic limitation because the instantaneous (pillanatnyi) correlation of the motion of electrons
is not taken into account, but only an average interaction.
The correlation energy is often defined as the difference between an ideal ‘exact’ energy and
the HF energy (if the basis set is infinitely large).
A more general wave function can be constructed by using more than one determinant by adding
so-called ‘excited determinants’. These are simply created by promoting electrons from
occupied orbitals to selected, originally empty virtual orbitals. Mathematically, more basis
functions improve the energy according to the variation principle. Physically, with the excited
determinants there may be electron distributions reachable that are just missing in the HF wave
function. All electron correlation methods are based on this idea of a multi-determinant wave
function.
(Electron correlation methods can be improved step by step, systematically, by including higher
and higher excitations. DFT, on the other hand, cannot be improved systematically.)
The Configuration Interaction (CI) Method
Excited determinants are formed by promotion of electrons from occupied orbitals to virtual
orbitals. These determinants form an orthonormal basis in the many-electron space. By
formulating the wave function as a linear combination of determinants, we obtain a
configuration interaction (CI) wave function:

If all possible excited determinants are included in the wave function, then the wave function
is restricted only by the finite one-electron basis. This method is called full CI, but it is
extremely expensive. More practical methods are defined by truncating the expansion after a
certain level of excitation: CIS, CISD (Configuration Interaction with Singles and Doubles),
CISDT, etc. After deciding on the level of approximation, the wave function parameters, i.e.,
the coefficients c of the expansion, should be obtained.

We can write the CI wave function in a general way: .

The best coefficients are found by minimizing the energy (var. princ.): .

A constraint is that the wave function should be normalized, which can be introduced using the
Lagrange multiplier Λ. The final eigenvalue equation that must be solved is: ,
where the matrix elements are: , and c is a vector of the cp coefficients.

9
The cost of the calculation grows exponentially with the number of one-electron basis functions
and consequently with the size of the molecule. Therefore, in practical calculations, limitation
on the excitation level must be imposed.
A meaningful and tractable approximation is CISD, which uses single and double excitations
only:

Doubles (doubly excited determinants) give the largest contribution to the energy, while singles
are needed for one-electron properties (plus their number is small).
CISD is not a satisfactory approximation if there are low-lying virtual orbitals, or at long
internuclear distances, or for calculating excited states.
Truncated CI methods are not size-consistent, which means that the wave function of two non-
interacting systems (infinite separation) is not the same for the two cases when the calculation
is performed on the supersystem, and when it is performed individually for the two systems
where the wave function of the supersystem is calculated as the product of the individual wave
functions. In the latter case, the wave function includes also higher-level excitations for the
supersystem. Consequently, the physical properties of the supersystem, including the energy,
will not be the same in the two cases, even though they should be. This error is also present in
the case of interacting systems, and this is referred to as the lack of size-extensitivity.
The Coupled-Cluster (CC) Method

10
These methods are very accurate for ground states of atoms and molecules. However, the cost
increases exponentially with the excitation level. The CCSD(T) method is very accurate at an
acceptable cost and is often referred to as “the golden standard of quantum chemistry”.

11
Perturbation Theory
Perturbation theory is a basic technique of (approximately) solving the Schrödinger equation
incases the Hamiltonian can be split up into a major component with known solution and a
smaller component as ‘perturbation’. In the present case, we treat “electron correlation” as a
perturbation on a known wave function taken from a lower-level approximation.

We partition the Hamiltonian: ̂0 the zeroth order part of the Hamiltonian

, with 𝐻
and 𝑉̂ the perturbation. It is assumed that the solution for 𝐻
̂0 is known: .
̂0 is the sum of
In quantum chemistry, we use the Møller–Plesset (MP) partitioning, where 𝐻
the one-electron Fock operators: . The zeroth order solution is the

“uncorrelated” system, and the perturbation 𝑉̂ is the difference between the full electron-
electron interaction and the averaged approximation. The unperturbed wave function is the
Hartree-Fock wave function: . The energy is the sum of the orbital energies:
, but this is not the Hartree-Fock energy. That is calculated as the expectation

̂ for the Hartree–Fock wave function, which is, in this context, the first order
value of the full 𝐻
energy:

Consequently, there is no 1st order PT (‘MP1’). Also, Hartree–Fock is a fairly good

approximation because the full Hamiltonian is used to calculate the energy. It is only the wave
function that ‘does not know’ about electron correlation.
Formulae for higher-order energies can also be obtained. In these, the perturbed wave function,
up to Ψⅈ−1 , also appears in the formula.
2nd order PT (MP2 or MBPT(2)) is very cheap and can be used to check if electron correlation
is important in the investigated system. MP3 is more expensive but not better. MP4 is the
highest level in routine PT calculations, but it is less accurate than CCSD(T). MP2 is popular,
the others are rarely used.
(The perturbation series may not converge.)

12
Basis Sets
All quantum chemistry methods make use of the LCAO-MO approximation: the molecular
orbitals are expanded in atomic orbitals. The set of functions chosen to represent atomic orbitals
is called the basis set. Both the type and the size may strongly influence the results of
calculations.
Types of Basis Functions
Slater-Type Orbitals (STOs) have equations similar to the eigenfunctions of the hydrogen
atom, with ζ (zeta) as a parameter of the function. Unfortunately, calculating the integrals of
these orbitals is difficult.
Gaussian-Type Orbitals (GTOs) have the following equation:

Here, i, j, k are integers, and these are used to form the necessary s, p, d functions (e.g., i = 1
and j = k = 0 for px). The only parameter is ζ. The product of two Gaussian functions is also a
Gaussian function, therefore their integral can be calculated easily. However, GTOs have a
different shape than STOs: their fall-off is too slow, and they are not pointy enough at the
position of the nucleus.
To overcome these flaws and improve accuracy, fixed combinations of so-called primitive
Gaussian functions gb can be used:

Such linear combinations are called contractions or Contracted Gaussian-Type Orbitals

(CGTOs). Parameters dab are obtained from preliminary calculations.
Number and Tpye of Basis Functions in Basis Sets
The quality of a basis set depends strongly on the number of functions per atom. This is
described by the so-called zeta number which gives the number of functions used to describe
an atomic shell.
Single-zeta (minimal) basis sets contain one set of functions for each occupied sub-shell. For
example, the O atom has one 1s, one 2s, and one set of 2p functions. Double-zeta (DZ) basis
sets have two sets of functions for each occupied subshell, while triple-zeta (TZ) basis sets have
three.
The addition of functions with higher angular momenta than the occupied ones, e.g., d functions
for a carbon atom, improves accuracy. These are called polarization functions. For example,
the DZ basis extended by polarization functions (DZP) includes:

13
In cases when the wave function extends to a larger space (like in the case of anions or when
intermolecular interactions need to be described), the inclusion of diffuse functions is necessary
for satisfactory accuracy. These functions have smaller ζ exponents and decay slower,
extending to distant spaces from the nucleus.
Gaussian Basis Sets
Split-valence basis sets use only one CGTO to describe atomic orbitals in the inner shell of
atoms, because the accurate description of the electronic structure of the valence shell is what
is important to chemists.
For example, the 6–31G basis set uses 1 CGTO made up of 6 GTOs for each AO in the inner
shell, while it uses 2 CGTOs, one made up of 3 GTOs and the other made up of 1 GTO, for
atoms in the valence shell. In this example, the C atom would have 1 CGTO (6 GTOs) for its
1s function, 2 CGTOs (3+1 GTOs) for its 2s function, and 3×2 CGTOs (3×(3+1) GTOs) for its
2p functions.
Polarization functions and diffuse functions can be added to these sets. For example, 6–31G*
includes one set of d polarization functions (5×1 CGTOs, 5×1 GTOs) for second-row atoms,
while 6–31G** also has polarization functions for H atom. 6–31+ has a diffuse function on
second-row atoms, 6–31++ has a diffuse function also on H atom.
Correlation consistent basis sets are designed specifically for calculations that include
electron correlation. These basis sets are abbreviated as cc-pVXZ. cc = correlation consistent;
p = polarized (they all include polarization functions); V = valence (they are also of split valence
type); XZ = x-ple zeta (DZ = double-zeta, TZ = triple-zeta), X is the cardinal number (written
as a number). Systematic improvement of the results can be achieved with increasing cardinal
number.

14