Problem with mapping atoms in custom LJ potential

Dear LAMMPS users,

I am trying to implement a version of LJ potential for coarse grained simulation which doesn’t have any corrugation artifact produced as a result of coarse graining a line into atoms. The basic algorithm is similar to pair_style line/lj but is implemented for atoms instead of lines.

In order to calculate lj potential between atoms i and j, I need orientation of each atom with respect to the adjacent atoms. It can be easily extracted if I can access the coordinates of adjacent atoms (i-1, i+1 & j-1, j+1).

I am implementing it as follows:

  1. get global id of local atom i by tag[i].
  2. calculate global ids of adjacent atoms as tag[i]+1 and tag[i]-1
  3. convert global ids to local ids, m and n, as atom->map(tag[i]+1) & atom->map(tag[i]-1).
  4. obtain images of m and n closest to i as domain->closest_image(i,m) & domain->closest_image(i,n).
  5. repeat same process for atom j.

The code works fine in serial version. But it crashes after running for some time step in parallel implementation with exit code 11. I tried changing several parameters. The only thing that worked was changing neighbour skin distance from 2 to 20 A. But it increases the computation cost significantly.

It will be very helpful if you can kindly let me know what might cause such error. I am using LAMMPS version dated 17Nov16. My system size is roughly 400,000 atoms. The cutoff for lj potential is 28 A and the style of atoms is angle with 10 A and 180 degrees as equilibrium bond length and angle.

Thanks & Regards,
Ankit Gupta

Dear LAMMPS users,

I am trying to implement a version of LJ potential for coarse grained
simulation which doesn't have any corrugation artifact produced as a result
of coarse graining a line into atoms. The basic algorithm is similar to
pair_style line/lj but is implemented for atoms instead of lines.

In order to calculate lj potential between atoms i and j, I need
orientation of each atom with respect to the adjacent atoms. It can be
easily extracted if I can access the coordinates of adjacent atoms (i-1,
i+1 & j-1, j+1).

I am implementing it as follows:
1) get global id of local atom i by tag[i].
2) calculate global ids of adjacent atoms as tag[i]+1 and tag[i]-1
3) convert global ids to local ids, m and n, as atom->map(tag[i]+1) &
atom->map(tag[i]-1).
4) obtain images of m and n closest to i as domain->closest_image(i,m) &
domain->closest_image(i,n).
5) repeat same process for atom j.

The code works fine in serial version. But it crashes after running for
some time step in parallel implementation with exit code 11. I tried
changing several parameters. The only thing that worked was changing
neighbour skin distance from 2 to 20 A. But it increases the computation
cost significantly.

It will be very helpful if you can kindly let me know what might cause
such error. I am using LAMMPS version dated 17Nov16. My system size is
roughly 400,000 atoms. The cutoff for lj potential is 28 A and the style of
atoms is angle with 10 A and 180 degrees as equilibrium bond length and
angle.

​you will get ​a segmentation fault, when the atom->map() function returns
a -1 because the requested atom is neither a local atom nor a ghost atom.
in serial execution all atom ids are local, so atom->map() will never
return -1 for a valid tag. in parallel execution

​i am confused by your use of the term "adjacent". your algorithm looks for
atoms with neighboring atom ids; what in your setup keeps those atoms​
close to each other? what happens for atoms that have tag 1 or the last
tag, i.e. where there is no atom->tag[i]-1 or atom->tag[i]+1?

​axel.​

Hi Axel,

Thanks for the prompt reply.

To address your second question first, I have put checks to make sure if atom tag is 1, I only add to it and if the atom tag is last, I only subtract from it.

Now coming to my definition of adjacent atoms. In my simulation the atoms with neighboring ids are connected by stiff springs and provided that their initial separation is 10 A, a cut-off distance of 30 A should capture those atoms even in their final displaced position.
As a safety check, in case the separation between neighboring atoms becomes greater than 30 A, I extended the ghost cutoff distance to 50 A using comm_modify command. I assume this should have communicated information of such atoms. Am i wrong in this assumption?

Regards,

Ankit Gupta

Hi Axel,

Thanks for the prompt reply.

To address your second question first, I have put checks to make sure if
atom tag is 1, I only add to it and if the atom tag is last, I only
subtract from it.

Now coming to my definition of adjacent atoms. In my simulation the atoms
with neighboring ids are connected by stiff springs and provided that their
initial separation is 10 A, a cut-off distance of 30 A should capture those
atoms even in their final displaced position.
As a safety check, in case the separation between neighboring atoms
becomes greater than 30 A, I extended the ghost cutoff distance to 50 A
using comm_modify command. I assume this should have communicated
information of such atoms. Am i wrong in this assumption?

​without knowing more details about what you are doing, this is difficult
to say.
on the other hand, this situation is easy to debug with adding tests into
your code that checks the output of atom->map() for being negative and then
print out which specific atoms are affected, their atom ids, positions and
distances. please note, that the most extreme case, you have to
accommodate, is that atom i is at the very border of a subdomain, and you
are looking at atom j as a ghost and need to have its two neighbors​ within
the ghost cutoff, too.

another possible source of a problem could be the domain->closest_image()
function. i am not sure, if this is supposed to work for the case of the
reference atom being a ghost atom instead of a local atom. an alternative
could be to compute the distance deltas and then apply
domain->minimum_image() instead.

axel.

Hi Axel,

Yes you were right. Atom->map() was generating a negative value for the case when atom i and atom j are at opposite border of sub-domain. Apparently increasing the ghost-cut off should have taken care of it but I was using comm_modify command wrong. After fixing this mistake, my code is working just fine on parallel machine as well.

Thanks for your help.

Regards,
Ankit Gupta