47001-AC6
Molecular Dynamics Study of Stability and Meta-Stability of Methane Hydrates
In the first year of research under the support of the PRF award, the research group led by the PI has concentrated on developing the adaptive force matching (AFM) method as outlined in the proposal. A code implementing the AFM procedure has been written by post-doctoral fellow, Omololu Akin-ojo. Omololu and graduate student Yang Song have made significant contribution to this project.
The AFM method is composed of three steps. The first step is the MD step, in which MD simulation is carried out using a guess force field to sample the configuration space. After equilibrium has been reached in the MD step, QM/MM calculations are performed in constant intervals (e.g., every 10fs) by treating part of the system at the QM level, while the rest of the system by the guess force field. This is constitutes the QM/MM step. After the QM/MM step, the force-matching (FM) technique is applied to reparametrize the guess force field to minimize the difference between the AFM and the QM forces in the FM step. The reparameterized force field will be used as the new guess force field and the AFM will start from the MD step again. The procedure will be repeated until convergence is reached. Taking a guess force field, the AFM method is capable of improving it to the quality of high level ab initio calculations in a few iterations (TOC figure). We have been successful in creating water forces field following the AFM procedure using several different ab initio methods as reference and the work has been published recently on the Journal of Chemical Physics. The paper was ranked the 4th most downloaded article on JCP for the month of August, 2008. In addition to the development of water force fields, significant progress has been made in developing a high quality potential for the methane clathrate system following the AFM approach. Recently, we found the water potential can be significantly improved by incorporating a new hydrogen bond term in the water model following the AFM procedure. We are writing up a second publication reporting our recent finding in the importance of this hydrogen bond term.
Another important direction outlined in the proposal is the investigation of the free energy of formation of methane clathrate using the non-equilibrium pulling approach with Jarzynski’s identity. Graduate student Gerrick Lindberg has incorporated the non-equilibrium pulling method into the software package DL_POLY. We are planning to first reproduce the published melting temperature of TIP4P water as a test of our implementation. In doing so, we realized the commonly accepted procedure for creating ice configurations is slightly biased. In the commonly accepted procedure for generating ice configurations, the total dipole moment of the box is constraint to zero. This leaves out other important configurations as configurations with a finite dipole moment contribute to the equilibrium ensemble according to the fluctuation-dissipation theorem. In order to properly create ice configurations sampling the proton reorientation degree of freedom, we developed the electro-static switching method to quickly generate equilibrium ice configuration, the method can properly sample the configuration space of ice. The dielectric constant of the ice has been investigated by measuring the dipole fluctuation of the ensemble of ice generated using our procedure. The result has been compared with previously published work. Not only a good agreement is achieved, the error bar produced by our method is significantly smaller. The method has been published as a journal article on the Journal of Physical Chemistry B.
Currently, we are proceeding to obtain the melting free energy of ice. When we are investigating this problem, we realized the pulling function in the non-equilibrium pulling calculations can be optimized to provide a tighter bound to the free energy difference for a given pulling time. The optimal pulling function has been investigated by a few groups. However, as it was formalized previously, the determination of the pulling function takes about the same effort as determining the free energy difference between the two states. We are working on an on-the-fly approach that can approximate this optimal pulling function with minimal additional computational effort. Very good results have been achieved for our new on-the-fly approach. Gerrick and I are preparing a manuscript to be submitted to the Journal of Chemical Physics reporting our effort.
