Dear LAMMPS developers,
For pair styles which require a kspace solver, such as lj/cut/coul/long, the special neighbors 1-2, 1-3, and 1-4 are included in the neighbor list. The short-range term q[i]*q[j]erfc(alphar)/r is computed for all pairs with r <= rc, including these special neighbors. However, whenever factor_coul < 1, a calculation is done so that the correct interaction is recovered afterwards when kspace includes the long-range term. Thus,
prefactor short-range special neighbor kspace
------^------ ------^------ ---------^------- ------^-----
(q[i]q[j]/r)[ erfc(alphar) - (1 - factor_coul) + erf(alphar) ]
This procedure works because erfc(x) + erf(x) = 1. For instance, interaction will fade out if factor_coul = 0.
What I do not understand is why the same procedure is done for lj/cut/coul/dsf and related pair styles, which do not involve a kspace solver and the term prefactorerf(alphar) will not be added afterwards. In the compute() routine of lj/cut/coul/dsf, the line if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*prefactor serves to exclude those interactions which should not have been considered initially. However, they are subtracted without performing any shifting/damping. I’m afraid there is an inconsistency there. Am I right?.
Moreover, I wonder if it is correct to apply any shifting/damping to the special-bond coulombic interactions. To the best of my knowledge, none of the force fields commonly employed in molecular dynamics includes shifting or damping in their parameterization procedures. The 1-4 coulombic interaction is an important part of a dihedral model, so that the fitted parameter values often reflect the presence of such interaction.
Thank you,
Ana Silveira
Dear LAMMPS developers,
For pair styles which require a kspace solver, such as lj/cut/coul/long,
the special neighbors 1-2, 1-3, and 1-4 are included in the neighbor list.
The short-range term q[i]*q[j]*erfc(alpha*r)/r is computed for all pairs
with r <= rc, including these special neighbors. However, whenever
factor_coul < 1, a calculation is done so that the correct interaction is
recovered afterwards when kspace includes the long-range term. Thus,
prefactor short-range special neighbor kspace
------^------ ------^------ ---------^-------
------^-----
(q[i]*q[j]/r)*[ erfc(alpha*r) - (1 - factor_coul) + erf(alpha*r) ]
This procedure works because erfc(x) + erf(x) = 1. For instance,
interaction will fade out if factor_coul = 0.
What I do not understand is why the same procedure is done for
lj/cut/coul/dsf and related pair styles, which do not involve a kspace
solver and the term prefactor*erf(alpha*r) will not be added afterwards.
In the compute() routine of lj/cut/coul/dsf, the line if (factor_coul <
1.0) forcecoul -= (1.0-factor_coul)*prefactor serves to exclude those
interactions which should not have been considered initially. However, they
are subtracted without performing any shifting/damping. I’m afraid there is
an inconsistency there. Am I right?.
no. the kspace term is implicit in coul/dsf or coul/wolf. the purpose of
this special treatment for special neighbors is that in those cases you
must exclude the *full* coulomb interaction to be consistent since you
cannot exclude them correctly in the kspace calculation, regardless of
whether you do that one explicitly or implicitly.
Moreover, I wonder if it is correct to apply any shifting/damping to the
special-bond coulombic interactions. To the best of my knowledge, none of
the force fields commonly employed in molecular dynamics includes shifting
or damping in their parameterization procedures. The 1-4 coulombic
interaction is an important part of a dihedral model, so that the fitted
parameter values often reflect the presence of such interaction.
of course not. that is why there is the special treatment in those pair
styles. with the code how it is currently in LAMMPS, you always exclude (or
scale) the same interaction strength, i.e. the undamped coulomb.
axel.