# How do I perform a partial optimisation in GAMESS?

Sometimes it is necessary to perform a partial optimisation of a molecule, i.e. keeping certain variables constant. One example might be to pre-optimise a transition state, taking advantage of not having to compute the force constants in the initial set-up.

For example, let's consider the transition state of the Diels-Alder Reaction of ethene and butadiene, where we want to freeze the distance between the two molecules.

How can I achieve this in the software suite GAMESS?

Unfortunately GAMESS is not as user friendly as other quantum chemistry packages, so when it comes to running non-standard calculations a little bit more effort is necessary.

First of all I recommend optimising the reactants by themselves, so that you have a somewhat valid idea what the bond lengths on your level of theory are going to be. After that you should assemble your best guess for the structure. In this example I have chosen a separation between the molecules of about 220 pm.

The coordinates of this structure are:

! DIENE
C       -0.486294808     -0.431915730     -1.395523790
C       -1.166676574      0.550508040     -0.705577944
C       -1.166676574      0.550508040      0.705577944
C       -0.486294808     -0.431915730      1.395523790
H       -0.352988862     -0.347708199      2.465648783
H       -0.084717264     -1.289523730      0.864198908
H       -1.712490163      1.338624758      1.223514576
H       -1.712490163      1.338624758     -1.223514576
H       -0.352988862     -0.347708199     -2.465648783
H       -0.084717264     -1.289523730     -0.864198908
! DIENOPHILE
C        1.559305232     -0.085874351     -0.690651129
C        1.559305232     -0.085874351      0.690651129
H        1.514911897     -1.022640981      1.232100569
H        1.603704743      0.851022589      1.232175888
H        1.603704743      0.851022589     -1.232175888
H        1.514911897     -1.022640981     -1.232100569


There are approaches, where you would simply fix the Cartesian coordinates. Here we only want to fix two bond lengths, 1(C)-11(C) and 4(C)-12(C). If we were to fix the coordinates of these atoms, we would also fix 1(C)-4(C) and 11(C)-12(C). Therefore a better approach is to use internal coordinates, in GAMESS the keyword group $ZMAT is used. Since there are a lot, it is tricky to set them up. We will therefore use the implemented automated procedure. Starting with a dry run, i.e. EXETYP=CHECK, so that we can check what the program actually does. We need to define$3\times N-6=42$internal coordinates, thus NZVAR=42. Also note, that we do not use symmetry in this case. Because of the nature of the structure we need to add some bonds already manually, NONVDW(1)=1,11, 4,12. This is the resulting input: $CONTRL
SCFTYP=RHF
RUNTYP=ENERGY
EXETYP=CHECK
COORD=UNIQUE
ICHARG=0
MULT=1
NZVAR=42
$END$BASIS
GBASIS=STO
NGAUSS=3
$END$ZMAT
AUTO=.TRUE.
DLC=.TRUE.
NONVDW(1)=1,11, 4,12
$END$GUESS
GUESS=HUCKEL
$END$DATA
C6H10
C1
CARBON      6.0     -0.486294808        -0.431915730        -1.395523790
CARBON      6.0     -1.166676574         0.550508040        -0.705577944
CARBON      6.0     -1.166676574         0.550508040         0.705577944
CARBON      6.0     -0.486294808        -0.431915730         1.395523790
HYDROGEN    1.0     -0.352988862        -0.347708199         2.465648783
HYDROGEN    1.0     -0.084717264        -1.289523730         0.864198908
HYDROGEN    1.0     -1.712490163         1.338624758         1.223514576
HYDROGEN    1.0     -1.712490163         1.338624758        -1.223514576
HYDROGEN    1.0     -0.352988862        -0.347708199        -2.465648783
HYDROGEN    1.0     -0.084717264        -1.289523730        -0.864198908
CARBON      6.0      1.559305232        -0.085874351        -0.690651129
CARBON      6.0      1.559305232        -0.085874351         0.690651129
HYDROGEN    1.0      1.514911897        -1.022640981         1.232100569
HYDROGEN    1.0      1.603704743         0.851022589         1.232175888
HYDROGEN    1.0      1.603704743         0.851022589        -1.232175888
HYDROGEN    1.0      1.514911897        -1.022640981        -1.232100569
$END  The resulting auxiliary outputfile .dat has now the group that we will use from here on in. And the .log file will give us the exact value of the internal coordinate which we want to freeze.  -------------------- INTERNAL COORDINATES -------------------- - - ATOMS - - COORDINATE COORDINATE NO. TYPE I J K L M N (BOHR,RAD) (ANG,DEG) ---------------------------------------------------------------- 1 STRETCH 1 11 4.1406435 2.1911343 ... 10 STRETCH 4 12 4.1406435 2.1911343 ...  We can now modify our $ZMAT group accordingly. Unfortunately the input gets very long, if you would like to retain some readability. Note that we create redundant internal coordinates, so we need to increase the number of NZVAR, too, to read all of them in.
Technically we do not have to set the generation of coordinates to false, but it helps remembering what you have done if you come back in a month or so, therefore AUTO=.FALSE.. Now we add the coordinates we would like to freeze with IFZMAT(1)=<list> and the corresponding values FVALUE(1)=<list>.

 $CONTRL SCFTYP=RHF RUNTYP=OPTIMIZE COORD=UNIQUE ICHARG=0 MULT=1 NZVAR=89$END
$BASIS GBASIS=STO NGAUSS=3$END
$ZMAT AUTO=.FALSE. DLC=.TRUE. IZMAT(1)= 1, 1, 11, 1, 1, 10, 1, 1, 9, 1, 1, 2, 1, 2, 8, 1, 2, 3, 1, 3, 7, 1, 3, 4, 1, 4, 6, 1, 4, 12, 1, 4, 5, 1, 11, 16, 1, 11, 15, 1, 11, 12, 1, 12, 14, 1, 12, 13, 2, 1, 11, 15, 2, 1, 11, 16, 2, 1, 2, 3, 2, 1, 2, 8, 2, 1, 11, 12, 2, 2, 3, 7, 2, 2, 1, 10, 2, 2, 1, 11, 2, 2, 1, 9, 2, 2, 3, 4, 2, 3, 2, 8, 2, 3, 4, 5, 2, 3, 4, 12, 2, 3, 4, 6, 2, 4, 3, 7, 2, 4, 12, 13, 2, 4, 12, 14, 2, 4, 12, 11, 2, 5, 4, 12, 2, 5, 4, 6, 2, 6, 4, 12, 2, 9, 1, 10, 2, 9, 1, 11, 2, 10, 1, 11, 2, 11, 12, 13, 2, 11, 12, 14, 2, 12, 11, 15, 2, 12, 11, 16, 2, 13, 12, 14, 2, 15, 11, 16, 3, 1, 2, 3, 4, 3, 1, 2, 3, 7, 3, 1, 11, 12, 14, 3, 1, 11, 12, 4, 3, 1, 11, 12, 13, 3, 2, 1, 11, 15, 3, 2, 1, 11, 16, 3, 2, 3, 4, 5, 3, 2, 3, 4, 12, 3, 2, 3, 4, 6, 3, 2, 1, 11, 12, 3, 3, 2, 1, 11, 3, 3, 4, 12, 13, 3, 3, 4, 12, 14, 3, 3, 4, 12, 11, 3, 3, 2, 1, 9, 3, 3, 2, 1, 10, 3, 4, 12, 11, 15, 3, 4, 12, 11, 16, 3, 4, 3, 2, 8, 3, 5, 4, 12, 11, 3, 5, 4, 12, 14, 3, 5, 4, 12, 13, 3, 5, 4, 3, 7, 3, 6, 4, 12, 11, 3, 6, 4, 12, 14, 3, 6, 4, 12, 13, 3, 6, 4, 3, 7, 3, 7, 3, 2, 8, 3, 7, 3, 4, 12, 3, 8, 2, 1, 9, 3, 8, 2, 1, 10, 3, 8, 2, 1, 11, 3, 9, 1, 11, 16, 3, 9, 1, 11, 15, 3, 9, 1, 11, 12, 3, 10, 1, 11, 12, 3, 10, 1, 11, 15, 3, 10, 1, 11, 16, 3, 13, 12, 11, 15, 3, 13, 12, 11, 16, 3, 14, 12, 11, 15, 3, 14, 12, 11, 16, IFZMAT(1)= 1, 1, 11, 1, 4, 12, FVALUE(1)= 2.1911343, 2.1911343,$END
$GUESS GUESS=HUCKEL$END
$DATA C6H10 C1 CARBON 6.0 -0.486294808 -0.431915730 -1.395523790 CARBON 6.0 -1.166676574 0.550508040 -0.705577944 CARBON 6.0 -1.166676574 0.550508040 0.705577944 CARBON 6.0 -0.486294808 -0.431915730 1.395523790 HYDROGEN 1.0 -0.352988862 -0.347708199 2.465648783 HYDROGEN 1.0 -0.084717264 -1.289523730 0.864198908 HYDROGEN 1.0 -1.712490163 1.338624758 1.223514576 HYDROGEN 1.0 -1.712490163 1.338624758 -1.223514576 HYDROGEN 1.0 -0.352988862 -0.347708199 -2.465648783 HYDROGEN 1.0 -0.084717264 -1.289523730 -0.864198908 CARBON 6.0 1.559305232 -0.085874351 -0.690651129 CARBON 6.0 1.559305232 -0.085874351 0.690651129 HYDROGEN 1.0 1.514911897 -1.022640981 1.232100569 HYDROGEN 1.0 1.603704743 0.851022589 1.232175888 HYDROGEN 1.0 1.603704743 0.851022589 -1.232175888 HYDROGEN 1.0 1.514911897 -1.022640981 -1.232100569$END


The optimisation should look something like the following animation; and took 19 steps in my example.

Now you could use the results of this calculation to run a full fledged transition state search. Remember, that you have all necessary information to restart a run in the auxiliary output, i.e. the .dat file. You should always keep this one until you are completely done with your calculations.

If you want to calculate along and wonder what the optimised geometry is, here are the results of the final step:

----- RESULTS FROM SUCCESSFUL RHF      GEOMETRY SEARCH -----
----- COORDS, ORBS, GRADIENT, AND APPROX. HESSIAN -----
COORDINATES OF SYMMETRY UNIQUE ATOMS (ANGS)
ATOM   CHARGE       X              Y              Z
------------------------------------------------------------
CARBON      6.0  -0.4847330011  -0.4122786646  -1.4086312155
CARBON      6.0  -1.1902020846   0.5288839653  -0.7019679480
CARBON      6.0  -1.1902020983   0.5288839622   0.7019679173
CARBON      6.0  -0.4847329728  -0.4122786820   1.4086312270
HYDROGEN    1.0  -0.4076405265  -0.3274824211   2.4860530248
HYDROGEN    1.0  -0.4019303863  -1.4244097276   1.0396872383
HYDROGEN    1.0  -1.5467259917   1.4139143889   1.2140730563
HYDROGEN    1.0  -1.5467259527   1.4139143906  -1.2140730901
HYDROGEN    1.0  -0.4076405381  -0.3274824126  -2.4860529903
HYDROGEN    1.0  -0.4019304494  -1.4244096950  -1.0396872193
CARBON      6.0   1.5587691569  -0.0858507824  -0.6884831136
CARBON      6.0   1.5587691834  -0.0858507803   0.6884831288
HYDROGEN    1.0   1.9167366148  -0.9561957837   1.2206602417
HYDROGEN    1.0   1.6758814772   0.8473264236   1.2212666939
HYDROGEN    1.0   1.6758814263   0.8473264124  -1.2212666800
HYDROGEN    1.0   1.9167365785  -0.9561957779  -1.2206602275
--- OPTIMIZED RHF      MO-S --- GENERATED AT Thu Dec 14 23:42:52 2017
E=     -230.0376739583, E(NUC)=  227.2020852580
\$VEC
...


Disclaimer: There is another way of freezing coordinates. I might add that at another point as another answer. It is a weeny bit more complicated, however offers also some more flexibility.