Advanced Kinetic Monte Carlo Simulations
PDFSimulation of Electrosis for the Production of Dihydrogen
Long Term ObjectivesOur objectives is to simulate the electrolysis of water for the production of hydrogen. To be able to simulate a large number of atoms, the simulation has be devided into two parts. One part is the simulation of the molecules which are near the electrode. These molecules may split the form dihydrogen. This process is a chemical reaction and has to be modelled by quantum mechanics. To simulate these molecules (near the electrode), we will use DFT calculation (Density Functional Theory). But, this is someone else project. My job is to work on the simulation of the water which is away from the electrode. The molecules in this area do not split so no quantum mechanics is required. Those molecules will be treated classically. Monte Carlo Simulation: Simulation of a very long Molecular Dynamics Simulation.We wish our simulation to be as close as possible to reality. For this we need to be able to simulate the water over long time scale. Using molecular dynamics simulation on todays computer capacities, we would probably be able to simulate at most 100 ns. We wish to simulate more. An alternative is to use Monte Carlo simulations. The idea of the Monte Carlo method is to simulate the Molecular Dynamics by making the system jump from minimum to minimum. The path followed is random but the probability for a minimum to be chosen as the next step is based on the barriers connecting to this minimum. The time spend by the system around one minimum is also chosen randomly based on the barriers surrounding the minimum. The monte carlo simulation requires to know the minima and saddle points of the system. The scanning of the landspace is presented in the next section. The Monte Carlo method we are using now is plain and we intend to program a more advanced algorithm. We hope that this algorithm will improve the efficiency of the simulation for systems which present a mix of low and high energy barriers.
|
 |
Scanning the energy landscape: to find the saddle points and minima
Starting from a minimum we need to climp up the potential in such a way that we shall end up to a saddle point. The idea is to follow the valley. A way to do this is to compute the Hessian (2nd derivatives of the potential), then compute the eigenvalues and eigenvectors of this matrix, then make a small step up the potential in the direction of the eigenvector which has the lowest eigenvalue. By repeating this process many times, we can expect to end up at a saddle point.
Calculating the hessian for large system is impossible. Instead we use other techniques to find the lowest mode.
Finding the Lowest Mode: Lanczos Iterative Method
The method I use to find the lowest mode is the Lanczos Iterative Method. The algorithm allows to build a matrix which approximates the Hessian but is much smaller. It is a iterative process which means that we start with a single scalar (i.e. 1×1 matrix) and step by step increases the size of the matrix. For very large system it is not possible to build the full matrix and consequently, it is not possible to obtain all the eigenvalues of full Hessian matrix. But the strength of the method is that the iterative process converges faster for the highest and smallest eigenvalues. As we are interested only in the lowest eigenvalue (and its eigenvector), it is possible to obtain a good approximation for the latter within a few iterations. For example, the system I am working on now is ice Ih. It has over 5000 degrees of freedom. It takes on average 3 iterations to obtain a sufficient approximation of the eigenvalue.
Another technique to find the lowest mode which is also used in our group is the Dimer Method.
Forcefield
The simulation is being tested with a simple forcefield TIP4P. But we plan in the future to use a more accurate forcefield called SCME. This forcefield gives water an electric field and a polarisabilty which match the those of DFT calculation up to the hexadecapole (4th degree).References
E.R. Batista and H. Jónsson, Diffusion and Island formation on the ice Ih basal plane surface, Computational Materials science, 20 (2001) 325-336.
G. Henkelman and H. Jónsson, A Dimer Method for Finding Saddle Points on High Dimensional Potential Surfaces Using Only First Derivatives, Journal of Chemical Physics, 111 (1999) 7010.
R.A. Olsen and G.J. Kroes, et al., Comparison of methods for finding saddle points without knoledge of final states, Journal of Chemical Physics, 121 (2004) 20, 9776.
F.D. Vila, E.R. Batista, et al., SCME: a Fast, Transferable H2O Interaction Potential Based on a Single Center Multipole Expansion, not published yet.
J.L.F. Abascal and C. Vega, A general purpose model for the condensed phases of water: TIP4P/2005, The Journal of Chemical Physisc, 123 (2005) 234505.
M.A. Novotny, A Tutorial on Advanced Dynamic Monte Carlo Methods for Systems with Discrete State Spaces, Annual Reviews of Computational Physics IX', editor D. Stauffer, (World Scientific, Singapore, 2001), pages 153-210 or http://front.math.ucdavis.edu/0109.0182.