QM/MM tutorial


QM/MM calculations on a Diels-Alder antibody catalyst.


V. Conclusions, Discussion and Outlook

We have performed a number of QMMM calculations on a Diels-Alder cycloaddition. We have seen that the reaction barrier decreases when going from vacuum to water and from water to the catalytic antibody. The conclusion so far are that the protein catalyzes the cycloaddition by stabilizing the transition state with respect to the reactant state. See table 5 and figure 9 for details. However, there a few critical points that have not been addressed.
Table 5. Relative Energies of the
reactant, transition state and product
geometries in vacuo, water and the
protein at the PM3/GROMOS96
QM/MM level. The energies are
reported in kJ/mol and are taken
relative to the potential energy of the
reactants.

vacuowaterprotein
Reactant0.00.00.0
Trans. St.192.902139.5113.5
Product-45.45-69.3-144.4

Figure 9. Energy profiles for the Diels Alder cycloaddition in
vacuum, in water and in the active site of the catalytic antibody,
calculated at the PM3/GROMOS QM/MM level.

First of all, we have used a Steepest Descent algorithm to do the energy minimizations. This algorithm, however, is known not to converge very well near a minimum. This was one of the reasons for increasing the convergence criteria. A better algortithm for minimization would be the BFGS algorithm, which employs higer order derivatives, but is much more expensive in terms of computational costs. It would give more accurate results though, which is the primary reason for using it.

Second, the results happen to be very sensitive with respect to the initial conditions. Slighly different starting configurations would result in different potential energy curves. This is another consequence of using the simple Steepest Descent method for minimization.

Third, real reaction barriers are Free Energy barriers. With the current setup for the Linear Transit is is possible to do a free energy calculation. We need to define a State A and a State B, where A and B represent the reactant and product configuration respectively. The constraint between the Dummy atoms in state A would be 0.14 and in state B 0.3. To actually perform a QM/MM Free Energy calculation, one needs to specify both the A and B state parameters in the constraint section of the topology file:

[ constraints ]
; atom1 atom2 type stateA stateB
dummy1 dummy2 2 0.14 0.3
Furthermore, one needs to 'tell' gromacs it is supposed to do a free energy perturbation calculation by addig the lines (see
gmx manual).
free_energy = yes
init_lambda = 0
delta_lambda = 0.01
to the mdp file. The calculations will be more time consuming, but the result is the free energy curve of the reaction, which is easier to relate to experimental data than the potential energy curves we computed thus far. But the tutorial was aimed to give a more qualitative picture of the reaction rather than a quantitative one.

We hope you enjoyed doing this tutorial and that you found it useful for your own work now or in the future. Of course, the system used in the tutorial was a rather easy starting point as the transition state analogue was known, but the techniques you learned should work without such knowledge as well. If, while performing QM/MM calculations you run into trouble, or would like to discuss something, you can always contact me.

Previous:IV. Optimization of product, reactant and transition state geometries in the fully solvated protein, using Linear Transit in gromacs

back to top


updated 07/09/04