Hi Lammps User.
I am trying to compute, during my simulation, energy differences betwen steps where the only difference is that I’ve added a particle to the simulation box, all the other particles remain in their original coordinates. In this situation the difference in energy between the system with N particles and the one with N+1 particle is only the interaction energy between the particle that I’m adding and the rest of the particles. As Lammps has a compute pe / atom command, I compared the two alternatives, computing energy differences of the systems and computing just the energy of the new particle. Unfortunately I’ve found differences.
In order to find what was going on I simulated a simple system of two water separated by their hydrogen bond distance. The model of water used was a coarse grain model based on a modify Stillinger-Weber potential
these are the results
PotEng (default thermo output) -5.4790852
compute all pe -5.4790852
compute all pe/atom (water 1) -2.7395426 !!!
compute all pe/atom (water 2) -2.7395426 !!!
sum of water 1 and water 2 -5.4790852
if I look the PotEng, I got a quite nice hydrogen bond energy around 5 Kcal/mol, but here is what I don’t understand, As I get from the manual the pe/atom is the energy of each particular atom interacting with the rest of the atoms in the simulation so as there are only two atoms here the energy of the water1-water2 should be the hydrogen bond energy -5.47 ,the same for the water2-water1 and that should be the total PotEng. So why does the pe/atom energy appears divided by two ?
I recall, I’m using a Stillinger-Weber with parameters for water molecules, these results where obtained by doing a minimization.