Chemistry 8003: Computational Chemistry

Problem Set 3

Winter Quarter 1996

( Due 3 / 11 / 96 )

This problem set will be very difficult to finish with the time you have remaining on your Cray accounts. I recommend that you work entirely on one problem before beginning the other (both will be worth 50 points). Once your account is closed, you will not be able to log in to see your files, so copy out the information you need as it becomes available.

A practical reminder: To run a job, you're input file should be called myfile.dat. Simply issue the command qrun myfile.dat and the script will take care of the rest. If you want to look at an output deck while a job is running, you can either vi myfile.out or more myfile.out.

A nomenclature reminder: the notation x/y//w/z means level of theory x using basis set y at a geometry optimized at level of theory w with basis set z. E.g., MP4/6-311G**//HF/6-31G* means the geometry was optimized at the HF/6-31G* level but the energy (and/or other properties) are being calculated at the MP4/6-311G** level.

    Some quick notes/reminders with respect to Gaussian94:
  1. Always include the keyword: scf=direct.
  2. Always include z-matrix as an option to fopt.
  3. Always specify freq=noraman for frequency calculations.
  4. To find transition states, use at least the keywords fopt=(ef,ts,z-matrix)
  5. You can save a lot of time by using useful information from previous calculations stored in the checkpoint file. Plan your calculations to try to save time.
    1. Keywords guess=read and geom=checkpoint get the wave function and the geometry, respectively, from the last completed calculation. So, if you have just done an optimization, and want to follow-up with a frequency calculation, you will certainly want to use these keywords. Note: Frequencies must be calculated using the same level of theory with which the geometry was optimized!
    2. If you just optimized a geometry at one level of theory, and want to repeat that optimization at a different level, include readfc in the fopt=() keyword. E.g., fopt=(ef,ts,readfc,z-matrix). This causes the program to read the force constants from the previous calculation, which speeds up the geometry optimization.
  6. When you do a CCSD(T) calculation, you also get the CCSD and MP2 energies--don't waste time running these as separate jobs. When you do a CCSD you also get MP2.
  7. MP4 as used in the problem set implies "full" MP4. In the G94 output that is often distinguished as MP4SDTQ (for singles, doubles, triples, quadruples). The keyword mp4 implies this level. When you do an MP4, you also get MP2 and MP3.

The Problems: (note that there is an appendix to simplify the compilation of some answers--report absolute energies in au to 5 decimal places, relative energies in kcal/mol to 1 decimal place, bond lengths in Ångstroms to 3 decimal places, and bond angles in degrees to 1 decimal place--written comments should still be provided separately)

1. Last problem set, we used AM1 and PM3 to calculate the "heats of formation" of the following molecules: F2, Cl2, ClF, and ClF3. Results were awful. We will now continue at the ab initio level. Since the first three are small molecules, we will attempt to reach basis set convergence with respect to various properties. To do that, we will use the correlation-consistent polarized valence double zeta (cc-pVDZ is the abbreviation used in the Gaussian94 keyword line), triple zeta (cc-pVTZ), and quadruple zeta (cc-pVQZ) basis sets of Dunning. The QZ set, besides being valence quadruple zeta, includes 2 sets of f functions and 3 sets of d functions on all atoms, i.e., it's big. The QZ calculations are so expensive that I have done them for you, and they are in the appendix-I have also filled in the CCSD(T)/cc-pVTZ energies there, so you only need to do up to a CCSD calculation to fill in the rest of the cc-pVTZ rows. Fill in all empty positions in the appendix tables. Note that all of the correlated calculations are done at the HF/cc-pVQZ geometries! You must use the right bond lengths from my calculations to get the right energies-bond lengths are: F2, 1.32749788 Å; Cl2, 1.97826267 Å; ClF, 1.58547516 Å.

To fill in the thermal correction block (because we want to compare to an enthalpy of reaction), perform a frequency calculation for each molecule at the HF/cc-pVDZ level. In the output file, take the value listed as Thermal correction to Enthalpy=.

The "best estimate" block at the bottom involves taking the best calculation that can be afforded with the cc-pVQZ basis, and correcting it for higher orders of correlation using results obtained with a smaller basis set, and finally adding the enthalpic correction. That is, Best estimate = CCSD/cc-pVQZ + [CCSD(T)/cc-pVTZ - CCSD/cc-pVTZ] + DDH. (The term in brackets thus analyzes the effects of triple excitations on the correlation energy using a basis set where the calculation can be afforded.)

Questions: How well converged are the RHF bond lengths by the cc-pVQZ level? Based on your results, estimate how much remaining error there might be relative to a "saturated" basis set. What about the convergence of the absolute energies? Are the correlated energies better converged or less well converged with respect to basis set size than the RHF energies? How does the best estimate for the heat of reaction to form ClF compare with the experimental value of -26.6 kcal/mol?

Other than using a still larger basis set, how could this calculation be further improved?

Recall that neither AM1 nor PM3 could reproduce the correct geometry for ClF3. Optimize the geometry of ClF3 at the RHF/cc-pVDZ level and report it with all necessary degrees of freedom. Using all of your data so far, make a best guess at the heat of formation of this molecule, and explain your reasoning. Don't do any additional calculations to answer this-you can't afford it!


2. Recall the ammonia N-oxide example used in class (if you took notes, you've already got a starting Z-matrix, although you'll recall it wasn't a great initial guess). This molecule is a tautomer of hydroxylamine, NH2OH. The two tautomers are interconverted via a 1,2-proton shift. Calculate the relative stabilities of the two tautomers and the barrier to tautomerization at the levels indicated in the appendix table. Also, provide the bond distances rNO, rNH, and rOH for the transition state structure at all levels of theory (where H is the proton in flight). It will be more efficient if you first find all of the structures at the RHF level, and then use these to start calculations at other levels. Verify your RHF transition state structure with a frequency calculation.

Questions: What is the magnitude of the imaginary frequency in your RHF TS structure? Looking at the associated displacement vector, what atom is moving the most in the imaginary vibration (which, of course, corresponds to the reaction coordinate)? Does the correlation energy appear converged by the MP4 level with respect to tautomer energies and barrier height? How do the RHF and the DFT results compare to the correlated levels? If you were going to assign a best estimate with an error bar for the barrier, what would it be (i.e., pick the numbers in the sentence: "The barrier is ##.# ± ##.# kcal/mol"). Explain your reasoning.