Hello,
In the pair_style I am working to develop I need to run a loop over all atoms connected to a given atom by bonds. I tried to use bond_atom property but it does not work for atom ID greater than natoms. Can anyone advise?
Thank you,
Mark Chaimovich.
This statement doesn’t make sense in two ways. First, the “bond_atom” array is a per-atom property and thus when are looping over local atoms you can get the “i” atom’s atom ID from atom->tag[i] and then loop over the bonded atoms and get their atom IDs. Second, the list is only complete when you require newton off for bonds. Otherwise bonds are only stored with one atom, not both bond partners.
That said, If you need to loop over atoms that are bonded neighbors of a given (local) atom, you can use the atom->special[][] and atom->nspecial[][] arrays. The entry atom->nspecial[i][0] says how many 1-2 neighbors atom “i” has and atom->special[i][j] contains their atom IDs. The atom->special[][] array also contains 1-3 and 1-4 neighbors (connected to atom “i” through two or three bonds, respectively).
The compute fragment/atom command for example, uses this information.
Axel,
Thank you for the quick response.
I tried it and special and special lists work good for atoms in the jlist that are less than natoms.
I have natoms = nlocal =1728 (I work with one processor).
This is the some output from the loop that runs over atoms in jlist (here j is the index taken from jlist: j = jlist[jj];).
j 42 molecule[j] 18 tag[j] 70 nspecial[j][0] 2
bonded to atom->special[j] 71 map(atom->special[j]) 43
atom->special[j] 72 map(atom->special[j]) 44
j 46 molecule[j] 77 tag[j] 305 nspecial[j][0] 4
bonded to atom->special[j] 306 map(atom->special[j]) 463
atom->special[j] 307 map(atom->special[j]) 47
atom->special[j] 308 map(atom->special[j]) 48
atom->special[j] 309 map(atom->special[j]) 49
j 3099 molecule[j] 84 tag[j] 334 nspecial[j][0] 0
Thank you,
Mark Chaimovich.
I can’t read this.
1 Like