[gmx-developers] Scaling velocities in Velocity Verlet
David van der Spoel
spoel at xray.bmc.uu.se
Sat Dec 3 21:22:33 CET 2011
I'm implementing an algorithm where I change the potential energy Epot
on time t by an amount DE. In order to conserve total energy I would
like to change the kinetic energy Ekin by an amount -DE. This is done by
scaling the velocities at time t by lambda = sqrt(1-DE/Ekin).
The problem is that the change in Ekin comes one step too late. I tried
to deduce the algorithm from the code (4.5 code base) and deduced the
following MD algorithm:
v(t) = v(t-dt/2) + a(t) dt/2
v'(t) = lambda v(t) <======= my addition
v(t+dt/2) = v'(t) + a(t) dt/2
r(t+dt) = r(t) + v(t+dt/2) dt
However, the Ekin that is printed is the uncorrected one even though I
scale both ekind->ekin and ekind->tcstat[x].ekinf .
Is the algorithm correct as I write it above? Is the Ekin(t) already
stored at the beginning of the time step?
David van der Spoel, Ph.D., Professor of Biology
Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone: +46184714205.
spoel at xray.bmc.uu.se http://folding.bmc.uu.se
More information about the gromacs.org_gmx-developers