compute rotational energy of flexible water molecules

Hi all,

I am doing simulations of water molecules with the flexible TIP3P model. I am interested in calculating the total translational and rotational kinetic energies of the system.

With that purpose I have defined a chunk for every water molecule, and with the commands "omega/chunk" and "inertia/chunk", angular velocity and inertia tensor can be computed for every molecule.

So my question is: How can I compute translational and rotational kinetic energies for every molecule and sum them to obtain total quantities? Should I use a compute or a variable to achieve this?

Thanks in advance,


Compute reduce/chunk could be used with compute ke/atom as
input to sum the ke per chunk. But that’s probably not what
you are looking for. That is just total KE of the system. I presume you want the partitioning of
that total KE into a rotational term and a translation of the COM
term. I think 0.5 I omega^2 is the rotational energy you want, you
could post-process is from compute omega/chunk and inertia/chunk
values, and subtract it from total KE to get translational KE of the COMs
I think. Or 2 new computes could be added: compute ke/chunk (for COMs)
and compute erotate/chunk.


Thanks a lot for your answer.

I have many molecules in my simulations, so I would like to avoid post-processing because it would be a very time-demanding task and the dynamics would be much slower.

Do you think is possible to define a variable to compute I omega^2 and sum it over chunks (molecules) within the existing commands, or the only way to obtain the total energies would be adding new computes to the code?


I think the variable formula would be long and messy since I is a tensor (3x3 matrix),
so its probably something like w_transpose I_matrix w_vector. You could try it.

Not clear on why post-processing would slow down the dynamics, unless you
think outputting what is needed would be prohibitive. I suggest doing a few
snapshots that way, to see if its less painful than you think. Also you certainly
want to know the right answer to compare to, if you start writing a compilcated
variable or C++ code (new compute) to do it.