# [gmx-users] non isotropic kinetic energy

David van der Spoel spoel at xray.bmc.uu.se
Tue Sep 15 11:27:54 CEST 2009

```herz wrote:
> Hi,
>
> I have done exactly what you suggest in your first paragraph, seperating
> translation and rotation for each water molecule.
> If you want I can send you the new ekin_vs_t_and_erot_vs_t.dat where you
> can see that the ekin is actually equidistributed
> but the erot isn't!
>
> So in this rigid body perspective each water molecule has 6 dof (3
> translat 3 rot) and the 3 rot dof are not
> equipartitioned! Now the rigid body is pretty much the text book example
> for the equipartition theorem.
>
> So you are right in so far as that the first measurements of ekin we did
> while ago only are insufficient for the rigid body.
> Now we measure ekin and erot seperatly for all 3 components per rigid
> body and they are still not equipartitioned!

Have you tested the equipartitioning per atom? This must be correct.

To make water move as a rigid body one would probably need to determine
the force and torque on the center of mass and use that to integrate the
equations of motion. A quaternion description of rigid water would do
this I presume. Did you analyze the EkRot around the center of mass of
the molecule?
Whether or not a quaternion is equivalent to a constrained atomic system
I do not know, but Berk's answer below seems to indicate it is not.
>
> Alex
>
> Berk Hess schrieb:
>> If you would consider a water molecule a generic rigid body,
>> you would have to calculate the kinetic energy in a different way,
>> you would have to separate the rotational and translational kinetic energy.
>>
>> I'll try to explain it with an example:
>> say you have a water molecule with all three atoms in the x-y plane.
>> Non-constrained there will be 3*2=6 degrees of freedom in the x-y plane.
>> If you constrain the two bonds and the angle you remove 3 degrees
>> of freedom in the x-y plane. You are left with 3: x and y translation
>> of the center of mass and rotation in the x-y plane.
>> Out of this plane (the z direction) you still have 3 degrees of freedom.
>> So in z you have the full Ekin, whereas in x and y you each have half
>> the degrees of freedom and half the Ekin of the unconstrained water.
>>
>> Berk
>>
>>
>>> Date: Mon, 14 Sep 2009 19:58:20 +0200
>>> From: alexander.herz at mytum.de
>>> To: gmx-users at gromacs.org
>>> Subject: Re: [gmx-users] non isotropic kinetic energy
>>>
>>> Hey,
>>>
>>> rigid bodies.
>>> The rigidity constraint is implemented via constraints (settle or so)
>>> but could also (in theory) be implemented
>>> using some euler angles or what ever scheme you prefer. It should not
>>> matter how you implement
>>> the rigidity, the result of a rigid body in some potential should be
>>> identical (if implemented correctly).
>>>
>>> The equipartition theorem must hold for a rigid body (it has 3
>>> translational and 3 (for water) rotational degrees
>>> of freedom which must contribute 1/2*k_b*T per molecule, irrespective
>>> how you implement rigidity.
>>> Also I wouldn't see how you could explain the energy difference if the
>>> molecules were actually implemented
>>> as rigid bodies with e.g. euler angles. Why should the euler angles care
>>> about the water orientation at the interface boundary?
>>>
>>> I ran some more analysis assuming the water molecules are rigid bodies.
>>> This shows that the translational energy of the molecules is actually
>>> correct (so 1/2k T) but the average rotational energy is incorrect (with
>>> settle the z part is too high and x/y ok,with lincs x/y is too low and z
>>> too high).
>>>
>>> Summarizing, I fail to see why choosing one method over the other for
>>> performing rigid body dynamics should have an influence over
>>> equipartition of energy in the system. I'm still not convinced the
>>> ensemble is correct.
>>>
>>> Alex
>>>
>>>
>>> Berk Hess schrieb:
>>>
>>>> Equipartitioning still holds.
>>>> The issue here is that constraints remove degrees of freedom.
>>>> If you remove a degree of freedom, there is not potential or kinetic
>>>> energy
>>>> anymore in that generalized coordinate.
>>>> This does not introduce any problems with ensembles or something like
>>>> that.
>>>>
>>>> I did not have time yet to implement pull code with pbc.
>>>> I will do it this week.
>>>>
>>>> Berk
>>>>
>>>>
>>>>> Date: Mon, 14 Sep 2009 16:08:55 +0200
>>>>> From: alexander.herz at mytum.de
>>>>> To: gmx-users at gromacs.org
>>>>> Subject: Re: [gmx-users] non isotropic kinetic energy
>>>>>
>>>>> Hm,
>>>>>
>>>>> apparently there are constraint algorithms which don't change the
>>>>>
>>>> ensemble
>>>>
>>>>> (e.g.
>>>>>
>>>> http://www.informaworld.com/smpp/content~db=all~content=a750934359)
>>>>
>>>> ------------------------------------------------------------------------
>>>>
>>>> _______________________________________________
>>>> gmx-users mailing list    gmx-users at gromacs.org
>>>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>>>> Please search the archive at http://www.gromacs.org/search before posting!
>>>> Please don't post (un)subscribe requests to the list. Use the
>>>> www interface or send it to gmx-users-request at gromacs.org.
>
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/search before posting!
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.