Finding the Center of Mass of the Density Profile of Benzene

Hi LAMMPS community,

I am trying to determine the density profile of benzene in my LAMMPS simulation by calculating the center of mass. I’ve written some commands, but I’m not entirely sure if they are correctly configured for the task. Could someone review my code and provide any corrections or advice on how to optimize it? Here is what I’ve implemented:

  1. group benzene type 4 3
  2. compute benzen_den benzene chunk/atom bin/1d x lower 1.0
  3. compute benzen_com benzene com/chunk benzen_den
  4. fix dbenzen benzene ave/chunk 1000 1 1000 benzen_com density/mass density/number ave running file


Posting in a forum for having your code double-checked is a terrible idea. If you want to properly test a code, its better to apply it to simple situations for which you know the answer, such as a single molecule at a fixed position in the box.

This looks more like something that would be constructed by an AI than a logically thinking human being. Let’s review what you do:

  • You define a group of all atoms in benzene molecules by picking there atom types
  • You then define chunks for a 1d binning
  • You then compute the center of mass for all atoms in each bin
  • The last line makes no sense, since you are passing the per-chunk data where you should be passing the compute id for defining the chunks.

So this doesn’t match what you are asking for and is in itself not logical.

That said, your description asks about a density profile computed from the center of mass of the molecule. That makes little sense as well, since that will different from the actual density profile which would depend on the individual atoms. So you need to explain why you want to do that

Thank you for your feedback and insights. Could you explain why LAMMPS cannot calculate the density profile using the center of mass of the molecule? My research suggests that using individual atoms is preferred, but I’d like to understand the specific limitations or reasons behind this.


I happen to be in a lecturely mood now …

As documented, LAMMPS maintains internal lists of the co-ordinates and masses of each atom. Therefore LAMMPS can easily assign each atom to a chunk (using compute chunk/atom), and then average (or sum) the masses of all atoms in a chunk (using fix ave/chunk).

By extension, fix ave/chunk can also aggregate other atomic properties tracked by LAMMPS, as documented.

LAMMPS does not automatically track “the centre of mass of the molecule this atom is a member of” for each atom. Why would it? The molecular CoM does not generally affect each atom’s trajectory (other than indirectly through the positions of other atoms – which LAMMPS already tracks and directly calculates with), and thus LAMMPS ignores it for efficiency.

You can, nonetheless, “tag” each atom with its molecular CoM (again, as documented), by first using a compute chunk/atom to assign molecules as chunks, then calculating the CoM with compute com/chunk and “spreading” the information back to each atom with compute chunk/spread/atom.

But should you? Molecules exist primarily because Nature is too lazy to break chemical bonds. Any experiment probing the structure of benzene, as far as I know, will either be too coarse to “see” individual molecules and their centers-of-mass, or fine enough to “see” individual atoms (certainly the carbons at least), making the molecular CoMs physically redundant.

In any case – you now know how to track molecular CoMs, as documented. All the best!

Before asking a followup question, wouldn’t it be the polite thing to first address what I was asking you to do (explain why you want to do this)?

The answer to your followup question is simple: it hasn’t been programmed. The majority of features in LAMMPS are available because somebody had a need for it, programmed it and then contributed the code to the distribution. The fact that some feature doesn’t exist in turn means, that nobody has (yet) felt the need for it.

I appreciate your lecturely mood and the effort you put into sharing this information.

I apologize for not addressing your question earlier. The reason I want to calculate the density profile based on the molecular center of mass is to gain a more accurate understanding of the spatial distribution of entire molecules, rather than just individual atoms.

I fail to see how this can be more accurate. I would say the opposite is the case, and you have more noisy data since you are no longer taking into account the orientation of the molecules.