# [gmx-users] non isotropic kinetic energy

Berk Hess gmx3 at hotmail.com
Tue Sep 15 09:15:56 CEST 2009

```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,
>
> I thought a little more about this. Eventually, the water molecules are
> 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)
```