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?