Hello all,
I am trying to perform simulations of a triclinic molecular crystal TATB. I am using a non-polarizable, flexible molecule force field. Non-bonded interactions are modeled using the Buckingham potential and Coulomb interactions with fixed partial charges. The potential is parameterized so that all intra-molecular pairwise interactions are excluded. My LAMMPS version is 10 Feb 2015 and given below are the relevant lines of my input script:
boundary p p p
units real
pair_style buck/coul/long 11.0 11.0
neigh_modify exclude molecule all
kspace_style ewald/disp 1e-7
I have attached a plot of ecoul+elong as reported by LAMMPS thermo output (run 0) for different values of the cut-off (the largest value of cut-off is in accordance with the minimum image criterion) for a TATB crystal. One would expect ecoul+elong to be independent of the cut-off and this is evidently not so in the plot. My understanding is that, this results from using “neighbor_modify exclude molecule all” with the k-space solver (the relevant LAMMPS warning is present in the output; special_bonds is not adequate because there will be 1-5, 1-6, 1-7 etc… interactions). My questions are:
-
From the discussions in the mailing list, I understand that with the “neighbor_modify exclude” command, the interactions are properly excluded for the ‘ecoul’ calculation but not for the ‘elong’ calculation (interactions with images) and this is the source of the inconsistency. Is this correct? (I see that the elong doesn’t change irrespective of whether I use ‘neighbor_modify exclude molecule all’ or not).
-
Is the inconsistency mentioned above present in the force calculations also?
-
Is there a way to correct for the inconsistent elong calculation when all intra-molecular pairwise interactions are excluded so that an accurate (or at least approximate) calculation of the inter-molecular pairwise energy (ecoul+elong+evdwl) can be obtained from LAMMPS?
Thank you in advance for any help/suggestions.
ecoul-elong.pdf (3.42 KB)