Molecular rotational correlation in LAMMPS

Currently I am calculating rotational correlation by dumping coordinates to text file followed by a simple python script to calculate molecule wise vectors and their dot product. But results are quite “noisy” as compared to some in our lab who did is long time ago by using Material studio.

Is there any way i can do it in lammps at run time itself? i.e.:

  1. define a vector for a molecule, say water, from midpoint of 2 hydrogen towards oxygen or from one hydrogen atom to another.
  2. calculate its dot product with t0 value and average all such products in simulation

if i read correctly fix ave/correlate can do latter part but I am not sure on defining molecule wise vectors

Interesting Q. Can’t think of an easy current way to

do this in LAMMPS. You could write a compute to

tally the vector for each molecule (chunk in LAMMPS-speak).

But that’s not trivial. One atom doesn’t “know” which other

atoms are in the molecule without searching for them.

How do you codify “vector from midpoint of 2 H to O” in a general

way for any molecule of any size. If the size is big, the proc

won’t even know all the atoms in the molecule.

There is a compute inertia/chunk which will give you info

on the shape of each molecule, which are eigenvalues. I suppose it could

be amended to compute the eigen vector(s) which for water

one of them would be close to the vector you are asking for.

Once you have a vector per-molecule other LAMMPS commands

can output or average them. You could also probably convert

the vector to an orientation angle and histogram them by angle

or spatially bin them, etc.


Although I dont think following is true always, but in all my input files (created from topo tools) usually there is a pattern in molecules. Like for example, suppose if oxygen is of type 4 and hydrogen is of type 1, input data file will have data as first for both hydrogens then oxygen then both hydrogens etc, so if order 1 1 4 1 1 4 1 1 4 and so on. Hence I always provide input to analysis scripts for arbitrary molecule as all types of atoms in molecule and total number of atoms. This creates a list of connected atoms from which i can pick and chose any atom etc i want. Like currently am defining is as " in every molecule take mid point of 2 and 3rd atom and subtract it from 1st atom" this way i can define any vector on the basis of atomic coordinates, COM etc.
I am brushing up my c++ (always used 98) after which I would be more than happy to try writing it for lammps (was looking for some thing to write in lammps just to get hang of it), but given slightly tight schedule, wont count on me! (only thing i ever tried in lammps was following a blog ( verbatim!)