Non-zero energy above cutoff for coul/dsf potential

LAMMPS LAMMPS General Discussion
1

Hi,
I constructed a simple test case with two charged atoms in a large box to check the working principle of the potential (pair_style) coul/dsf, which uses the complimentary error function (erfc) to damp the coulomb potential. However, the potential energy of the system is not zero, when the distance between the atoms is larger than the cutoff radius. My input file looks as follows:

###########################################
units metal
atom_style charge
pair_style coul/dsf 0.02 3.5
read_data structure.inp
pair_coeff * *
thermo 1
thermo_style custom step pe ecoul
run 0
#######################################

with the structure.inp:

test system
2 atoms
2 atom types

0.0 40.000000 xlo xhi
0.0 40.000000 ylo yhi
0.0 40.000000 zlo zhi
0.000000 0.000000 0.000000 xy xz yz

Masses

1 91.224000
2 15.999400

Atoms

1 1 1.000000 0.000000 0.000000 0.000000
2 2 -1.000000 0.000000 0.000000 4.000000

The distance between the atoms is 4 Angstroem. Cutoff is 3.5 A. But LAMMPS returns a potential energy of -8.2278401 eV. How can that be?

Also, I could not reproduce the potential energy returned by LAMMPS for systems where the atoms are closer together than r_cutoff when I insert the values into the formula given in the documentation of coul/dsf (the formula given there is correct, I checked it with the original paper).

Am I understanding something wrong or is the problem in the implementation of the potential?

3

Since this is - like the Wolf summation - an approximation to the Ewald sum, the pair style also computes a “self-energy” and that term is independent of the cutoff. See lines 97-100 of the pair_coul_dsf.cpp file.

4

Thanks for your quick reply.
I made some more tests: keeping the distance between the atoms constant at 4 A and changing r_cutoff to lower values (e.g. 3.0 A, 2.5 A, …) leads to a change of the potential energy. That means that there is an energy shift deltaE depending on r_cutoff (and also alpha). For a fixed r_cutoff, deltaE stays constant when r is increased (r > r_cutoff). If r is decreased to values r < r_cutoff, the potential energy computed by LAMMPS is identical to the sum of the energy shift deltaE and the value from the formula of the potential given in the manual.