Dear Lammps users,
I read in the manual (https://docs.lammps.org/compute_pe_atom.html) that
Note that the energy of each atom is due to its interaction with all other atoms in the simulation, not just with other atoms in the group.
For my problem focusing on twisted bilayer graphene (where my groups are defined as atoms in either the top or bottom layer), I wonder if a work-around to this limitation would be to relax the bilayer system, and then, remove one of the layers by hand, to then calculate the potential and/or kinetic energy for the atoms in the relaxed, possibly corrugated left-over single layer graphene layer, without doing any further relaxation of the system.
But before heading in that direction, I wonder if there is a ready-to-use solution/trick to only include the interaction with atoms that are in the same group?
Regards,
Nicolas
Dear Nicolas,
If your interlayer potential is pairwise and there are no bonds between layers, you can use “compute pe/tally” to calculate interlayer energy per atom and then subtract it from energy per atom.
Another approach is to use “neigh_modify exclude”, to remove interactions between layers.
Cheers,
Michał
What can be done depends on the kind of pair style and/or force field you are using and whether you need the sum of the potential energy or the per-atom info.
If you need per-atom info, the best approach would be to first run you simulation without computing pe/atom data and just record a suitable trajectory.
Then you can use rerun and then only define the interactions you want to include.
Otherwise you could use compute pe/mol/tally, or compute group/group or give your different groups different atom types and then use multiple substyles with pair style hybrid and extract the energy via compute pair. All of these options have limitations and thus may not be suitable or may require some post-processing or creativity.
The approach using rerun is the most flexible and striaghtforward, but also more time consuming.
Axel.
Thank you Axel and Michal for useful answers.
The pe/tally and group/group options don’t seem suitable in my case as I am using the (many-body) rebo potential for intralayer and KC type drip potential for interlayer. Both corresponding manual pages seem to heed against using it with such many-body potential. Also got a warning when testing pe/tally.
So instead of using pe/tally, I am just using pe/atom now, for the single layer (intralayer) and the double layer (total energy) system, after which I subtract their respective contributions to get the corresponding interlayer energy per atom. My first test gave zero interlayer energy, so will confirm my steps.
The rerun seems quite useful indeed for ke/atom info.
Regards,
Nicolas