The usage of all three scripts is: scriptname.scr dist1 dist2 steps where dist1 and
dist2 are the first and last points along the reaction coordinate
respectively and steps is the number of points we want to have on the
reaction coordinate.

We also need to write a QM/MM topology file for this
compound. See for example topol_A.itp. Have a look at this file
befroe we proceed. For the linear transit we will need a new topology
with a different constraint distance for each point of the Linear
transit. In our case the scripts mentioned above will take care of
this, which is why we can skip this part.

Furthermore we need an index file, with entries for the
atoms that are supposed to be part of the MM and QM subsystems. We can
use make_ndx for that, but we will skip this part as well, and use the
prepared index file instead. index

Finally, we create an mdp file with the options for the
QM/MM linear transit calculation. We need to specify that we want to
energy minimize, how we want to minimize, that we want to do QM/MM,
what QM method we want to use. The mdp file we are going to use is LT.mdp.

With the input files and the scripts at hand, we can now
perfom the actual linear transit calculation. We will let the reaction
coordinate constraint vary between 0.12 and 0.5 nm in 200 steps.

execute the create_tops.scr scripts to create 201 subdirectories
(called step_0, ..., step_200) and create a topol_A.itp with different
constraints lengths in the range 0.12 to 0.5 nm.

localhost:~>./create_tops.scr 0.12 0.5 200

execute the run.scr scripts to run the minimizations in the
different subdirectories. To speed up the convergence, the script
takes the output coordinates of the previous Linear Transit point as
input in the current minimization

localhost:~>./run.scr 0.12 0.5 200

This will take a while, depending on the speed of your computer
system.

execute the get_ener.scr script to retrieve the energies from the
individual linear transit points. Redirect the output to a file, let's
call it eqmmm.xvg.

localhost:~>./get_ener.scr 0.12 0.5 200 >
eqmmm.xvg

Let's have a look at the results

Import eqmmm.xvg in xmgr to see the energy as a function of the
reaction coordinate.

localhost:~>xmgr eqmmm.xvg

The maximum is at 0.2093 nm (-149.924 kJ/mol). This corresponds to
point 47. Let's have a look at that structure. We simply go into the
directory step_47 and use editconf to convert the confout.gro into
confout.pdb

localhost:~>editconf -f confout.gro -o confout.pdb

We then use rasmol to visualize the structure:

localhost:~>rasmol confout.pdb

It looks verly much like the one we optimized in the previous
step. You can check if it is a TS by perfoming a frequencies
caclulaiton on this structure, using gaussian, just like we did before.

The minima are at 0.4468 nm (-342.826kJ/mol) and 0.1485 (-388.276
kJ/mol) and correspond to the reactant and product state respectively,
and are shown here, along with the transition state again.

Figure 5. Reactant (a), Transition State (b),
and Product (c) geometries in vacuo, found with the Linear
Transit method. Theenergies of these structures are listed in table
2.