MOLLY Integrators for Ensembles other than Constant Energy

To facilitate the comparison of calculated results with experiments, it is desirable to perform simulations in canonical (constant temperature T) and isothermal-isobaric (constant temperature T and pressure P) ensembles. For example, in studying dilute solutions, it is worth while to simulate both the pure solvent and the dilute solution at the same T and P. This corresponds to the usual experimental situation in which partial molar quantities of the solute are measured. Examples of such systems which particularly interest us include protein membrane simulations, solvated proteins, and quantum-dot molecules. In this part of the project we propose to develop and implement MOLLY integrators for constant temperature and/or pressure molecular dynamics. The main objective of this part is two fold: (i) to increase the step size used for integration during simulations and (ii) to implement ensembles.

For constant T and/or P simulations, Anderson [7] proposed a hybrid of MD and MC methods. In his approach, the particles change their velocities by stochastic collisions. The distribution of the velocities of the particles that collide is chosen to reproduce the canonical ensemble. Because of the sudden change of the velocities by collisions, the trajectory in the phase space becomes discontinuous. Nosé [67,68,116] achieved a major advance by showing that the canonical distribution can be generated with smooth deterministic and time reversible trajectories. He introduced an additional degree of freedom which acts as an external system for the physical system and thus, the system is simulated using an extended system Hamiltonian. Recently, Nosé's derivation was corrected to make it symplectic [23]. Klein and coworkers [102,103,105,145,148] by applying the Verlet method to the Trotter factorization of the Liouville operator implemented canonical and isobaric-isothermal ensembles to extended system Hamiltonian for MTS. This Verlet based MTS method is limited by numerical considerations and the longest time step is limited by resonance. They used a time step of 12fs to simulate liquid pseudoatom tetradecane.

In our earlier work [76], we introduced two new MTS integrators, LM and BBK-M. Our calculations showed that long time step of 14fs are possible by using LM with mild damping. When using mild damping, LM was found to be superior to other methods because it is symplectic in the zero damping limit. When larger damping coefficients were used, LN [15,14] was superior, with LM closely following. In this part of the project we will further improve the LM method. The main questions we will address here are: what is the longest time step possible for LM method and what improvements can be made to further lengthen the time step without sacrificing accuracy and stability.

We also propose to implement extended system Hamiltonian to simulate canonical and isothermal-isobaric ensembles using multiscale MOLLY methods that have been shown to use longer time step than 12fs [76]. This project will be implemented in three steps. The first step will be the development of equations of motion for extended system Hamiltonian that uses MOLLY integrators for constant T, in the second step, equations of motion for constant P will be developed and finally, we will combine equations of motion for constant T and P. At each step the algorithm, thus developed will be implemented in PROTOMOL. This effort will help simulate canonical and isobaric-isothermal ensembles. We will also compare the quality of dynamics generated using both Langevin and extended system molecular dynamics. Another possibility is to combine extended system Hamiltonian and stochastic damping to further increase the length of the time step used during simulations.