# [lammps-users] computing electric field

Dear Lammps-Users,

I am interested in computing the value of the electric field on any given atom in my simulation. I am looking at a simple system of 2 spc/e water molecules located near the center of a cubic box with length = 200 Angstroms.

I compare the following two numbers:
(1) Force and Efield computed from the definition of force and efield due to point charges. (Thus a sum over all charges in the system divided by either r^2 or r, respectively).
(2) Force and Efield computed from the force calculation based on the pair potential. I keep track of the forces due to Lennard-Jones and SHAKE and subtract those from the total force to give me (what I presume to be) the Coulombic part of the force: forceCOUL = forceTOTAL - forceLJ - forceSHAKE
(NOTE: special_bonds is set to special_bonds 1 1 1 1 1 1.)

When I compare values from (1) and (2) using pair lj/cut/coul/cut with cutoff 20, I get the same values. Number (2) is in fact not sensitive to the cutoff length as long as it's larger than about 3.

When I compare the values from (1) and (2) using pair lj/cut/coul/long with cutoff 20, the numbers differ. I assumed this was because of the long-ranged forces, but when I keep track of the long-ranged forces (in ewald.cpp) and get forceCoulomb from forceCoul = forceTOTAL - forceLJ - forceSHAKE - forceLONG, the results still do not match up.

My questions are
1. Is there an additional contribution to the force in lj/cut/coul/long that I am neglecting? What is it, and where is it coming from?

2. Units. Forces in LAMMPS are in kcal/mol-Angstrom. I am interested in the electric field in units Volts/Angstrom = J/C-Angstrom. Since charge is basically a number in LAMMPS (multiples of electron charge), I need a multiplicative factor alpha so that the units work out

Efield = alpha * forceCoul / q
alpha = 0.24 / e Na = 2.49e-6

where e = electron charge in Coulombs, Na = Avogadro's number, q is the charge of the atom of interest. This gives me values that are about 1.0e4 off in magnitude from what I expect.
I've gone through the unit conversions and I don't see my mistake. Any help here would be much appreciated as well!

Thanks in advance for any help!

Joyce

My questions are
1. Is there an additional contribution to the force in
lj/cut/coul/long that I am neglecting? What is it, and where is
it coming from?

lj/cut/coul/long does include a damping term (erfc()) on the
Coulombic pair interactions, to match up with PPPM long-range.

2. Units. Forces in LAMMPS are in kcal/mol-Angstrom. I am
interested in the electric field in units Volts/Angstrom =
J/C-Angstrom. Since charge is basically a number in LAMMPS
(multiples of electron charge), I need a multiplicative factor
alpha so that the units work out

Efield = alpha * forceCoul / q
alpha = 0.24 / e Na = 2.49e-6

where e = electron charge in Coulombs, Na = Avogadro's number, q
is the charge of the atom of interest. This gives me values that
are about 1.0e4 off in magnitude from what I expect.
I've gone through the unit conversions and I don't see my
mistake. Any help here would be much appreciated as well!

The conversion from qq/r to energy units used the factor qqr2e in
update.cpp, which is 332 for real units. If you can figure out
that one, you should be able to get yours.

Steve