Dear Axel Kohlmeyer
Thank you for your fruitful reply to my questions and also about your advices. I am indebted to you. I think our discussion about the external electric field in the reactive MD simulations can help to clear the subject and also find good solutions to the probably problems.
About the energy conservation along with the external electric field (or force) I think when we are imposing an external electric field (or force) continually to the simulation box (if we have charged particles in the system, this is for the efield), actually we are injecting some energy continually to the system and it is clear that the total energy of the system is not conserved during the simulation. In the other words, the system is not isolated and it is closed. I questioned this subject (in my last email) because of this reasons:

We had written an article (entitled: a simulation with ReaxFF by applying an external efield), but the referee made a major revise due to two reasons that one of them was about the problem with conservation energy along with the ReaxFF simulation when we use external efield. He/she told “*The most important thing in reactive molecular dynamics simulation is that the total energy of the system is conserved. If the ReaxFF program uses the Verlet method for the integrator, it is impossible to preserve the total energy of the system, and the calculation results shown in this paper cannot be verified.*The author need to use the sixth symplectic integrator with a time step of 0.1 fs. the symplectic integrator is the numerical integraton scheme for Hamiltonian system, which conserve the symplectic twoform exactly, so that (q(0), p(0)) → (q(τ), p(τ)) is canonical transformation. This algorithm is accurate and has no accumulation of numerical errors for total energy in constrast to the other common algorithm to solve the Hamiltonial equation of motion. ” We know that ReaxFF uses the velocity Verlet, not the ‘normal’ Verlet integrator.

In reactive MD simulations by applying continually constant electric field, when charged particles are crossing the PBC I seen a jump in kinetic energy. I think this is because the opposite sign of the efield at the in front wall in direction to the efield. I did this with ReaxFF Fortran source code.

In a simple system I checked the differential total energy in alternatively similar delta t during the simulation with continually constant efield but I didn’t gain the equal total energy differences. I think this is due to continually charge equilibration (and changing the atomic partial charges) in each time step and also change the position of the particles in the box. So, we have different forces and energies due to the external electric field and also different delta total energy in each delta t. I did this with ReaxFF Fortran source code.
It is not clear for me that exactly/truly what happens for the energy when we use external electric field in the Reactive MD simulation! If you have any suggestion or criticism to the causes and effects raised above, I gladly welcomed.
Thanks in advance for your help and feedback.
Best regard
Mohammad Ebrahim izadi,
Department of Chemistry,
Tehran University,
Islamic Republic of Iran,
Phone : +98 – 21 – 61113358
Fax : +98 – 21 – 66409348