#1

# CalculatingPotential EnergySurfaces withGaussian

9 January 2018

An important part of my doctoral research was using quantum chemical calculations to simulate organic molecules and gather theoretical information about them. One of these studies involved calculating the potential energy surface for rotation around a dihedral angle of a molecule. Googling this provided a good deal of helpful websites explaining how to code the input file for these types of jobs, but I couldn’t find a decent walk-though using GaussView. Thanks to a very helpful colleague, I was introduced to GaussView’s redundant coordinate editor, and that’s what I’m going to walk you through today.

Queue flashbacks of your endless lectures on alkane stereochemistry and conformational isomerism. For this system it should seem completely logical to expect the molecule’s most stable conformation to be where the large electronegative chlorine atoms are anti to each other, minimising the destabilizing electrostatic interactions between them (though have you heard of the Gauche effect?). So lets prove that computationally.

First things first, draw DCE with the chlorine atoms in anti conformation and run an optimisation to get the molecule optimised before we start the scan.

`# opt b3lyp/6-31g(d)`

### Setting up the scan

The tool for scanning dihedrals in GaussView is called the redundant coordinate editor. Open up your output file containing your geometry optimised structure of DCE. To set up redundant coordinates. Go to Edit > Redundant Coordinates…

• Click Add in the top left corner.
• Click in one of the four boxes with a ? in, then in your molecule window, click on each of the four atoms of the bond you want to scan. It’s important to click them in succession, don’t click them in a random order! So in this case it is 8, 1, 4, 7, but the reverse of that is analogous.
• Change the dropdown from Unidentified to Dihedral
• Change the second dropdown from Add to Scan Coordinate
• In the new area that comes up, change the values for Take # Step(s) of size # degrees to whatever you deem responsible for your molecule.

That’s it! You’re done, click OK. There is also a section in that window that lets you set the initial dihedral angle if you so wish.

### Checking the input file before running

Lets examine a typical dihedral scan input file:

```%mem=2GB
%nprocshared=16
%chk=dce_pe1.chk

# opt=modredundant b3lyp/6-31g(d)

Dihedral PE scan of dichloroethane

0 1
C                 -2.06409138   -0.62360123   -1.25740497
H                 -2.42234231    0.38464245   -1.25838371
H                 -2.41916602   -1.12912802   -2.13105517
C                 -0.52409393   -0.62118176   -1.25599770
H                 -0.16584419   -1.62942543   -1.25464200
H                 -0.16901940   -0.11532905   -0.38253613
Cl                 0.06258402    0.20794606   -2.69334109
Cl                -2.65076821   -1.45326554    0.17962928

D 8 1 4 7 S 360 1.000000```

The two key lines to notice here are

`# opt=modredundant b3lyp/6-31g(d)`

which tells Gaussian to expect redundant internal coordinate definitions (used for scan or to constrain elements) – this site has an approachable in-depth explanation if you’re interested – and…

`D 8 1 4 7 S 360 1.000000`

which is the key line and contains all the information we entered into the redundant coordinate editor:

• D : Dihedral – we want to look at a dihedral angle
• 8, 1, 4, 7 : The number (coordinates) of the atoms in question – this relates to the C–Cl–Cl–C dihedral
• S : Scan – we want to perform a scan of the dihedral
• 360 : How many steps to take during the scan
• 1.0 : How many degrees to rotate between each step

Doing a 360 step scan is very extreme, for most systems fewer steps would give you a decent sized dataset (remember that after each rotation Gaussian will perform an optimisation of all degrees of freedom apart from the dihedral which remains fixed, then does a single point energy calculation, so the bigger your molecule the longer the calculation time). Most dihedral potential energy surfaces should also be symmetrical around an given 180 degree scan, but I don’t want to suggest you just scan 180 degrees in the event yours isn’t (and there are reasons why it might not be).

### Results

Once your job has run, open up your output file in GaussView. You can also download mine if you want to play around without actually running the calculation. In GaussView, go to Results > Scan. You should see a neat little graph like this.

showing the results of the scan. You can play the animation of your dihedral rotating and watch the little red circle zoom along the curve. Gaussian considers energy in units of hartrees but to better understand the energetics its helpful to convert these to eV or kJ/mol. If you right click on the graph window you can save the curve data as a text file, which you can use for plotting your own graph.

```# Scan of Total Energy
# X-Axis: Scan Coordinate
# Y-Axis: Total Energy (Hartree)
# X Y
179.9786056287 -999.0189993100
180.9786056287 -999.0189949430
181.9786056287 -999.0189809930
182.9786056287 -999.0189576350
183.9786056287 -999.0189249060
184.9786056287 -999.0188828640
...```

The total energy of each single point energy calculation is given in hartrees. We’re less interested in the actual values of the total energy of each conformer, but do want to know the difference in energy between each conformer and the lowest energy conformer. Simply find the smallest total energy – the most stable conformation – and take it from all the other values to get ΔE in hartrees, and then convert this energy to eV using

`1 eV = 0.0367493 hartree`

to get the ΔE (eV).

If we call the most stable conformation 0°, we can plot the relative energies of the conformations of dichloroethane in eV, shown below with an image showing the Newman projections of the molecule at the minima and maxima of the scan. Considering that at 25 °C thermal energy kT = 0.257 eV, and the above potential energy surface of DCE at the B3LYP/6-31G(d) level of theory, we can suggest that, at 25 °C, 1,2-dichloroethane has free rotation between -150° → +150°, however the molecule does not possess sufficient energy to rotate the chlorine atoms past each other.

Next Post →

February 2018