After delete_atoms command is executed and some atoms are deleted,
(which is done by overwriting unneeded atom info -- avec->copy(nlocal-1,i,1) )
it's very likely that some bonds will point to global-atom-ID which is
not there anymore.
I have two questions:
After the atom_map is reset in delete_atoms.cpp with no compression,
(atom->map_init(); atom->map_set(); )
1) when I access a bond via
deletedID = atom->bond_atom[localID][k],
atom->map(deletedID) should give me "-1",
so I can tell that this atom was deleted, right?
2) With regards to usual calculation of bond forces (e.g. bond_harmonic.cpp)
via
int **bondlist = neighbor->bondlist;
int nbondlist = neighbor->nbondlist;
for (n = 0; n < nbondlist; n++) {
i1 = bondlist[n][0];
i2 = bondlist[n][1];
type = bondlist[n][2];
...
}
Does bondlist have updated info on deleted bonds and will not include
bonds with atoms which were deleted? So calculation will go ok and
deleted atoms will not spoil force calculation?
Probably neighbor->anglelist / neighbor->nanglelist
behave exactly the same as bondlist/nbondlist? And should also be ok?
p/s/
If not, what would be a simplest way to fix bonds after deleting atoms?
It's straight forward to go through the bonds of local atoms (in delete_atoms.cpp) and
change atom->num_bond[localID] and rearrange atom->bond_type[][] and atom->bond_atom
appropriately, but i'm not sure if it's enough or do I need to do it at all?
If not, what would be a simplest way to fix bonds after deleting atoms?
The answer is not to delete atoms that are in bonds. If you do,
that's a user error. There is a delete_bonds command. Use that
to delete the bonds you want, then delete the unbonded atoms.
I looked at delete_bonds command, and I don't see how can it help me.
Say, I want to make a hole (delete_atoms inside a sphere). So, I
make a group which corresponds to that sphere and do:
delete_bonds mySphere bond 1 remove
But it will not delete bonds between
atoms inside the sphere and the ones which are outside, right?
"For all styles, an interaction is only turned off (or on) if all the atoms involved are in the specified group"
Is there a simple way to accomplish what I need?
I can also modify delete_bonds.cpp, more precisely the following condition:
if (mask[i] & groupbit && mask[atom1] & groupbit) ->
if (mask[i] & groupbit || mask[atom1] & groupbit)
Is there a reason why it's "&&" only?
Going through the code I think it's also possible to have "||" there.
Your comments are correct. The code in both
delete_bonds and delete_atoms is trying to
protect you from dong the bad thing of deleting
only some of the atoms in an interaction (bond, angle, etc).
But maybe there needs to be some additional options
to enable carving out voids in systems with bond (angles, etc)
that cross the boundary.