Lab 4
Geometry Optimizations

4.1 Overview

In this lab the process of obtaining molecular structures by calculation will be considered. We will begin with an alternative coordinate system that is useful for chemical applications. We will then consider some of the basic theories of molecular structure such how a stable structure is defined in terms of the molecule's potential energy surface. Following this, we will deal with how a geometry optimization converges and how we can follow the progress of a calculation. Finally, we will make use of geometry optimizations in several mini-projects.

4.2 Coordinate Systems

In the second lab we encountered a few ways of describing molecular structures in the context of reading a structure from a file. One of the structure formats described there, XYZ coordinates, is used extensively in structure optimizations. An alternative method that we haven't considered yet, which is perhaps more relevant for chemists because it deals in bond lengths, angles and dihedral angles is the system of internal coordinates or the z-matrix. For further information on internal coordinates visit http://chemistry.umeche.maine.edu/Modeling/GGZmat.html, http://cbc.arizona.edu/courseweb/074/CHEM5181/Zmatrix/ZmatrixDef.htm and http://en.wikipedia.org/wiki/Internal_coordinates. It should be noted that these web pages outline one way that a z-matrix can be written but not the only way.

Exercise 4.1: After consulting some or all of the web pages listed above write down a z-matrix to describe the geometry of hydrogen peroxide H₂O₂. Assume that the O-O bond length is 1.5 Å, that the O-H bond length is 1.0 Å, that the O-O-H angle is 109.5 degrees and that the H-O-O-H dihedral angle is 180 degrees.

The web page http://www.shodor.org/chemviz/zmatrices/babel.html is a form that will convert a z-matrix into xyz-coordinates.

Chemical reactions are generally thought of in terms of making and breaking bonds so internal coordinates are often a good choice when studying reactions.

Figure 4.1: The Model:Coordinates option showing both Cartesian and internal coordinates

In order to look at the Cartesian coordinates or internal coordinates that correspond to the geometry that you are investigating with ADFINPUT choose the Model:Coordinates option. The window on the right will now contain information about your molecule in terms of Cartesian coordinates (figure 4.1). The menu above the coordinates definition (default: Cartesian) allows you to choose between Cartesian, internal and reduced coordinate systems. We will not be considering reduced coordinates in this course.

Exercise 4.2: Build and clean up the geometry of the methanal molecule CH₂O. What z-matrix does ADFINPUT obtain for this molecule?

4.3 Geometry Optimizations: The Basics

Implicit in the idea of a molecule is the concept of molecular structure. That is, a molecule is not only a collection of atoms but is a collection of atoms in a particular set of locations in space.

Often, the structure of a molecule is of interest. Even if structure is less important, most important molecular properties are structurally dependent. Determining molecular structures is therefore a very important aspect of chemistry. In the area of computational chemistry, structure determination is achieved through a geometry optimization.

The energy of a collection of atoms as a function of their 3N-6 internal degrees of freedom is known as a potential energy surface (PES). Molecular structures correspond to minima on this surface. The function of a geometry optimization is therefore to find a minimum in the potential energy surface.

4.3.1 Measuring Geometrical Parameters

Once a molecular structure is obtained, it is often necessary to derive quantitative values of important parameters such as certain bond lengths and angles from the structure. Sometimes the required parameters can be read from the text output or output summary files. If the molecule is displayed in the build window of ADFINPUT then these parameters can also be obtained simply in the following ways.

Figure 4.2: Measuring the distance between two atoms.

Bond Length: Select two atoms. The distance between the two atoms is printed at the bottom right of the build window (figure 4.2).
Bond Angles: Select three atoms. The angle between the three atoms is printed at the bottom right of the build window. Note that the order in which the atoms are selected is important.
Dihedral Angles: Select four atoms. The dihedral angle made up by the four atoms is printed at the bottom right of the build window. Note that the order in which the atoms are selected is important.

4.3.2 Modifying Geometrical Parameters

It is often useful to be able to modify certain properties of the molecular structure. For example, if a dihedral angle can be modified then butane can be changed from the anti to the gauche conformer. The parameters described above can be modified using ADFINPUT. When a particular parameter is displayed by selecting the appropriate atoms, that parameter can be modified by editing the value displayed.

Alternatively, if the Model: Coordinates menu option has been chosen then the Cartesian coordinates or Z-matrix can be edited directly. Note that if an atom or atoms are selected in the build window then the geometrical parameters of those atoms are highlighted (figure 4.2). The selected atoms can be moved up or down in the overall list which may be important when creating a Z-matrix.

Exercise 4.3: Build and clean up the water molecule. What are the O-H bond lengths and the H-O-H bond angle obtained?

4.3.3 A Simple Geometry Optimization

Figure 4.3: Choose the Geometry Optimization option from the Main Options.

Unless an experimental geometry is available, the first step in a computational study is usually a geometry optimization to obtain an accurate structure of the molecule(s) of interest. The simple building and cleanup process described in lab 2 is unlikely to give a geometry accurate enough for use in further studies especially if heavier elements are involved.

All of the calculations carried out so far have been the "single point" type meaning that they considered only a single point on the PES. To instead choose to perform a geometry optimization, click on the menu next to preset on the right of the window (below title and above Basis Set) and choose the option Geometry Optimization (figure 4.3). If the options to the right don't seem correct, check that Main Options is at the top (figure 4.3). If this is not the case click on the button to the far top right. This should give a long list of options. The one that you want (Main Options) is at the top.

If you now save your work and run the calculation, ADF will look for a minimum on the potential energy surface starting from the geometry that you have specified.

Exercise 4.4: Build an N₂ molecule and run a single point calculation on it. Choose the Geometry Optimization option and run a second calculation. Use the default DZ basis set and the LDA functional in both cases. Compare the final energy obtained in the geometry optimization calculation with that from the single point calculationn. Would you expect the energy of the optimized geometry to always be lower than that of the geometry you start from? Would you expect it always to be higher? Would you expect it to sometimes be higher and sometime be lower? Why?

4.3.4 Obtaining the Results of a Geometry Optimization

Figure 4.4: An optimized geometry. The picture on the left shows where the results are found in a logfile and on the right in the detailed output.

There are a number of ways to determine what the results of the geometry optimization calculation were. You can:

Clear previous work using the File: New command and then using File: Import Coordinates import the results of your calculation from the TAPE21 file resulting from your calculation. Make sure that you select the file type ADF results files (*.t21) . Once the molecule has appeared in the build window, bond lengths, bond angles and torsion angles can be measured as described above.
The final geometry can be found in the summary output (logfile) of the calculation under the text Geometry Converged near the end of the file (figure 4.4). See section 1.5.1 in Lab 1 for instructions on how to open a logfile.
The final geometry can also be obtained from the detailed text output file. This file is opened using the Output command from the SCM menu. Once the file is open, the final geometry is reached though the Section:Final Geometry command (figure 4.4).

Further Information: The TAPE21 file saves a wide range of results from your calculation and not just the final optimized geometry. This file is called X.t21 where X is the filename that you gave for your calculation.
Exercise 4.5: What is your final optimized bond length for N₂? How does it compare with the starting geometry? Find the experimental geometry of N₂ at http://webbook.nist.gov/chemistry/. Which is closer to the experimental value, the starting guess or the final optimized value?

4.4 Convergence

In a geometry optimization the calculation doesn't just go immediately to the optimized geometry. The program takes information from the starting geometry and then evaluates a new geometry that should be closer to a minimum in the PES than the previous geometry. The information that is used includes the energy, the first derivative of the energy with respect to changes in the positions of the atoms and sometimes the second derivative of the energy with respect to atomic position. At the new geometry, the energy and derivatives are obtained again and another geometry that should be even closer to a minimum is obtained. This process is repeated until the maximum number of steps is reached or the calculation is satisfied that it is close enough to a minimum. If the maximum number of steps is reached the calculation is "not converged." Usually, the next step is to run another geometry optimization this time starting from where the last one left off. If the program is satisfied that a minimum has been found then the geometry is said to be "converged." A minimum in a PES has zero first derivatives in the energy with respect to movement of the atoms. Thus, a calculation is "close enough to a minimum" when the calculated gradients at a particular geometry are closer to zero than a predefined value. In order to be more sure that a minimum has been reached other criteria are often used such as the change in energy or bond lengths and angles between steps.

Figure 4.5: Geometry iteration data.

Like most programs, ADF writes information about the convergence of a geometry optimization after every step of a calculation.

Exercise 4.6: Consider your geometry optimization calculation of N₂. Look at the detailed output file using the SCM:Output program. To look at the geometry steps of the calculation select the Geometry Cycles option from the Iterations menu of the OUTPUT program. Doing so will put you at the results of the first geometry step. The information presented includes the gradients of the energy with respect to nuclear movements, the energy of the current geometry and the previous geometry, the new geometry predicted from the information known and the convergence information (see figure 4.5). The convergence information includes "Item" (what we are considering), "Value" (the value of that item), "Criterion" (what Value needs to be less than for this item to be converged) "Conv" (Whether Value is less than Criterion or not) and "Ratio" (the ratio of Value from the current iteration divided by Value from the previous iteration. This gives an idea of whether you are getting closer to convergence or not). When all of the Conv are YES your calculation is converged.
In this exercise you should make up a table listing two of the items of convergence data (two out of "change in energy", gradient max and so on) for each geometry iteration. You can proceed to the next geometry cycle by again choosing the Iterations: Geometry Cycles menu option. Plot the two sets of the geometry convergence data against geometry cycle number. How do each of the criteria behave as the calculation proceeds?

4.4.1 Visualizing Geometry Optimizations

Figure 4.6: The MOVIE program.

An alternative way to follow the progress of is to visualize the structure of the molecule at each geometry iteration. ADF provides a simple method for doing this with the MOVIE program. Following a geometry optimization, open MOVIE with the SCM: Movie menu option. A new window will open displaying a movie of how the molecular structure changes as the optimization proceeds (see figure 4.6). The control buttons of the MOVIE program can be used to move through the frames of the movie. Most of the menu options of the MOVIE program are reasonably standard. The most interesting new option is Files: Save Movie Frames.... Selecting this option will save a picture for each geometry iteration. These pictures can be used to make an .mpeg or .gif movie for later use. We won't do this here because it will take a long time.

Exercise 4.7: Choose one of your geometry optimization calculations. Examine the progress of the optimization with the MOVIE program. Save two of the frames of the optimization with the Files: Save Picture... command of the MOVIE program. One of the frames should be near the beginning of the calculation and one near the end.

4.5 Projects

Do project 4.5.1, either project 4.5.2 or 4.5.3 and project 4.5.4 or 4.5.5 for a total of three projects.

4.5.1 Calculations with Different Basis Sets and Functionals

Calculate the geometry of H₂O with the SZ, DZP, TZ2P and QZ4P basis sets. Use the LDA functional. Prepare a table comparing the calculated geometries of the molecules using different basis sets. Your table should include the O-H bond lengths and the H-O-H bond angle. Include in your table how long each calculation took. Is there any particular trend in the optimized geometries going from SZ to DZP to TZ2P to QZ4P basis sets? Which basis set gives the best agreement between theory and experiment? The experimental O-H bond length for water is 0.957 Å and the experimental bond angle is 104.5 degrees.
Calculate the geometry of water with three different functionals of your choice. Choose from LDA or any of the GGAs. It is not possible to do geometry optimizations with the other functionals with ADF at the present time. Use a TZ2P basis set. Prepare a table comparing the calculated geometries of the molecules using different functionals. Your table should include the O-H bond lengths and the H-O-H bond angle. Include in your table how long each calculation took. Is there any particular trend in the optimized geometries going one functional to another? Which functional gives the best agreement between theory and experiment?

Isomerization

The concept of isomerization is an important one in chemistry. Now that we can optimize geometries, we are in a position to make predictions about isomers.

If anything unusual happens in the course of the calculations, describe the events and give an explanation if you have one.

Be sure to check that your geometries are converged.

4.5.2 Conformational Isomers: Gauche and Staggered n-butane

Build n-butane. Before running any calculations, predict whether the gauche or staggered conformation of n-butane will be lower in energy. Optimize the geometry of the molecule in the gauche and staggered conformers. In order to achieve the two conformers, you will have to use different starting geometries for your optimizations. Use a DZP basis set. Be sure to check your optimized geometries to ensure that they are what you expect. From the results, evaluate the energy difference between the gauche and staggered forms.

4.5.3 Conformational Isomers: Chair and Boat forms of Cyclohexane

Calculate the difference in energy between the chair and boat forms of cyclohexane. Use a DZP basis set.

4.5.4 Structural Isomers: Hard and Soft Acids and Thiocyanate

A ligand that can bind to a metal through two different atoms is called ambidentate. An example of such a ligand is the thiocyanate ion (SCN^-). It is not always clear ahead of time which way around an ambidentate ligand will be bound in a given complex. Calculations obviously can be used to make predictions about how a ligand will bind or to support whatever experimental evidence that can be obtained once a compound is synthesized.

Some idea about the likely orientation of a SCN^- ligand can usually be obtained using the concepts of hard and soft acids (http://www.meta-synthesis.com/webbook/43_hsab/HSAB.html) since the N atom in the ligand is "harder" than the S atom.

In these calculations we will try and predict whether SCN^- will be bound in a series of complexes through its S atom or its N atom and whether the results of these calculations agree with what is known about the hardness and softness of the central metals of interest.

We will be considering three complexes: [Sc(CO)₅SCN]²⁺, Cu(CO)₅SCN and [Zn(CO)₅SCN]⁺. Calculate the geometries of two forms of each of these complexes, one with the thiocyanate ligand bound through the S atom and one with it bound through the N atom. Use a DZP basis set. Evaluate the difference in energy between the two isomers of each complex and decide which you would expect to be made in a synthesis. Do the complexes follow the trend that would be expected from hard and soft acid arguments? Explain.

Note: All of these are closed-shell molecules so a restricted calculation is required. Some of them are cations so you may need to set the charge.

4.5.5 Clusters: Six atoms of Be

Sometimes very little information about the structure of a molecule exists beyond its formula and what might be guessed from chemical intuition. This is often the case with small clusters. With the rise of nanotechnology small atomic clusters, a type of nanoparticle, are increasingly being synthesized and investigated. It is not always obvious what the structure of a cluster with a given formula might be. Calculations can help by giving a clear prediction of a cluster structure or at least narrowing down the list of candidates.

Consider a cluster of six Be atoms. Propose as many structures with the formula Be₆ as you can. Build these structures and then optimize their geometry. How many stable structures do you find? Which is the most stable?

File translated from T_EX by T_TH, version 3.85.
On 26 Sep 2009, 16:37.

Lab 4 Geometry Optimizations