Force Field Calculations in Gaussian 03

Despite the fact that Gaussian is a program for quantum mechanical calculations, it also offers some possibilities for performing force field calculations with the AMBER, Dreiding, and UFF force fields. In contrast to quantum mechanical calculations, the input files for these types of calculations are slightly more complex as the force field calculations require clear specifications of atom type and connectivity. The following is an input file for the optimization of methanol with the AMBER force field:
#P amber opt=Z-Matrix geom=connectivity

AMBER opt methanol, explicit connectivity

0  1
H-H1-
C-CT-  1  r2
O-OH-  2  r3  1  a3
H-HO-  3  r4  2  a4  1  180.0
H-H1-  2  r5  3  a5  1  d5
H-H1-  2  r5  3  a5  1  -d5

r2=1.08105957
r3=1.39956997
r4=0.94629421
r5=1.08744287
a3=107.170012
a4=109.447010
a5=112.035787
d5=118.773010

1 2 1.0
2 3 1.00 5 1.0 6 1.0
3 4 1.0
4
5
6
      

The definition of the geometry of the system differs from that used in quantum mechanics calculations in two points:

1) The atom types are defined through atomic symbols as well as a force field atom type. The first line of the above input:

H-H1-

defines a hydrogen atom which is connected to a carbon atom carying one electronegative substituent. The same atom type is chosen for hydrogen atoms H5 and H6. A different atom type -HO is selected for hydrogen atom 4 connected to the oxygen atom. The atom type definitions are those defined for the AMBER force field. A full list of possible atom types can be found here.

2) The input file now also includes information on the atom connectivity. Prompted by the geom=connectivity keyword the program reads after the end of the Z-Matrix definition additional lines (one per atom) that define the bonds present in the system. The line:

1 2 1.0

defines a bond of bond order 1.0 between centers 1 and 2. This type of information is given once for every bond contained in the system. Based on a proximity criterion the program can also define bonds by itself. This does, of course, only work properly if the starting geometry is chosen well enough.
#P amber opt=Z-Matrix

AMBER opt methanol, no connectivity, no charges

0  1
H-H1-
C-CT-  1  r2
O-OH-  2  r3  1  a3
H-HO-  3  r4  2  a4  1  180.0
H-H1-  2  r5  3  a5  1  d5
H-H1-  2  r5  3  a5  1  -d5

r2=1.08105957
r3=1.39956997
r4=0.94629421
r5=1.08744287
a3=107.170012
a4=109.447010
a5=112.035787
d5=118.773010
      

The Gaussian output file contains rather limited information on the geometries generated during the geometry optimzation if only the opt keyword is being used. This has to do with the original purpose of the implementation of force field methods in Gaussian for combination with quantum mechanical methods. In these types of QM/MM methods several microiterations are performed for the force field subsystem before one geometry optimization step is performed for the quantum mechanical system. Only after the latter a new geometry is written to the output file. A full reflection of geometry changes during optimization can be achieved using either the opt=nomicro keyword or the definition of a particular coordinate system with the opt=Z-Matrix or opt=Redundant keywords.

For each structure along the optimization pathway Gaussian reports the energy of the system as:
(Enter /scr1/g98/l402.exe)
 AMBER calculation of energy and first derivatives.
 MO and density RWFs will be updated.
 Recover connectivity from rwf.
 EGHMM:  NBAlg=1 ICut=1 CutNB= 1.43D+03
 Energy=    0.000793358514 NIter=   0.
The energy listed here is that defined through the bonded and non-bonded interactions present in the system, the hypothetical "unstrained" state serving as the zero point of energy.

One additional aspect of force field calculations concerns the calculation of Coulomb interactions between non-bonded centers. If not specified in the input file, no charges are assigned to the atoms of the system and no Coulomb interactions are calculated. If, however, charges are included in the input file, the Coulomb interactions are calculated. This also changes the definition of the force field itself and the reported energy of the system is not directly comparable anymore to that without Coulomb interactions. The following input file for methanol contains charges derived with the CHELPG scheme at the HF/6-31G(d) level of theory:

#P amber opt=nomicro

AMBER opt methanol, no connectivity, CHELPG charges

0  1
H-H1-0.04
C-CT-0.36   1  r2
O-OH--0.49  2  r3  1  a3
H-HO-0.17   3  r4  2  a4  1  180.0
H-H1--0.04  2  r5  3  a5  1  d5
H-H1--0.04  2  r5  3  a5  1  -d5

r2=1.08105957
r3=1.39956997
r4=0.94629421
r5=1.08744287
a3=107.170012
a4=109.447010
a5=112.035787
d5=118.773010
      

Please observe that the charge parameters directly follow the atom type specification without a blank.

The process of building an appropriate input file for force field calculations can be vastly facilitated using GaussView. For this purpose we will use GaussViews facility to generate complete input files for any type of calculation. After designing the system of interest using, for example, the fragment library, select Gaussian... from the Calculate pull down menu. In the popup window select Mechanics.. in the Method: panel and Amber as the force field model. After checking the Write Connectivity button in the General menue, you can generate the input file selecting the Edit.. and the Save buttons and provide a file name. The generated file then appears in a new window, which can be closed without any changes. The subsequently appearing submit prompt can be canceled. This type of aborted job submission generates the file with all necessary connectivity information and atom type declarations. If contained in the fragment library, atomic charges are also included.

last changes: 09.11.2006, HZ
questions & comments to: zipse@cup.uni-muenchen.de