Chapter 5. Quantum chemistry in Molecular Modeling


5.3 Solving the electronic Schrödinger equation : Hartree-Fock Self-Consistent Field theory.

5.3b Molecular Orbitals : the LCAO/MO method

Sofar we have not discussed the functional form of the wavefunction. In almost all practical calculations the orbitals are Molecular Orbitals, written as linear combinations of basis orbitals which resemble the orbitals of atoms ("Atomic Orbitals"). We will come back to those basis orbitals later. Let us first consider the implication of the LCAO-MO method, as it is often called.

MO is a linear combination of N basis functions :


Note that the expansion of the wavefunction in terms of basis functions leads to a limitation of the accuracy of the ab initio Hartree-Fock approach only because there is a limited number of basis functions available.
The greater the number of basis functions (provided they are well chosen) the better the wavefunction, the lower the energy. The limit of an infinite basis set is known as the Hartree-Fock limit. This energy is still greater than the exact energy that follows from the Hamiltonian because of the independent particle approximation, see section 5.4.

The HF equations can be rewritten as follows :

in which the overlap is

while the elements of the Fock matrix are :


In this equation P is the density matrix :

The significance of the density matrix is that it describes the electron density of the molecule. Note that the criterion for judging convergence of the SFC procedure can refer to the density as well as to the energy, because both have to be stationary at self-consistence.

The two-electron integrals over atomic basis functions give rise to a major practical problem in the application of the ab initio HF method because for a reasonable quality calculation a large number of basis functions is required, resulting in a very large number of two-electron integrals, approximately N4/8 for N basis functions. Not only is the calculation of the integrals time-consuming, also their storage on disk is practically impossible for large systems.
Recently, direct SCF methods have become available which reduce the problem dramatically. In this approach the two-electron integrals are not stored but recalculated as required. This makes sense because the cpu of modern computers is very fast, while i/o operations take quite a lot of time. Secondly, only those integrals that are expected to have a significant value are actually calculated.
Furthermore, selection and convergence criteria can be adjusted dynamically, leading to tighter criteria as the calculation approaches completion.
With these tricks built into modern programs the direct algorithms are actually faster than the conventional ones for systems of more than about 100 basis functions (depending on the particular computer). On small workstations direct SCF methods are the only practical option, even for small systems.

Basis functions

Above it was stated that the basis functions in the LCAO-MO method are atomic orbitals. Indeed, the basis orbitals used in practical calculations mostly are atom-centered functions that resemble orbitals as they can be found for isolated atoms.
The radial part of such orbitals is an exponentially decaying function. Basis orbitals of this type are called Slater-type orbitals (STO).
For practical calculations they have the disadvantage that evaluation of integrals involving such functions is time-consuming. Therefore these orbitals are approximated by a linear combination of gaussian basis functions (GTO):

As can be seen in the figure below, the functional form of such a GTO (Gaussian-type orbital) is different especially in the vicinity of the nucleus. Therefore at least three GTO's with different exponents alpha are necessary to give a reasonable basis orbital.

A sample STO and two different GTO's

The exponents alpha and contraction coefficients k can be determined in different ways, e.g. by fitting to an STO or by optimizing the energy in ab initio calculations on atoms and small molecules. Once these values are determined they define a standard basis set, and are not changed in further calculations on other molecules.

Standard basis sets

Several standard basis sets are nowadays commonly used [4,6,7,8].
The most popular are the basis sets designed by J.A. Pople and his coworkers. These basis sets are not necessarily the best, but they have the advantage that they are in widespread use so that there is much experience with their performance in practical calculations.
The so-called minimal basis sets have one basis orbital per two inner shell electrons and one basis orbital for each valence atomic orbital.
Thus, for first-row elements there are basis functions resembling 1s, 2s, 2px, 2py, 2pz atomic orbitals. The STO's are replaced by n GTO's (STO-nG). The common minimal basis set is the STO-3G basis.

In order to increase the flexibility of the SCF wavefunction one can increase the number of basis functions per atom.
In Double Zeta basis sets there are two functions for each AO of the minimal basis,
one which is closer to the nucleus,
the other allowing for electron density to move away from the nucleus.
For first row elements this gives 1s, 1s', 2s, 2s', 2px, 2px', 2py, 2py', 2pz, 2pz' basis functions.
When this doubling of the minimal basis is done only for the valence orbitals one gets the Split Valence Basis Sets. This results for first row atoms in 1s, 2s, 2s', 2px, 2px', 2py, 2py', 2pz, 2pz basis orbitals.
Some well-known examples of split-valence basis functions are 3-21G and 6-31G.

In the 3-21G basis the 1s AO of a first-row element is represented by a fixed combination of 3 GTO's, the 2s (2px etc.) are approximated by a fixed combination of 2 GTO's and the extra valence orbitals 2s'(2px' etc.) are just one GTO. The 3-21G basis set offers a reasonable compromise between computational burden and quality of results, and is at the moment the most often used starting basis set for ab initio work.
Its "bigger brothers" such as 6-31G (1s : 6 GTO's; 2s (2px etc.): 3 GTO's; 2s'(2px' etc.) : one GTO) give a much lower energy but at the expense of having more integrals over GTO's to be calculated, while for chemical phenomena (energy differences) the advantage is only small.

Thus, after 3-21G the usual next step is to go to basis sets which have polarization functions, e.g. d-orbitals for first row elements which allow for a lower symmetry of the electron distribution in a molecule compared to an atom.
Note that in this case the exponents cannot be optimized by calculations on atoms because the polarization functions are not occupied. The 6-31G* (6-31G + d-functions for first row atoms; also written 6-31G(d)) or 6-31G** (6-31G + d for first row, + p for H; 6-31G(d,p)) basis sets are at present considered to give good quality geometries for most molecules, at least when the Hartree-Fock method is applicable (see section 5.6).

In practice it has been shown that good geometries can often be obtained with simple basis sets such as 3-21G, but that relative energies are better described using more extended basis sets.
A commonly used notation has emerged for this procedure, e.g. : 6-31G**//3-21G means a calculation with the 6-31G** basis, at a geometry optimized with the 3-21G basis set.

Example

To get some idea of the consequences of the choice of the basis set let us look at a practical example in terms of the number of two-electron integrals, computation time and results. In a study of the conformational analysis of 3-methyl-1,3,5-hexatriene [9] the geometries were optimized for different isomers and conformers.

tZt conformer
tEt cEt equilibrium

The molecular formula is C7H10.
The 3-21G basis set has on each carbon atom 3 "shells" (1s, 2sp, 2sp'), 9 basis functions, on each hydrogen 2 shells (1s, 1s'), 2 basisfunctions. In total there are 83 basis functions.
In a conventional SCF calculation (using the HONDO program) there were 3 x 10^6 two-electron integrals kept (N^4/8 = 6 x 10^6).
The total energy was E = -269.33677 au. (tEt form).
The cpu-time needed for geometry optimization was ca. 90 minutes. With the 6-31G* basis each C has 4 shells (1s, 2sp, 2sp', d), 15 basis functions, each H 2 shells (1s, 1s'), 2 basis functions, resulting in a total of 125 basis functions. In this case ca. 15 x 10^6 two-electron integrals were kept (N^4/8 = 30 x 10^6).
The time for geometry optimization was ca. 360 minutes. With the bigger basis set the total energy is much lower : E = -270.84127 au. (tEt), but the chemically significant energy difference between the most stable tEt form and the cEt rotamer was 2.9 kcal/mol with both basis sets, and the predicted torsion of the s-cis part was also essentially the same (33 vs. 34 degr.).
Interestingly, this energy difference which is compatible with the experimental fact that the cEt form is not observed at room temperature, was predicted to be much smaller with the usually rather reliable MMP2 method.
Thus, in this example ab initio gives a correct result while molecular mechanics failed. The geometry optimization at the SCF level with the 6-31G* basis set already took a substantial amount of cpu time and disk space (15 x 10^6 x 8 bytes = 120 MB).

In section 5.6 we will discuss the quality of HF results in more detail.

References:

  1. Levine Quantum Chemistry, 4th ed., 1991
  2. Reviews in Computational Chemistry, Volumes 1 to 4, K.B. Lipkowitz and D.B. Boyd, eds., VCH publishers.
  3. G. Nàray-Szabò, P.R. Surjàn and J.G. Angyàn Applied Quantum Chemistry, Reidel, 1987
  4. W.J. Hehre, L. Radom, P. von R. Schleyer and J.A. Pople Ab initio molecular orbital theory, Wiley, 1986
  5. Szabo and Ostlund Modern Quantum Chemistry, McGraw-Hill, 1989.
  6. Spartan User's Guide, version 3.0, Wavefunction, Inc., 1993.
  7. Foresman, J.B.; Frisch, A. Exploring Chemistry with Electronic Structure Methods: A Guide to Using Gaussian, Gaussian Inc., 1993.
  8. D. Feller and E.R Davidson, in Reviews in Computational Chemistry, K.B. Lipkowitz and D.B. Boyd, eds., VCH, 1990, pp. 1 - 43.
  9. Brouwer, A.M.; Bezemer, L.; Jacobs, H.J.C., Recl. Trav. Chim. Pays-Bas 1992, 111, 138-143

Next paragraph, 5.4 Limitations of the HF method; Electron correlation
Previous paragraph 5.2 The Schrödinger equation
Chapter 5 MM Syllabus 1995 MODIFIED November 8, 1995
Fred Brouwer, Lab. of Organic Chemistry, University of Amsterdam.