[lammps-users] parallelization and bond lists

Hi. I am trying to simulate reversible (polymer) networks in LAMMPS as follows. Each chain end atom starts off unbound, but can form bonds with any other chain end. All ~ (# of chains)^2 _possible_ interchain bonds are included in the initial data file. I then use measures similar to those found in bond_quartic.cpp to turn interchain bonds on and off (in the sense that forces are imposed or not). This all works fine on one processor. However, in parallel LAMMPS crashes immediately with a "Bond atoms [i j] missing on proc [k] at step 0" message (generated from line 45 of neigh_bond.cpp). This is not terribly surprising, because some "bonded" atoms (with the bond inactive, i. e. forces are not computed) may be on distantly separated chains and belong to different processors. Question is, is this scheme (i. e. listing all _possible_ bonds) doomed to break parallelizability, or is there something I can do to make this work?

Thanks,
Rob

You could alter the logic in neigh_bond::bond_all() to ignore
the bond when it can't find it, rather than throw an error.

Steve

Hi, Steve. I tried that. The best way I could think to do it was to modify neigh_bond::bond_all() as follows:

   for (i = 0; i < nlocal; i++)
     for (m = 0; m < num_bond[i]; m++) {
       atom1 = atom->map(bond_atom[i][m]);
       if (atom1 == -1) {
// char str[128];
// sprintf(str,"Bond atoms %d %d missing on proc %d at step %d",
// tag[i],bond_atom[i][m],me,update->ntimestep);
// error->one(str);
       }
       if (newton_bond || i < atom1) {
         if (nbondlist == maxbond) {
           maxbond += BONDDELTA;
           bondlist = memory->grow_2d_int_array(bondlist,maxbond,3,
                                                "neighbor:bondlist");
         }
        if (atom1 != -1) {
            bondlist[nbondlist][0] = i;
            bondlist[nbondlist][1] = atom1;
            bondlist[nbondlist][2] = bond_type[i][m];
            nbondlist++; }
       }
     }

The last if statement is there because obviously atoms not owned by the executing processor shouldn't be added to the bond list. However, in parallel this caused a "Lost atoms" error (where the final number of atoms was GREATER than the true #) the first time bond_all() was called. If I remove the if (atom1 != -1) condition (which would, I think, cause other problems) there is still the same lost atoms error. How were you thinking I should modify this?

Thanks,
Rob

Quoting Steve Plimpton <[email protected]>: