I am doing MD simulations of triclinic molecular crystal TATB, in which the
non-bonded interactions are modeled using the Buckingham potential and
Coloumb interactions with fixed partial charges. Bonded interactions are
modeled using bonds, angles, dihedrals and impropers. All intra-molecular
non-bonded interactions are excluded. Long range Coloumbic interactions are
computed using the PPPM solver. The parameters are taken from J. Chem. Phys.
139, 074503 (2013). My LAMMPS version is Dec. 7, 2013.
There are 6 atom types in my molecule. I need to introduce an additional
pairwise, repulsive interaction only between the intra-molecular 1-2 pairs.
This repulsive interaction should be excluded for all inter-molecular 1-2
pairs. This is what I think I should do:
1. Tabulate the r^12 repulsion with non-zero values only for the 1-2
2. Use pair_style hybrid/overlay buck/coul/long table
3. use neigh_modify exclude molecule group36, where group36 is composed of
all atoms of type 3-6
4. modify pair_table.cpp with an 'if loop' to calculate the repulsion only
for the intra-molecular 1-2 pairs.
I'm not sure what an "if loop" is in this context. Can you clarify step 4?
I can think of two other ways:
If you know exactly which pairs of atoms of type 1 and 2 are likely to
be close contact, or if there aren't too many of them, you can try
connecting all of these pairs of atomswith explicit bonds and use
bond_style table (or hybrid harmonic table, etc..) to add specific
1/r^12 interactions between these atom pairs.
Alternately, there is a pair_style which allows you to use different
Lennard-Jones coefficients depending on whether the two atoms belong
to the same molecule. Check out:
You can download the code for that pair style here:
At least 3 different groups have used it successfully so far.
It was tested in late 2013. (Let me know if LAMMPS has changed and
this code no longer compiles.)
The only problem is that this pair_style assumes you want to use
short-range coulombics/electrostatics. If you want to use long-range
coulombics, you would have to turn off the short range electrostatics
forces, and then use
"pair_style hybrid/overlay" to combine these Lennard Jones forces with
There's a (long and confusing) discussion of how to combine this style
with long-range coulombics here:
Here's an excerpt from that post:
pair_style hybrid/overlay &
lj/charmm/coul/charmm/inter 29.5 30 0 0 &
pair_coeff 1 1 lj/charmm/coul/charmm/inter 0 0 0 0 0.0298 3.4
pair_coeff 1 2 lj/charmm/coul/charmm/inter 0 0 0 0 0.0298 3.5
pair_coeff 2 2 lj/charmm/coul/charmm/inter 0 0 0 0 0.0298 3.6
pair_coeff * * coul/long