Compute hybrid substyle energies for different groups


I wanted to ask if it is possible to extract energies of different substyles used in a hybrid potential for different groups of atoms. Specifically, if I have a system which is a hybrid of lj/cut and born/coul/dsf, and I have two groups of atoms A and B, is it possible to obtain the energy of the A-B interactions, split by the LJ component and the BMH component? As far as I can tell there is no out of the box way to do this. These are the options I have considered:

  1. Use compute group/group for the A-B energies - however this would return the net interaction between the two groups, and cannot be split by individual components of the hybrid potential energy function.
  2. Use compute pair, where the A-B interaction can be computed as U_AB = U_all - (U_AA + U_BB). This command allows extraction of energies of individual substyles (e.g. compute 1 all pair lj/cut), but not by group, since (as mentioned in the LAMMPS documentation) the group specified for this command is ignored.
  3. Use compute pe/atom for the three groups all, A and B. However here the energy components can only be split as bonded or (net) nonbonded interactions, and not by individual components of a hybrid pair_style.

It seems to me that the most appropriate way to do this would be the 2nd option, but where the group command in compute pair actually is taken into account. I wanted to ask if anyone could provide me guidance on what specifically to edit in the compute_pair .cpp and .h files to enable this (or if there are any other suggestions for how to obtain these energies).


There is a fourth option: record a trajectory without any specific subsystem output and then use the rerun command on the same trajectory multiple times, where only the interactions of interest for each subset are enabled for the individual runs.

Hi Axel,

Thanks for your reply. I have never used the rerun command before so I want to confirm if I understand its usage correctly. If I already have a dump file with a trajectory, I can then use a separate input file with the rerun command instead of a run command, which will read in positions from that dump file instead of actually doing time-integration. Is that correct? So then in this rerun input file I can still define groups, perform computations, get thermodynamic output in a log file etc like in a regular simulation?