Gaussian geometry optimisation

Geometry optimisation using Gaussian

This exercise follows the drylab Building Z matrices using Molden and uses the Z matrices for [PCl₄]⁺, PCl₅ and [PCl₆]^- saved in your unix filespace then as .com files. You will now use the modelling program Gaussian, running in unix, to optimise the geometry of these species, starting from these Z matrices. Besides writing a .log file of all the results it calculates, which can be read by human beings as well as by Molden, Gaussian also saves a checkpoint file .chk in machine-readable code, for its own future use. This contains much more detail of molecular orbitals etc.. In the final drylab you may do the last exercise in this series, in which Gaussian will read the checkpoint files it saves now, and carry out vibrational analysis calculations to predict stretching frequencies.

Gaussian is not designed for you to interact with it directly: it is always used in batch mode

You give it a ready-made .com command file to tell it what to do
It works silently in the background in a unix computer
It writes its results to a .log file
When it has finished, you read the .log file using more or a screen editor, or Molden can read it for you

You find out what Gaussian has done only after it has finished doing it

Because molecular orbital calculations can be very long, involving numerical calculation of millions of integrals, taking days of machine time, this is the obvious practical way: if you try to do too big a geometry optimisation in a PC program such as Arguslab, where you are watching it happen in real time, you can get very tired of waiting
The present exercises, involving only moderately small molecules and only the pm3 semiempirical method, take only a few seconds of time each in powerful unix computers, so are suitable for a drylab environment in which you wait for Gaussian to finish

Molden saves only the Z matrix for a species in its .com file, as you have seen. To make a full command file for input to Gaussian you need to add command lines to the file

For this you will use a screen text editor, which may also be used to examine the .log file of results

Screen editors in unix

Whereas in Microsoft Windows there is a standard text editor Notepad, there is no comparably universal standard screen editor in unix

This partly is because many unix computers are not used via their own screens, so a screen editor has to work via a terminal emulation in a distant computer
The most standard unix screen editor is called vi, which works via terminal emulators providing a limited range of basic services, such as moving the mouse cursor or clearing all or part of the screen as vi requires it
I find that vi is not very user friendly
Rather than teach you to use vi, which you may also find elsewhere, I have decided to teach the use of Curlew instead. Curlew is a local Newcastle product which you are unlikely to find in other institutions
Curlew is no longer supported by the Computing Service, but is still provided in their unix machines
Curlew works similarly to vi, but is, in my opinion, easier to use and is ideally suited to scientific work. In particular, it has very good facilities for moving columns of data and for duplicating blocks of text, as needed for making up input files for programs like Gaussian
The Curlew keyboard is completely programmable so that combinations of key presses can be used e.g. to search for known phrases in the voluminous Gaussian .log file, e.g. the tests of geometric convergence
If you are already an experienced user of vi, there is no reason why you should not go on using it. Obviously, the instructions below for using Curlew will not then apply to you

Before proceeding with these instructions, you should now read the separate document Curlew in unix via an Exceed xterm window. Do the setting up specified there, and refer back to the table of Curlew commands there whenever you need them in the work here.

Making a command file for Gaussian, to optimise the geometry of PCl₅

After setting up the Curlew keyboard in Exceed, open two xterm windows in aidan
In one window

cd ~nbwt3/work/listen

more templatepm3.com

This file contains template command lines for a geometry optimisation Gaussian job, followed by instructions on how to do it. I shall give more explicit instructions now instead, but you will copy the template lines from this source

In the other window

cd g98

cl pcl5.com

This starts a Curlew session editing the partial command file you saved from Molden, to which you are now going to add the command lines
The black cursor should be over the first character of line 1, which is the line specifying that p is the first atom in the Z matrix
Press control-o to split the line at this point
You now have a blank line into which to paste the command lines for Gaussian

In the first window, containing the template, drag the mouse cursor from the beginning of the % sign of the first line (starting %chk= ) down to immediately after the 1 in the line containing 0 1

Do not drag so far that the line blacks across to the right of the window
Let go of the mouse
This takes a little skill with the mouse: if you do not succeed, left click anywhere in the window, and try again

Left click on the title bar of the second window, containing your Z matrix, to raise it without changing the selection

Curlew's cursor should still be at the beginning of the new blank line
Centre click (e.g. press the mouse wheel inwards) to paste the lines from the template into your .com file
There should not be a blank line between the 0 1 line and the first line of your Z matrix

If there is, you dragged too far. The black cursor should be on the blank line. Hold down the control key and tap u once to delete the line

You now need to alter the command lines you have pasted in from the template
Using the arrow keys, move the cursor to the beginning of YOURUNIXID in the first line of the file
Type in your n.. unix id
Delete the remaining characters of YOURUNIXID using control-e or the Delete key
Move to the beginning of FILEFIRSTNAME and press control-y to delete to the end of the line
Type pcl5 , making sure you do not leave a space after the / after your id

This line now will tell Gaussian to put a checkpoint file called pcl5.chk in a directory called /usr/local/tmp1/n.. where n.. is your user id
You do not put the .chk file extension on the end of the filename, because Gaussian will do that for you
Do not be concerned if the directory /usr/local/tmp1/n.. does not yet exist, because the script which will run Gaussian for you will automatically find out whether it exists, and if it does not, it will create it
/usr/local/tmp1 is a directory into which you, or Gaussian acting on your behalf, can put directories and files which are too big to go in your file space

The files do not get charged to your account, but they are rubbed out from time to time by the operating staff
Very probably, the file(s) you are creating now will survive for at least a week, until you use them in the next drylab. If not, it is very quick to create them again, because you will still have the command file you are making now, to do it with
When this course is completely over, it would be suitably public spirited if you signed on to unix and gave the command

/usr/local/tmp1/n..

rm *

The second line of the file (%Mem= ) tells Gaussian how many 8-byte words of memory to use for its calculations: you see that it is set to 64 MByte
The third line is the command line to tell Gaussian to run a pm3 geometry optimisation
Now, separated by a blank line before and after it, comes your title for the job. You are allowed up to five consecutive lines, but one should be enough

Alter YOUR TITLE FOR THE MOLECULE to a sensible title for PCl₅, so that you would know what the resulting output was about if you found it later in your file space. Remember you are not allowed to use subscripts

The 0 1 line tells Gaussian that the molecule has no charge and is in a singlet state, i.e. all the electrons are paired
Use Esc then h to go to the end of the file

Make sure there are at least two blank lines after the variables list for your Z matrix
Press Enter to add another one to make sure

Save the file with Esc then x

Running Gaussian

For running a modelling job of any size, Gaussian is normally used in a 'batch' environment, in which your job queues for the large resources required. However, the current jobs are so small, because you are using a semi-empirical method on a small molecule, rather than doing ab initio calculations, that you can run Gaussian simply in the background from your present unix signon. The essential is that you should remember to put the & to cause a background process, on the end of the command line, before you press Enter

From your work/g98 directory, where you have stored the completed command file pcl5.com, give the command

g98 pcl5 &

Notice that you do not type the full name of your file pcl5.com, because Gaussian adds the .com ending for you

You should get a receipt of [1] (meaning background job 1) followed by a process number
If you have already run one of these jobs, you will also get a message 'n.. exists', which is the script telling you that it does not need to create a directory for you in /usr/local/tmp1. You can ignore this
Give the BWT command

recent

You will see that a file pcl5.log has been created by Gaussian
After probably a couple of minutes of elapsed time, when you press Enter you should see a receipt

This shows that the Gaussian run had ended
You should not have seen any error messages starting ***
If you did, seek help as to what is wrong with your .com file, repair it with Curlew, and try again

Looking at the output from Gaussian

Only when Gaussian has finished, look at the output file, using Curlew, by giving the command

cl pcl5.log

Use Esc then h to go to the bottom of the file
The last line of the file should be 'Normal termination of Gaussian 98.'
Before that, you will see one of Gaussian's large store of quotations, and how long the computer took to do the calculation
Use Esc then g to go to the top of the file, then use function key F2 to scroll down the file, to examine the following

Initial Parameters

These are the internal coordinates of the molecule, set up from your Z matrix, which are to be refined

Distance matrix

You can see all the distances, bonded and non-bonded, in the molecule
This is followed by the point group, which should be D3H for PCl₅

Standard basis

This will be VSTO-3G
This is a very simple basis set, with just one basis function of three Gaussian primitives for each valence shell orbital
A little further down you see the number of basis functions, the number of primitives, and the number of electron pairs. You should confirm that these are all correct for the valence shell of PCl₅
How many filled MOs and how many empty MOs should this produce?

Simple Huckel Guess

Gaussian has used the extended Huckel semiempirical method to get the starting guess at the electron distribution. At the end of the Leave Link line just below that, you see that it took only a fraction of a second of processor time to do that

Tap Esc then n. This is one of BWT's customisations of the Curlew setup to look for the next occurrence of the changes of parameters table from Gaussian

You will see that only R4 and R5, the bond lengths to the apical chlorine atoms, are changing much in this first iteration. They are getting shorter

Tap Esc then c. This BWT customisation looks for the next occurrence of the table of four tests for geometric convergence

All four tests have got a long way to go before passing the preset thresholds

Press function key F7 to repeat the search

Each time you do this to get to the results of the next iteration, the forces and displacements have got less
Finally, all four tests have a YES flag under Converged?

Use F3 to bring the top of the table of Optimized Parameters (American spelling!) to the top of the window

You see that Gaussian has maintained the symmetry. If it was D_3h then the three equatorial bonds should have stayed equal in length to each other, and the two apical bonds show a slightly greater length, as taught in Inorganic Chemistry courses for the known structure of PF₅

Search further down the file for 'Total atomic charges'

These are calculated by the Mulliken method
The charge on P is very positive, and is hard to believe. Probably there should be more pi bonding than allowed in this model, which would reduce the charge separation by Cl to P pi donation
The apical chlorine atoms carry more negative charge than the equatorial chlorine atoms, which is one of the principles used by VSEPR in predicting the shapes of 5-coordinate species

Exit from Curlew using Esc then z. If you are asked whether you want to discard changes, reply y
In your other xterm window, if you still have more running, exit from it by pressing q

Give the commands

cd work/g98

Now look at the same output file, using Molden
Give the command

molden

In Molden's control window, click Read, then select pcl5.log from the list of files. Close the File Select dialog box
Switch to Solid, Ball & Stick and Shade and put on Sticky Pointer rotation
Reduce the model to a reasonable size with Out, and rotate the model so that you can see the trigonal bipyramidal structure clearly
Click Movie to see an animation of the geometry optimisation

You just see the bond lengths change, which all that the symmetry allowed to happen
You can see the changes more slowly if you click First and then repeatedly click Next

Click Geom. conv. to see graphs of the changes in energy, step sizes and forces during the optimisation

This is the same information, presented graphically, which you saw as text when you examined the four convergence tests in the output using Curlew
Click Geom. conv. again to remove the graphs, and the skull and crossbones to exit from Molden

[PCl₆]^- and [PCl₄]⁺

Repeat the setting up and optimisation procedures above for the other two models
Instead of copying and pasting lines from the template file, you can use more to view your own file pcl5.com, and copy the command lines from there. You will then have less alterations to make. Do make sure you change the name of the checkpoint file in the first, %chk= , line, otherwise you will overwrite your checkpoint file for PCl₅, which you will need for doing its vibrational analysis
When you edit the .com files, remember to specify that [PCl₆]^- is an anion and that [PCl₄]⁺ is a cation, otherwise you will have an odd number of electrons and the jobs will fail because you are saying that the species are in the singlet spin state, which is impossible for an odd number of electrons

Do this by changing the 0 1 line to -1 1 or 1 1 respectively

For the tetrahedral species [PCl₄]⁺, do not be surprised that Gaussian says that the symmetry is only C_3v

This because Molden specified the tetrahedral angle, for three of the angles, with a rounded decimal part, or at least it is when translated into the binary number system
As a result, Gaussian calculates the remaining three angles as being slightly different, and forming another set of three. (There are six edges to a tetrahedron, and each one subtends a bond angle)
This is elongation along one of the threefold axes of the tetrahedron, so the symmetry is C_3v
This could be got around by specifying the ligand positions as alternating corners of a cube, using dummy atoms at the centres of the faces in a Z matrix, when only whole-number angles are involved, or, more simply, by specifying Cartesian coordinates directly in the .com file instead of a Z matrix
With a Z matrix made in the obvious way, as you are using, there is nothing you can do about it
Look at the Optimized Parameters for the final geometry: you will see that all four bond lengths round to the same number, but if you look down a further 20 lines, you will come to the Distance matrix, where distances are given to higher precision. You see that one of the P-Cl bonds is longer than the others by one in the last decimal place, in keeping with C_3v symmetry, which Gaussian has taken steps to retain
In fact, this difference in length is far less than the accuracy of the method you are using, so the molecule is really tetrahedral to within the error of the calculation