Chapter 1 Introduction 1 General Introduction

Finally, we demonstrate how the multiple time step QM/MM technique can accelerate the equilibration process and ultimately configurational sampling compared to a standard simulation in which there is no mass rescaling. Two simulations of normal undecane(C11) have been performed where the terminal methyl groups make up the QM region as illustrated in Figure 5.14.^◊ One simulation was performed with the standard Verlet propagator where "standard" masses for the whole system were utilized (12.0 amu for C and 1.5 amu for H). In the other simulation, the standard masses were rescaled in the MM region by a factor of 1/400 and the multiple time step algorithm was used with an over-sampling ratio of 20:1. The masses used in the multiple time-step simulation are detailed in Figure 5.14 for the carbon atoms of the backbone(hydrogen atoms were similarly rescaled). A high energy conformation of the hydrocarbon chain was selected as the initial structure.^§ For both simulations, the free, unthermostated dynamics was commenced from this structure with no initial velocities. Additionally, the constraint shown in Figure 5.14 was applied to the system. Here the distance between one of the terminal carbon atoms and a proton of the other terminal carbon atom was constrained to a distance of 4 Å throughout the dynamics (The reason for applying the constraint will be explained later)

The temperature evolution of the QM and MM subsystems for both the standard and the multiple time-step simulation is shown in Figure 5.15. Since the dynamics was initiated from a high energy conformation, there is a steep rise in the temperature of both the QM and MM systems in the first 100 fs of the simulation. When the whole system is properly equilibrated, both the QM and MM systems should be at the same temperature. For the multiple time-step simulation the temperature evolution plotted in Figure 5.14b reveals that the QM and MM systems both attain a steady state temperature of approximately 400 K. The temperature of the two systems is well equilibrated at about the 800 fs mark of the simulation. In contrast, for the standard Verlet simulation with no mass rescaling, equilibration does not occur in the first 1400 fs of the simulation. Even at the 1400 fs mark, the fluctuations are large and disparity between the QM and MM temperatures is still rather large. The continuation of the temperature evolution for the standard time step simulation is plotted in Figure 5.16. The plot reveals that equilibration does not begin to occur until the 5500 fs mark where an equilibrium temperature of 400 K is being established. However, even after 8000 fs temperature equilibration has not been fully established.

The distance constraint was imposed during the dynamics of undecane in order to evaluate if the multiple time step method would be effective in accelerating the convergence of the constraint force. The convergence of the force is important since the slow growth simulations we perform to map out reaction free energy profiles involves the integration of the constraint force. A faster convergence of the constraint force would then correspond to an acceleration in the configurational averaging in these simulations and also smaller errors. Although the undecane system doesn't correspond to a real reaction, the simulation does imitate the condition where the structure of the MM subsystem strongly influences the force on the reaction coordinate. In the undecane simulation, the constraint force involving the distance between the terminal methyl groups will be regulated by the structure of the hydrocarbon backbone. Thus, the faster the structure of the back bone equilibrates, the faster the constraint should converge. Since most of the backbone resides in the MM partition, we expect that the QM/MM multiple time-step procedure will accelerate the convergence of the force compared to the standard simulation.