MO Calculation for a Diatomic Molecule
Introduction
The properties of any molecular system can in principle be found by looking at the solutions to the corresponding time independent Schrodinger equation (TISE) for the molecule. For a diatomic molecule AB having N electrons, the Hamiltonian for the molecule may be written as
H = Tn + Te + Vne + Vnn + Vee
where
Tn = - [ (2/2mA) A2 + (2/2mB) B2]nuclear kinetic energy
Te = - (2/2me) i=1N i2 electron kinetic energy
Vne = - [(ZAe2/40) i=1N (1/rAi) + (ZBe2/40) i=1N (1/rBi)]nuclear-electron potential energy
Vnn = (ZAZBe2/40) (1/rAB) nuclear-nuclear potential energy
Vee = (e2/40) i=1N j>iN (1/rij) electron-electron potential energy
The above can be extended in a straightforward manner for use in polyatomic molecules.
While writing the Hamiltonian operator is simple, finding solutions to the corresponding TISE is not. Some simplification can be obtained by making the Born-Oppenheimer approximation. This approximation is based on the observation that electrons are much less massive than nuclei. As a result, one can, to a very good first approximation, assume that the nuclei are stationary relative to the electrons, which are, due to their light mass, able to adjust to changes in nuclear geometry "instantaneously". As a consequence of this, the nuclear kinetic energy term (Tn) can be removed, while the nuclear-nuclear potential energy term becomes a constant whose value depends on the geometry chosen for the molecule (which, for a diatomic molecule, means the bond length).
The more simple TISE obtained when making the Born-Oppenheimer approximation is still too difficult to solve exactly (except for a single electron system, such as in the H2+ ion). However, there are a number of approximation methods that have been developed over the past decades that make it possible to find solutions where, within the approximations made, the error is small.
Molecular orbital theory
The first approximation usually made is the molecular orbital approximation. Each electron in the system is assumed to occupy a separate molecular orbital. The wavefunction for the molecule is then the product of the individual one electron orbital functions. Satisfaction of the Pauli principle (that the overall wavefunction be antisymmetric with respect to exchange of any two electrons) is obtained by writing the solution in terms of a matrix (the Slater determinant).
The starting point in finding the best set of orbitals for a molecule is the Hartree-Fock (HF) approach. Molecular orbitals are found by numerically solving the Hartree-Fock equation, a relationship that is derived from our initial TISE. Solutions are often obtained in terms of a set of functions called Slater type orbitals (STOs). Since the numerical calculations required when using STOs are very time consuming, the STOs are themselves usually approximated by a set of Gaussian functions. This substitution leads to great improvements in the speed with which the required numerical calculations can be carried out.
The quality of a particular molecular orbital calculation generally depends on two factors:
(1) Basis set - The basis set is the set of initial functions from which the final molecular orbitals are constructed. In principle the best basis set would be a complete set, which would make it possible to (within the Born-Oppenheimer and orbital approximations) find an exact solution to the TISE for a molecule. However, the time required to carry out the calculations involved scales approximately as the square of the number of basis functions being used, so as a practical matter a complete basis set cannot be used.
There are two common basis sets used in molecular orbital calculations:
Minimal basis set. Uses all of the core and valence atomic orbitals for each atom in the molecule.
Augmented basis set. Uses all of the orbitals in the minimal basis set, along with additional unoccupied atomic orbitals, called virtual orbitals. Using more virtual orbitals improves the quality of the calculation but increases the time required to carry out the calculation.
As an example, a minimal basis set for a molecular orbital calculation for the HF molecule would use a 1s atomic orbital for hydrogen, and the 1s, 2s, 2px, 2py, and 2pz atomic orbitals for fluorine (usually represented as STOs). This minimal basis set could be augmented by including the 2pz atomic orbital for the hydrogen atom (the p orbital that points towards the fluorine atom in the HF molecule), or by including virtual atomic orbitals from both the hydrogen and the fluorine atom.
The basis set used in a particular molecular orbital calculation is usually indicated by the label used to identify the calculation. For example, minimal basis sets are often written as STO-nG, where STO indicates that Slater type orbitals are being used in the calculation, and n is the number of Gaussian functions being used to approximate the Slater type orbitals. So, for example, the STO-3G basis set uses Slater type orbitals, each of which has been approximated as a combination of three Gaussian functions. More complicated basis sets will be indicated with different types of symbols (which are often difficult to decipher for those who do not carry out molecular orbital calculations on a regular basis).
(2) Corrections to the Hartree-Fock method. The Hartree-Fock method makes some approximations to find the solutions to the TISE for a molecule. There are a variety of methods that have been developed to go beyond the basic Hartree-Fock approach to improve the results that are obtained.
Most corrections concentrate on accounting for the correlation of motion of the electrons in a molecule. Correlation refers to the fact that the presence of an electron at a particular point in space affects the locations of the other electrons in the molecule, due both to electron-electron repulsion and spin effects related to satisfying the Pauli principle. A common method for correcting for the effects of electron correlation is called configuration interaction (CI). Configuraton interaction calculations make use of virtual orbitals. In full-CI all possible assignments of electrons to virtual orbitals are made. However, because of the time involved for such calculations, CI is often restricted to configurations where only a few electrons are allowed to occupy virtual orbitals. For example, CI-SD restricts the CI calculations to states with one (single) or two (double) electron excitations.
A second method for correcting for electron correlation is called the coupled cluster method (CC). This method is based on the theory of many body interactions in multiparticle systems. As is the case for the CI method, different levels of approximation are indicated based on the number of excitations that have been included. So, for example, CCSD includes all terms with one or two electrons in excited states.
An additional correction procedure, called the Moller-Plesset method, makes use of perturbation theory to improve the quality of the molecular orbital calculation. One starts with the solution to the Hartree-Fock equation, and then treats corrections to this solution as small perturbations, which are addressed using standard techniques. The label used for this type of correction indicates the order to which corrections have been made. For example, in MP2, corrections up to second order in the energy have been carried out, while in MP3 corrections up to third order have been made.
Additional methods for improving the quality of a molecular orbital calculation have also been developed. As is the case with basis sets, abbreviations to represent the types of corrections being included are usually indicated when a particular calculation is carried out.
Density functional theory
As an alternative to the molecular orbital approach, a second method, called density functional theory (DFT), has been developed for theoretical calculations on molecules. Instead of attempting to solve the TISE for a molecule, the density functional approach attempts to find the probability density for electrons in the molecule, (x,y,z). In 1964, Hohenberg and Kohn proved that all of the properties of a molecule can be calculated from this probability density. Unfortunately, there is no method for finding an exact expression for (x,y,z) for a molecule. However, as is the case for molecular orbital theory, there are a variety of approximation methods that have been developed to implement the density functional approach. The quality of density functional calculations for molecules is now in many respects comparable to that found using molecular orbital theory, although there are still difficulties in modeling molecular dissociation and transition states. Density functional theory is generally the best method to use in solid state calculations.
The experiment
In the present experiment you will carry out several different molecular orbital calculations on the diatomic molecule you were assigned in the literature search in CHM 3410L (or a molecule assigned to you if you did not do the literature search experiment). The three calculations you will do are as follows:
hf/6-31g* - A Hartree-Fock calculation using an augmented minimal basis set. The STOs are approximated using six Gaussian functions. The symbol 6-31 indicates use of a "split" basis set, which gives more flexibility to the calculation, while the symbol * indicates the use of additional polarized functions, one way of augmenting a minimal basis set.
mp2/6-31g* - The same as the first calculation, with the addition of Moller-Plesset correction terms out to second order.
b3yp/6-31g* - A density functional calculation using the same basis set of functions used in the Hartree-Fock calculation. The symbol b3yp stands for the Becke three parameter approximation to the Lee-Yang-Parr exchange correlation functional, which is essentially a method for accounting for electron correlation in the calculation.
The calculations will be done using Gaussian-98, an older version of a commercial pckage of programs used to carry out molecular orbital calculations (the most recent version of this package is Gaussian-09), which is available on the computer in CP 375. To get to the screen where you will input the data for the calculation click FILE---New, and then enter the following information:
% Section
Leave blank
Route Section
Enter one of the following
#hf/6-31g* opt freqfor the Hartree-Fock calculation
#mp2/6-31g* opt freqfor the Moller-Plesset calculation
#b3yp/6-31g* opt freqfor the density functional calculation
(Note - Opt freq stands for "optimal frequency". This instruction tells the program to find the minimum in the potential energy for the molecule, and to then find the vibrational constant e.)
Title section
You need to put something here or the program will not run. So, for example, you could put your name and the kind of calculation being done (Joens-HF-O2, for example).
Charge and multiplicity
For a calculation on a molecule the charge will be zero, while for an ion it will be the ion charge. The multiplicity is equal to 2S+1, where S is the total spin quantum number for the ground electronic state of the molecule. For example, for O2, S = 1, and so the multiplicity is 3. Enter this information as the charge, a space, and the multiplicity (0 3 would be the entry for a calculation for the O2 molecule).
Molecule Specification
This enters the geometric information used as a starting point in the calculation, and is very simple for a diatomic molecule. Enter the data in the same way as shown below for a calculation on the O2 molecule.
O
O 1 R
R=1.30
This means the following
O - identifies the first atom in the diatomic molecule, using the symbol for the atom (in this case, an oxygen atom)
O 1 R - identifies the second atom in the diatomic (in this case, O), and tells the program to look for a starting value for R, the bond distance.
blank line
R=1.30 - The starting value for R to be used by the program (in Å, angstroms; 1 nm = 10 Å).
The results of the calculation are not sensitive to the starting value used for R, as long as a reasonable value is used. For your trial value for R you can use the experimental value for re from the literature search, rounded to two significant figures and converted to Å. The program will vary the geometry until it find re, the equilibrium bond distance.
After entering the data, click RUN (upper right side of the screen). The program will ask for a file name. I would suggest you use your name and an abbreviation for the type of calculation (joenshf-O2, for example).
There are three pieces of information that you need to obtain from the calculation:
re (equilibrium bond length)
The value for re will appear in a table labeled "optimized parameters" found approximately in the middle of the output file, and again near the bottom of the output file. It will have the following appearance.
Optimization completed.
-- Stationary point found.
------
! Optimized Parameters !
! (Angstroms and Degrees) !
------
! Name Definition Value Derivative Info. !
------
! R1 R(1,2) 1.1677 -DE/DX = 0. !
------
In the above example the value of re is 1.1677 Å.
e (vibrational constant)
The value for e will be given towards the bottom of the output file, and will have the following appearance.
Harmonic frequencies (cm**-1), IR intensities (KM/Mole),
Raman scattering activities (A**4/AMU), Raman depolarization ratios,
reduced masses (AMU), force constants (mDyne/A) and normal coordinates:
1
SGG
Frequencies -- 1998.4307
Red. masses -- 15.9949
Frc consts -- 37.6366
IR Inten -- 0.0000
Raman Activ -- 43.0051
Depolar -- 0.3103
Atom AN X Y Z
1 8 0.00 0.00 0.71
2 8 0.00 0.00 -0.71
In the above example the value of e is 1998.4 cm-1.
Computer time required for the calculation
This is the amount of CPU time required to carry out the calculation. The larger the value the longer the calculation takes. This value will appear at the bottom of the output file, and will have the following appearance.
Job cpu time: 0 days 0 hours 0 minutes 6.0 seconds.
In the above example the calculation took 6.0 seconds of CPU time.
Note that there is a large amount of additional information given in the output file, including the distribution of electrons, the rotational constant, the dipole moment, dissociation constant, an estimate of the infrared and Raman intensity, and thermodynamic information calculated from the molecular parameters.
You should make a table of results that includes the values for re, e, and the CPU time required for each calculation. Your values for re and e should be compared to two experimental and two theoretical values from the literature, either from your literature search project in CHM 3410L or obtained from a new literature search.
Briefly discuss your results. Focus on which method seems to best agree with experiment. You should also look at whether your calculations do a better or worse job than the literature calculations in finding re and e, and comment on the relative amount of CPU time required for each calculation.
1