Dear Lammps users,
I am interested in computing the electric field on a given atom in my system of spc/e water molecules. To do this, I introduced a new member to the Atom class to keep track of the forces due to the Lennard-Jones part of the potential. Then in a fix (executed at end_of_step), I use the relation
F_coul = F_total - F_LJ
to obtain the total electrostatics force on a given atom, and thus the electric field on that atom, since E_i = F_coul(at i)/q_i, where i indicates my atom of interest.
When I compare my results for this electric field with other published results, they do not match up. I've been through the code and as far as I could tell the above procedure is a reasonable way to compute F_coul and thus the E-field. I also compared my E-field value distributions for various cutoff values and different system sizes and I get the same results.
My guess is that I am misunderstanding some part of how the LAMMPS code works in dealing with the forces the way that I do. I would appreciate anyone being able to point out what's wrong with the above procedure.
much thanks,
Joyce