QM/MM Calculations on Solvated Molecules

Principal Investigator

Peter L. Cummins

Division of Biochemistry and Molecular Biology

John Curtin School of Medical Research

The theoretical study of complex chemical systems consisting of many hundreds or thousands of atoms requires the use of various approximations for the interaction potential between atomic nuclei. These approximations are necessary because the exact, i.e. quantum mechanical (QM), solution rapidly becomes intractable as the number of atoms is increased. Consequently, the simpler molecular mechanics (MM) potentials are widely used in the simulation of large molecular aggregates and macromolecules to calculate thermodynamic properties or describe phenomena associated with the classical dynamics of the system. It only becomes necessary to use one of the QM methodologies in the calculations where, for example, chemical bond formation/breaking or electronic excitation processes are of interest, i.e. any process which involves a large change in the electronic structure. Since many of these processes are comparatively localised events, e.g. chemical bond formation, an obvious approach to the computational bottle-neck is to treat only a relatively small part of the system quantum mechanically. The remaining larger part of the system would be treated using molecular mechanics potentials.

In our approach to the combined QM and MM (QM/MM) methodology, we take advantage of the semiempirical QM approximation which is computationally more efficient than either the density functional or ab initio methods. This approach allows us to explore the possibility of generating ensemble averaged quantities using molecular dynamics simulation, thus facilitating the direct comparison of theory with experiment. The new QM/MM methods are now being used to gain greater insight into the chemical processes occurring in complex biomolecular systems.


Jill E. Gready

Division of Biochemistry and Molecular Biology

John Curtin School of Medical Research


u51 - VPP, PC

What are the results to date and the future of the work?

The work on the parameterization of the QM/MM interaction terms for small molecules in aqueous solution is complete. These QM/MM force-fields have now been applied to a solvated

protein system, dihydrofolate reductase (DHFR), in order to study the hydride-ion transfer step between the NADPH cofactor and substrate molecules.

The free energy change between reactant and product states for the hydride-ion transfer has been computed using molecular dynamics (MD) and free energy perturbation (FEP) methods. The semiempirical AM1 hamiltonian has been used for the QM system and the AMBER/TIP3P force field for the DHFR/solvent system. In particular we have looked at how the free energy change for the reaction depends on the choice of QM fragments and the FEP pathway. The results we have obtained demonstrate the feasibility of using combined QM/MM force-fields in MD simulations on enzymes. We have also begun preliminary studies into the use of free energy gradients (calculated numerically using an FEP-based procedure) as a computationally efficient means of locating local minima (reactant and product) and transition states.

Total MD trajectories of up to nanosecond (1000 ps) time scales are being calculated in order to provide the necessary equilibrations and statistical averages for free energy determinations using FEP methods. On the DHFR systems studied to date, we achieve 0.5 to 1 ps per cpu hour, depending largely on the size of the QM part of the system.We will need to characterise transition states as well as reactants and products. In addition, it will be necessary to calculate various free energy components in order to gain a detailed understanding of the nature of enzyme catalysis.

What computational techniques are used?

We use standard molecular dynamics (MD) simulation and free energy perturbation (FEP) methods, as implemented in Molecular Orbital Programs for Simulation (MOPS). The semiempirical AM1, MNDO and PM3 molecular orbital methods can be used to model the QM part of the system. Reaction coordinates are defined by imposing distance constraints in the MD simulations using the SHAKE algorithm.


Cummins, P.L., Gready, J.E., A Coupled Semiempirical Molecular Orbital and Molecular Mechanics Model (QM/MM) for Organic Molecules in Aqueous Solution, Journal of Computational Chemistry, in press.

- Appendix A