Using a neighborlist in a custom fix

Dear LAMMPS users,

I am developing/modifying a fix which needs a separate neighborlist, since only a much smaller cutoff distance is needed (around 12A for short range Coulomb and vdW and around 5.5A for the custom fix). I request the neighborlist in my fix by the following lines of code in the init method of the fix:

 int irequest = neighbor\->request\(this,instance\_me\);
 neighbor\->requests\[irequest\]\->pair=0;
 neighbor\->requests\[irequest\]\->fix=1;
 neighbor\->requests\[irequest\]\->half=0;
 neighbor\->requests\[irequest\]\->full=1;
 neighbor\->requests\[irequest\]\->cut=1;
 neighbor\->requests\[irequest\]\->cutoff=cutoff\_radius; // usually around 5\.5A

To get it working properly, I have to configure the neighborlist to rebuild in every time step by invoking the following command:

 neigh\_modify every 1 delay 0 check no

else, I get a corrupted neighbor list in the sense that pairs are missing.

Do you have any idea how to get such a setup running without the need for rebuilding the neighborlist in every step?

Best

Johannes

Have you tried setting the flag

neighbor->requests[irequest]->occasional = 1;

and then use the following to build the neighbor list as needed?

neighbor->build_one(list);

axel.

That’s odd, the perpetual neighbor list should be reconstructed according to the same skin trigger as the default pair style list etc. Which neighbor list object pointer did you end up using in the implementation of your fix? theres also an empty init list method in the fix parent you can inherit and implement to get the neighbor list pointer for your full list or just simply dereference neighbor->lists[irequest] since “lists” has public access as I recall. As long as you didn’t set occasional to 1 on the request it should’ve led to a neigh list object thats considered perpetual and rebuilt automatically to be accurate.

Adrian Diaz

Judging from what I’m seeing in nbin full and npair full it seems like your custom cutoff might need to explicitly include the skin distance (default 2.0 set in lammps unless you altered it with the neighbor command in your input script); this might explain why its working with a rebuild every timestep. If you specify 7.5 as your cutoff it should work I think since that’s synchronized with the rebuild trigger by adding the skin distance.

Adrian Diaz

Dear all,

thanks for your suggestions. I tried already the hints from Axel, but they do not solve problem. I will now try what is happening when explicitly including the skin distance.

Best,

Johannes

You can actually just use neighbor->requests[irequest]->cutoff=cutoff_radius+neighbor->skin; to include the skin rather than worrying about what the skin is for your run.

Adrian Diaz.

Thank you very much for your help. It seems that adding the skin distance solves the problem!

I have another neighborlist related question:

I know that I can get via atom->map(atom->tag[i]) the real index of atom i in case that it is a ghost atom. Is there also the possibility to get the relative position of the ghost atom with respect to current domain or for a serial run with respect to the unitcell? For example if the original atom sits at (0.25, 0, 0) and one of its ghost atoms at (1.25, 0, 0) in fractional coordinates, then the relative position of this ghost atom with respect to the original atom would be (1,0,0).

Is this saved somewhere or do I need to compute this?

Best,

Johannes

That’s good news. Not sure what you mean at the moment, however the tags and positions of all ghosts are available in the tag and x arrays respectively. Not sure if this is known or not but you can access all the ghosts by considering the array access between atom->nlocal and atom->nlocal+atom->nghost to find all the ghost atoms. The neighborlist indices will also already lie within this range for ghosts and the tag and x arrays are allocated to accommodate that information.

Adrian Diaz

I have another neighborlist related question:

I know that I can get via atom->map(atom->tag[i]) the real index of atom i in case that it is a ghost atom. Is there also the possibility to get the relative position of the ghost atom with respect to current domain or for a serial run with respect to the unitcell? For example if the original atom sits at (0.25, 0, 0) and one of its ghost atoms at (1.25, 0, 0) in fractional coordinates, then the relative position of this ghost atom with respect to the original atom would be (1,0,0).

this is not really a neighborlist question, but a question about domain decomposition. In general, you cannot know whether a ghost atom is a periodic replica of some atom or just a copy from a neighboring domain. technically, it always is the latter. LAMMPS doesn’t pay any attention and is not subject to minimum image conventions (so there may be multiple copies in case of a small cell and a large cutoff).
you could compare its coordinates with domain->box{lo,hi} numbers to see whether it is inside or outside the principal box, but even that is not a safe measure, since atoms can move slightly outside the box between neighbor list rebuilds, so you can only trust this number if neighbor->ago is 0.
…and that only works for rectangular boxes, in the general case, you would have convert the position to lambda coordinates (aka fractional coordinate), but even then the previously mentioned caveat remains.

axel.

Is this saved somewhere or do I need to compute this?

it