[gmx-users] RE: Difference between kinetic energy from .xvg file and calculated kinetic energy from velocities

[Mark Abraham]

On 28/03/2010 12:25 AM, Samantha Sanders wrote:
>
> On 27/03/2010 1:57 PM, Samantha Sanders wrote:
>> Hello!
>>
>> I have run a simulation and using g_energy on the .edr file, I have been
>> looking at the kinetic energies in the resulting .xvg file. It is my
>> understanding that these are the instantaneous kinetic energies for
>> those times. However, when I try to calculate the kinetic energy on my
>> own using the frames from my trajectory file and being sure to compare
>> the same timestamps of the trajectory and the .edr file, I am not
>> getting the same values for the kinetic energy. My first question is
>> should I get the same values between the two sources: the .xvg file and
>> the velocities in the frames of my trajectory file? If so, any
>> suggestions for why I may not be getting the same values or if there are
>> any other adjustments I should make? I have simply been calculating the
>> kinetic energy as the sum of .5mv^2 for all particles. I have read the
>> manual and the various mailing list responses and I am still unsure
>> about this. There may be something simple I am missing. Any help is
>> appreciated. Thank you!
>
> See manual 3.4.4. The leap-frog integration algorithm never computes the
> velocities and positions at the same time point. Thus, your
> after-the-fact computation cannot agree with those reported by GROMACS.
>
> Mark
>
> Thank you for your reply.  I understand that the velocities and positions
> are not at the same time point.  I am just using the velocities reported
> from the trajectory to calculate the kinetic energy.  I am not using the
> positions at all.

Ah, quite right. Sorry, I confused some ideas in my head.

>  However, I did look at that section in the manual.  I
> think I understand what it is saying now, but let me restate it. Since the
> kinetic energy is calculated before velocities are updated in the second and
> third steps respectively, when the energy and velocities are output in the
> fourth step, the velocities are now different from the ones used to
> calculate the outputted energy for that step.  Is this correct?

If you're referring to Figure 3.3, that looks right.

> 	Also, if I understand the algorithm correctly when positions and
> velocities are output in the same step, the outputted positions are at t +
> deltat and the outputted velocities are at t + (deltat/2).  Correct?

That seems about right, but I won't speak definitively. If you use 
nstvout = 1 and nstenergy = 1, you should be able to compute the kinetic 
energy with one of two possible sets of velocities. Positions and 
energies will always match in t. If you really want a detailed 
understanding, set up a toy Lennard-Jones simulation in a spreadsheet 
that you can compute alongside GROMACS for a step or two.

Mark