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
cd /usr/local/tmp1/n..
followed by to delete all the files you have there
• 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
• 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
[1]   Done   g98 pcl5
• 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
cd work/g98
to get back to your own file space in that window as well
• Now look at the same output file, using Molden
• Give the command
• 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