Gaussian is a complicated Fortran program which can perform a variety of computational tasks such as computing energies, optimizing geometries, computing vibrational frequencies, etc., and with a variety of computational methods, including simple Molecular Mechanics (force-field) approaches, Semi-empirical wavefunction-type methods, Hartree-Fock Self-Consistent Field methods, so-called "post-Hartree-Fock" or "correlated" methods, and Density Functional Theory approaches. Here you will learn how to:

- Generate Gaussian input files.
- Issue Commands to Gaussian.
- Specify the Molecular Geometry.
- Run a Gaussian Job.
- Read Gaussian output files.

Gaussian input is always in the form of a *unformatted text*
or ASCII file, usually called 'jobname'.com.
For historical reasons (Gaussian has been around for more than 25 years,
and new versions try to be partly compatible with older ones,
so the input format has not changed much), the order of the different items
in the file has to be followed very rigorously, because otherwise Gaussian does
not understand the input (typically, it tries to find piece of information 'A'
where you have written 'B', and that makes no sense). When this happens, the
job stops, and you get a rude & incomprehensible message in the Gaussian Output
(or logfile). Whenever this happens, consult the following section!!:

There are two types of lines of input in Gaussian:

- Single-line Input. This is on a single, self-contained line. E.g. the 'Charge and Spin' data below.
- Multiple-line Input.
This is on any number of lines. For Gaussian to know where this
group stops, it needs to be followed by
**ONE**blank line. Example: the geometry specification.

Here is a sample Gaussian Input File:

The sections indicated in this file are the most common input groups. In the order in which they should appear in Gaussian, these are:

**'Link 0' commands**. This is the way Gaussian designates general commands concerning the way in which the program should run, what it does with files, etc. For most jobs, you will not need*any*Link 0 command. The line shown here, which starts with the % sign, tells Gaussian to save to your current directory the temporary 'Checkfile', and to call it ethene.chk. If you do not do this, the checkfile is kept in a temporary, 'work', directory, and deleted when the job ends. Link 0 commands are*Single-line Input*.- The
**Command Section**. This is the most important part of the input: It tells Gaussian what to do ! Although it is only one line long here, it is*Multiple-line Input*and needs to be followed by one blank line. As with other Gaussian input, this section is case-insensitive. The first line of the Command Section needs to start with the '#' sign. Here, the program is being instructed to perform a geometry optimisation (opt) of the specified molecule, at the Hartree-Fock level of theory, and using the 3-21G standard basis set on all atoms, which it will read from its basis set library. - The
**Job Title**. This is free for you to choose. It can be as brief or as long as you like, and may span several lines. It is multiple-line input, so it needs to be finished by a blank line. As a more verbose alternative to the dull job title I have put here, the following would be fine, too: - The
**Molecular Charge and Spin**setting. This begins your description of the molecule you wish Gaussian to calculate. You need to specify the*overall molecular charge*and the*spin multiplicity*(as two integers). Here, ethene is a neutral species (charge=0), and it is a spin-Singlet, with S=0 and the multiplicity, 2*S + 1, equal to 1. - The
**Geometry Specification**. This is multiple-line input, so it needs to be followed by a blank line. Because it is often the last part of a Gaussian input file, you may forget to add a blank line, especially as with pico, you cannot*see*whether a blank line is there or not ! The format can be Cartesian, as here, or in Z-matrix format (see below).

This file is the first Gaussian Input File I look at.
Running it will cause Gaussian to perform
an ab initio Hartree-Fock Geometry
Optimisation of the ethene molecule!!!!

As mentioned above, the command line tells Gaussian what to do. It will generally need to specify the following things:

- The theoretical
*method*to use, including - The
*basis set*, in most cases. - The
*type of job*to be done, and any options chosen.

This is often put as the first item in the Command Section. It is written as:
*method/basis set*.
The methods you may meet include Hartree-Fock (written more simply as "HF"),
MP2, the density functional theory approach called "B3LYP",
or the semi-empirical "PM3" method.
You will always need to specifiy a basis set, *except for semiempirical
computations*, where no basis is needed. Basis sets you might need to use
include STO-3G, 3-21G, 6-31G, 6-31+G, 6-31+G**, aug-cc-pVDZ, ...
You can read about the methods and basis sets in the
*Ab Initio*
Background Section.

If you do not specifiy anything else than the method and basis set,
Gaussian will perform the default job-type, a
*single-point energy computation* at the geometry specified.
This can also be done by adding the keyword "SP". Here is a list of
jobtypes you will need, and some indications on options you will need:

- Single-Point Job. Keyword:
**SP**. - Geometry Optimisation. Keyword:
**Opt**. By default, this leads to optimisation of the nearest minimum (stable structure). This is generally trouble free. Optimising the geometry of a transition state is harder, requires you to specifiy several extra options. Frequency Calculation. Keyword:

>To avoid this problem, it is probably best to use**Freq**. This will calculate the matrix of second derivatives of the energy with respect to the position of the nuclei, and diagonalise it to generate vibrational frequencies.**WARNING:**This type of job only yields meaningful results if you have previously optimised the geometry using the*same theoretical method*and the*same basis set*. If you have*not*done this, you will*not*get a warning or error message from Gaussian, and you will apparently be able to read the frequencies in GaussView. The results will however be**meaningless**, as is mentioned (rather discreetly, though), in the output.*both keywords together*, as in: Opt Freq. Whereas you should normally only specify one job-type, these two can be given together. Gaussian first optimises the geometry, then performs a frequency job at the optimised geometry.

There are two main ways of specifying the molecular geometry: Cartesian
coordinates, and internal (or Z-matrix) coordinates. The former are more compact,
but they cannot readily be set up by hand. The latter, whilst they are longer, can
be used to specify any molecule's geometry in a text file. In the past, this is how the
geometry of most (all !) molecules was entered into programs such as Gaussian. Nowadays,
most large molecules are built using graphical programs such as
GaussView, and are then saved in Cartesian or
(more rarely) Z-matrix format. This is what you will do in the project about the Diels-Alder
reaction. However, for the triatomic system F + H_{2}, Z-matrix format is simple
to use and more straightforward, so it is described here.

Geometries in Cartesian coordinates are entered simply as follows:

label1 x1 y1 z1
label2 x2 y2 z2
label3 x3 y3 z3
...

(As shown in the ethene example, above). Note that the unit of distance is Angstrom.

A Z-matrix describes a geometry in terms of atom-atom distances, atom-atom-atom angles, and dihedral angles between two planes formed by four atoms. Graphically, this looks like this:

The image on the right shows the meaning of the dihedral angle d_{1234}: looking down
the central 2-3 bond, it is the angle between bonds 1-2 and 3-4 (this is also the angle between the
planes defined respectively by atoms 1, 2 and 3 and 2, 3 and 4).

The corresponding Z-matrix would be:

atom1
atom2 1 r12
atom3 2 r23 1 a123
atom4 3 r34 2 a234 1 d1234
(atom5 4 r45 3 a345 2 d2345
...)

The atoms referred to in the Z-matrix do not have to be in order as shown here. To make
this all a bit less abstract, here is a valid Z-matrix for H_{2}O:

O
H 1 0.96
H 1 0.96 2 105.

Note that the decimal point is included in the angle specification (which is in degrees).
This is **very important**. Without the decimal point, Gaussian interprets that number as
an integer, where it is expecting a real number. This will lead to rude (and incomprehensible)
error messages... (yet again !)

With a dihedral angle of 0 or 180 degrees, the four atoms are coplanar. Here is an example
Z-matrix for formaldehyde CH_{2}O:

C
O 1 1.35
H 1 1.1 2 122.
H 1 1.1 2 122. 3 180.

The numerical values for lengths and angles in the two previous examples can be replaced
by *variables*. The value that these variables take must then be defined in a *
separate, multiple-lined input section* situated just after the z-matrix. So the
formaldehyde example above could be rendered instead as:

C
O 1 rCO
H 1 rCH 2 aHCO
H 1 rCH 2 aHCO 3 d
rCO 1.35
rCH 1.1
aHCO 122.
d 180.

This is what Z-matrices generated by GaussView will look like.

Back to Top
To be able to rund a Gaussian job on Minerva, **you should be a member of
the Gaussian user group.**

Before running a Gaussian job you have to issue the following commands (once per session):

- If your default shell is bash then
`csh`

- To make Gaussian available for your shell:
`module load gaussian`

- followed by either

`source $g09root/g09/bsd/g09.login`

(when your login shell is C-shell)

or by

`source $g09root/g09/bsd/g09.profile`

(when your login shell is bash)

When you have created your input file using your favorite text editor, you should save it under the name 'jobname'.com. Then you can run the calculations, by typing:

`g09 < input-file > output-file &`

Gaussian writes its output to a file called the *logfile*, which has the same name
as the input file, but the '.log' extension. It writes to this file progressively as the job
advances. When you open the file with a text editor, you actually open a *copy* of that
file, so you can open the output file * while the job is running* to inspect the job
progress, without that affecting the way the output is written. You should be very careful
*not to save the file* if you have opened it in this way, however !!

The output is quite verbose, and the present discussion will focus on only a small part of the information you can find in the logfile.

The first few tens of lines of the output are taken up by system details and copyright information ... not very interesting !

Then you will see the following 'header':

***************************************
Gaussian 98: SGI-G98RevA.6 19-Oct-1998
6-Sep-2000
***************************************

Shortly after that, your Command Section and Job Title will be echoed.

Then, you will see lots of details about your geometry. In the following block, for example, you will see the Cartesian Coordinate representation of your molecule (even if you input it as a Z-matrix):

Standard orientation:
---------------------------------------------------------------------
Center Atomic Atomic Coordinates (Angstroms)
Number Number Type X Y Z
---------------------------------------------------------------------
1 6 0 0.000000 -0.677600 0.000000
2 1 0 -0.924521 -1.216266 0.000000
3 1 0 0.927705 -1.210763 0.000000
4 6 0 0.000000 0.677600 0.000000
5 1 0 0.924521 1.216266 0.000000
6 1 0 -0.927705 1.210763 0.000000
---------------------------------------------------------------------

Note that the atoms are represented by their * atomic number*, not
by their symbol.

Then, the calculation will start. You will see details of the initial 'guess' wavefunction which Gaussian generates, then, if all goes well, a message stating that your wavefunction has converged. This will look different for different types of job.

For a semiempirical (e.g. PM3) job, the relevant section will look as shown below. The total energy at this geometry is given by 0.174514138371 atomic units. 1 atomic unit (or Hartree) equals 2625.5 kJ/mol.

(...) It= 7 PL= 4.07D-04 DiagD=F ESCF= 47.482071 Diff=-9.42D-04 RMSDP= 3.93D-05. It= 8 PL= 1.76D-04 DiagD=F ESCF= 47.486848 Diff= 4.78D-04 RMSDP= 2.02D-05. It= 9 PL= 1.17D-04 DiagD=F ESCF= 47.486816 Diff=-3.18D-06 RMSDP= 3.07D-05. It= 10 PL= 3.47D-05 DiagD=F ESCF= 47.486769 Diff=-4.72D-06 RMSDP= 4.50D-06. It= 11 PL= 1.90D-05 DiagD=F ESCF= 47.486792 Diff= 2.31D-06 RMSDP= 2.94D-06. It= 12 PL= 1.37D-05 DiagD=F ESCF= 47.486792 Diff=-6.30D-08 RMSDP= 5.79D-06. It= 13 PL= 2.30D-06 DiagD=F ESCF= 47.486790 Diff=-1.52D-07 RMSDP= 5.45D-07. Energy= 0.174514138371 NIter= 14. Dipole moment= 1.182830 0.500037 0.671478 (...)

When you run frequencies, it will look a bit different, such as here (where the total energy is 0.103689468239 a.u.:

(...) SE2nd: IAtom= 19 IXYZ=1 IS=2. SE2nd: IAtom= 19 IXYZ=2 IS=1. SE2nd: IAtom= 19 IXYZ=2 IS=2. SE2nd: IAtom= 19 IXYZ=3 IS=1. SE2nd: IAtom= 19 IXYZ=3 IS=2. Maximum difference in off-diagonal FC elements: I= 16 J= 3 Difference= 1.3800844072D-04 Max difference between analytic and numerical forces: I= 17 Difference= 6.1396859930D-05 Energy= 0.103689468239 NIter= 17. Dipole moment= 0.078083 -0.081208 -0.007392 (...)

- For a Hartree-Fock job, there will be a section such as the following (the energy is obvious in this case):
The output for an MP2 job is unusual, in that these jobs have two main steps: first, the Hartree-Fock wavefunction is determined, then the MP2 correction is computed. You will accordingly see output similar to the above, followed by the following, which gives the energy:

(...) Spin components of T(2) and E(2): alpha-alpha T2 = 0.5486066740D-02 E2= -0.2458707205D-01 alpha-beta T2 = 0.3482713231D-01 E2= -0.1361644265D+00 beta-beta T2 = 0.6522273303D-02 E2= -0.2267776322D-01 (S**2,0)= 0.83038D+00 (S**2,1)= 0.79224D+00 E(PUHF)= -0.10047163912D+03 E(PMP2)= -0.10065203695D+03 ANorm= 0.1023149780D+01 E2 = -0.1834292618D+00 EUMP2 = -0.10064225142696D+03 (...)

In this example, the energy (EUMP2) is -100.64225142696.

**Be careful**not to be confused ! It is given in scientific notation, it is**NOT**the same as the HF energy, and there is a third energy E(PMP2) which is similar but it too is not the one you want.

(...)
Requested convergence on energy=5.00D-05.
Keep R1 integrals in memory in canonical form, NReq= 468559.
SCF Done: E(RHF) = -77.5985276906 A.U. after 5 cycles
Convg = 0.7567D-04 -V/T = 2.0034
S**2 = 0.0000
(...)

After this section, the output will differ depending on whether you have asked for a geometry optimisation, a single-point, or a frequency job.

For a single-point job, or for the first step of a geometry optimisation, you will have some information on the wavefunction, e.g. the Mulliken charges, the overall dipole, etc. For the work on the Diels-Alder reaction, you will be interested in orbital energies. These are given in atomic units (Hartrees) and appear in the output in this region, in the following form:

(...)
The electronic state is 1-A'.
Alpha occ. eigenvalues -- -1.27277 -1.09035 -0.82094 -0.65970 -0.60769
Alpha occ. eigenvalues -- -0.49878 -0.49658 -0.46789 -0.45542 -0.40803
Alpha virt. eigenvalues -- 0.00072 0.05192 0.05756 0.07565 0.10793
Alpha virt. eigenvalues -- 0.12388 0.16303 0.18013 0.24643
(...)

In this example, the energy of the Highest Occupied MO is -0.40803, whilst that of the Lowest Unoccupied (=virtual) MO is 0.00072.

If the job is not a geometry optimisation, this will more or less be the end of the output. If all has gone well, the ultimate lines of the logfile will be something like this:

(...)
Job cpu time: 0 days 0 hours 0 minutes 13.3 seconds.
File lengths (MBytes): RWF= 10 Int= 0 D2E= 0 Chk= 5 Scr= 1
Normal termination of Gaussian 98.

Otherwise, some or other error has occurred.

For a geometry optimisation, the energy gradient will be computed, and the geometry updated. There will be a convergence test, the results of which appear in the output as follows:

(...)
Item Value Threshold Converged?
Maximum Force 0.080863 0.000450 NO
RMS Force 0.048414 0.000300 NO
Maximum Displacement 0.181555 0.001800 NO
RMS Displacement 0.173872 0.001200 NO
(...)

All the four criteria need to be converged for the job to end. If they are not, there will be a new energy computation at the new geometry, and so on. As the optimisation proceeds, the Values should all decrease, and should start to drop below the Thresholds. Some of the 'NO's will be replaced by 'YES's. After the geometry has converged, there will be a message such as:

(...)
Maximum Force 0.000048 0.000450 YES
RMS Force 0.000013 0.000300 YES
Maximum Displacement 0.000999 0.001800 YES
RMS Displacement 0.000184 0.001200 YES
Predicted change in Energy=-1.137258D-08
Optimization completed.
-- Stationary point found.
----------------------------
! Optimized Parameters !
! (Angstroms and Degrees) !
------------------------ -------------------------
! Name Definition Value Derivative Info. !
-----------------------------------------------------------------------------
! R1 R(1,6) 1.3724 -DE/DX = 0. !
! R2 R(1,9) 1.087 -DE/DX = 0. !
! R3 R(1,14) 1.0957 -DE/DX = 0. !
! R4 R(2,3) 1.3856 -DE/DX = 0. !
(...)

The final geometry will be presented in several ways, and will be followed by the wavefunction information as mentioned above.

For a frequency job, or an optimisation + frequency job as discussed earlier, there will be a further section to the output, summarising the results of the frequency run. This will look like this:

(...)
Full mass-weighted force constant matrix:
Low frequencies --- -0.0010 -0.0008 0.0006 114.6950 189.0486 193.7753
Low frequencies --- 929.2244 1074.4697 1116.6974
Harmonic frequencies (cm**-1), IR intensities (KM/Mole),
Raman scattering activities (A**4/AMU), Raman depolarization ratios,
reduced masses (AMU), force constants (mDyne/A) and normal coordinates:
1 2 3
BU AU BG
Frequencies -- 929.2244 1074.4697 1114.7150
Red. masses -- 1.0406 1.1607 1.4837
Frc consts -- 0.5294 0.7895 1.0862
IR Inten -- 1.0153 142.6862 0.0000
Raman Activ -- 0.0000 0.0000 7.4672
Depolar -- 0.0000 0.0000 0.7500
Atom AN X Y Z X Y Z X Y Z
1 6 -0.04 0.00 0.00 0.00 0.00 -0.08 0.00 0.00 0.15
2 1 0.23 -0.44 0.00 0.00 0.00 0.50 0.00 0.00 -0.49
3 1 0.23 0.44 0.00 0.00 0.00 0.49 0.00 0.00 -0.49
4 6 -0.04 0.00 0.00 0.00 0.00 -0.08 0.00 0.00 -0.15
5 1 0.23 -0.44 0.00 0.00 0.00 0.50 0.00 0.00 0.49
6 1 0.23 0.44 0.00 0.00 0.00 0.49 0.00 0.00 0.49
4 5 6
AU AG AG
Frequencies -- 1140.4661 1385.7888 1434.9169
Red. masses -- 1.0078 1.4914 1.7488
(...)
VIBRATIONAL TEMPERATURES: 1336.94 1545.91 1603.82 1640.87 1993.83
(KELVIN) 2064.51 2305.24 2471.81 4795.66 4817.39
4933.10 4974.48
Zero-point correction= 0.054602 (Hartree/Particle)
Thermal correction to Energy= 0.057575
Thermal correction to Enthalpy= 0.058519
Thermal correction to Gibbs Free Energy= 0.033080
(...)

This gives you the symmetry assignment of each vibrational mode, its frequency (in wavenumbers), its IR and Raman intensity, the motions of the atoms to which it corresponds, etc. It also gives you the total Zero-point energy, in atomic units, and lots of other information besides...

Last modified: 10/01/2014 | Back to Top | Go to Mezei Lab Home page |