Quantumchemistry, NWI-MOL406, Computer exercises, week 5

A severe test for DFT methods is the calculation of atomization energy (the energy needed to break up a molecule into separate neutral atoms). We'll calculate the atomization energy of H2, and next that of CH4.

Beware: use the online manual!

The hydrogen atom

As a pre-requisite, we first calculate the H atom. As basis set take: 6-311++G** (overkill, but we don't want to change it later). Do an energy calculation for 3 DFT functionals and Hartree-Fock (HF):
  1. a plain LDA: Slater exchange and Perdew-Zunger (1981) correlation,
  2. a good gradient corrected functional: PBE (Perdew, Burke and Ernzerhof),
  3. a standard hybrid functional (DFT and HF mixed): B3LYP.
  4. standard HF (you can do that from the "dft" mode, using the keyword "XC HFexch")
Note that this atom is an open shell system (look in the manual how to set multiplicity). Make a list of the energies and a list with the eigenvalues... amazed?

The atomization energy of H2

Now we're going to calculate the atomization energy of H2. Build the molecule, optimize the structure, and calculate the frequencies. Do it for the 3 different DFT functionals above. Make a list of total energies, bond lengths and frequencies. The experimental bond length is 0.74 Angstrom, The experimental atomization energy at 0 Kelvin is 103.5 kcal/mole. Why did you have to calculate the frequencies? Again reflect why you might be amazed...

The atomization energy of CH4

Now you have done this for H2, do it for a slightly more challenging molecule: CH4. This example ought to show the superiority of functionals 2 and 3 over 1. Note that CH4 is a tetrahedral molecule. Think about the spin-multiplicity of C. Inspect the eigenvalues of the alpha and beta orbitals. Note that in LSDA (local spin density functional) theory we have two densities, one for the alpha and one for the beta spin.
The experimental atomization energy is 392.5 kcal/mole.

A chemical reaction

Now we'll do a simple chemical reaction. We'll use CH4, H2 and C2H4. So first repeat the exercise of CH4 for C2H4 (this should be routine now). Again you'll need an experimental number (531.9 kcal/mole). Next calculate the (T=0K) enthalpy change of: 2 CH4 → C2H4 + 2 H2, both experimental and for the 3 functionals. What do you observe? Do things get better for the LDA? Why?

The water dimer

Use B3LYP and optimize. Start from the model below [from Halkier et al., Theor. Chem. Acc. (1997) 97:150-157] that is already good.
O -1.41084  0.063712  0
H -1.74764 -0.415688 -0.756952
H -1.74764 -0.415688  0.756952
O  1.48172 -0.066612  0
H  0.53235  0.055606  0
H  1.83806  0.82179   0
The two molecules are linked via a hydrogen bond. Its strength (at 0 K) is approximately 3.5 kcal/mole. Calculate the strength in two ways: just by subtracting the energy (forget about zero-point contributions) of the fragments from that of the complex, and by doing the same with BSSE (Basis Set Superposition Error, see note) correction (see dft part of manual). What is happening here? How do energies change with this so-called counterpoise correction method? (you don't need to further optimize the isolated molecules, that effect should be minute, but you can check).

Adding and removing an electron

If time is left, you can do this. We're going to look a bit into ionization potentials and electron affinities. First revist you results on H. In EXACT DFT the highest occupied KS orbital should give you minus the IP. How large is the error here? Again you're confronted with a large self-interaction error.
Let's try to do it the "Delta-SCF" way: You do separate calculations on Li+, Li and Li-, and obtain the IP and EA as total energy differences. You can use the `CHARGE' keyword to make ions. Just use plain LDA. Do a series of calculations for various basis sets: 3-21G 6-311G* 6-311G** 6-311+G* Calculate both IP and EA as a function of basis set, and also keep track of the KS eigenvalue of the highest occupied KS orbital. You'll note something specials for the negative ion. Compare to the experimental numbers (IP = 5.39 eV, and EA = 0.62 eV). You'll observe good agreement, both IP and EA, but there is something fishy here... what? how do you explain the variation of the EA with basis set?


Back to [ Quantum Chemistry ]
Valid HTML 4.01! Last updated: 14-Nov-2014