Replicate causes differences in per-atom electrostatics

Hi all,

In the attached input and data file (derived from Data for: "Assessment of GAFF and OPLS force fields for urea: crystal and aqueous solution properties" - University of Strathclyde), I am simply doing a run 0 before and after a replicate. The per-atom electrostatic energy and the max force are quite different before and after the replicate 5 5 5 command (themo_modify norm yes has been issued):

     Fmax          E_mol          E_pair         E_vdwl         E_tail         E_coul         E_long    
 1.18487        1.4489072     -25.078458     -0.80193461    -0.059558078    1.1583971     -25.43492     
     Fmax          E_mol          E_pair         E_vdwl         E_tail         E_coul         E_long    
 6.7109147      1.4489072     -25.582828     -0.80662136    -0.059558078   -1.0421779     -23.734029    

What gives? Is this expected? I tried several things to overcome this (tightening kspace_style tolerance, kspace_modify force), and nothing made a difference.
data.minCELL (2.3 KB)
test_replicate.lmp (1.7 KB)

Did you see this warning?

WARNING: Bond/angle/dihedral extent > half of periodic box length (src/domain.cpp:1174)

This is a problem because then the identification of the nearest neighbors and exclusions can go wrong (which appears to be the case here).

This warning goes away when you do replicate 2 2 2 and you have fairly consistent results afterwards with further replications.

Those force/energy values will be even more consistent when you do the following 3 steps that reduce the degree of approximation and thus noise in the calculation.

  • tighten the kspace convergence to 1.0e-6
  • disable coulomb tabulation with pair_modify table 0
  • check on the FFT grid and use kspace_modify mesh to adjust it so it is exactly replicated from the smaller systems. LAMMPS will use the smallest possible grid that can be factorized into small prime numbers for the FFT part of PPPM. For larger systems, typically there are more possible factorizations and LAMMPS will use the smallest grid that is sufficient for the requested accuracy.
1 Like

I did see the warning, but I did not think anything of it because the error was not in the bonded terms, but I guess this because the accounting in special_bonds affects nonbonded terms?

You are correct that if I start with 2x2x2, further replication makes a much smaller difference. Thank you, this was very frustrating!

While LAMMPS technically does not require minimum image conventions (with a few notable exceptions), there are some ambiguities that arise when you have bonded interactions that extend beyond half the box. Then molecules interact with their own periodic replicas. This should be avoided in general.

2 Likes