neighbor list problem

Dear all. Following Axel suggestions, I've implemented a preliminary part of mesvb method by using different replica.
Now I am doing a preliminary test that consists in running two replica of a single bond system.
In my system I have one bond(O1-H) associated with the first replica while the second bond is O2-H.
However at the beginning both replicas have O1-H bond while O2 need to be found by a specific routine and then O2-H can be built.
After finding the atom O2 I change the topologies

      swap_topologies(atom1,atom2);

and rebuild the neighbor list

    neighbor->build();

The routine swap_topologies will swap the properties of the two atoms, i.e. q, type, bond_type, num_bond, bond_atom, bondlist,molecule,
tag, sametage. Thes is something missing?

In order to test if this approach works I did two test. The first test consist in running two lammps simulation, one with the bond O1-H and the other
with the bond statically built in O2-H. The second test consist in running the following: mpirun -np 2 lmp -partition 2x1 -in lmp.in. I ran for 1 timestep with timestep=0.0.
So the energies should be the same.
With the first test the bond are built at the beginning of the simulation while with the partition run the second bond is defined after the routine swap_topologies is called.

I printed out bonds, kspace, eng_vdwl, eng_coul energies. Bonds and kspace are computed correctly while eng_vdwl and eng_coul are not after the swapping.
However if I loop over the neighbor list and I change the entries of jlist the energies are correct, e.g.

    for(int ii=0; ii<inum; ii++){
      i=ilist[ii];
      jlist = firstneigh[i];
      jnum = numneigh[i];
      for (int jj = 0; jj < jnum; jj++) {
        j = jlist[jj];
        t = j&NEIGHMASK;
      if(i==oxygens_id_local[nh][k][0] && t==hydrogens_id[nh])
        {
          jlist[jj] = j_atom2[k];
        }
     if(ii==oxygens_id_local[nh][k][1] && t==hydrogens_id[nh]){
          jlist[jj] = j_atom1[k];
        }
}

This approach looks ugly and dangerous to me but it seems to work. This suggest me that the operation neighbor->build does not rebuild the topology correctly
when the configurations are swapped. Any suggestion on how rebuild correctly the neighbor list in the swapped configuration?

Thanks very much

Nicola

[...]

The routine swap_topologies will swap the properties of the two atoms,
i.e. q, type, bond_type, num_bond, bond_atom, bondlist,molecule,
tag, sametage. Thes is something missing?

the exclusion settings in special[][]?

[...]

I printed out bonds, kspace, eng_vdwl, eng_coul energies. Bonds and
kspace are computed correctly while eng_vdwl and eng_coul are not after
the swapping.
However if I loop over the neighbor list and I change the entries of
jlist the l are correct, e.g.

that confirms that your exclusion settings are not updated.

    for(int ii=0; ii<inum; ii++){
      i=ilist[ii];
      jlist = firstneigh[i];
      jnum = numneigh[i];
      for (int jj = 0; jj < jnum; jj++) {
        j = jlist[jj];
        t = j&NEIGHMASK;
      if(i==oxygens_id_local[nh][k][0] && t==hydrogens_id[nh])
        {
          jlist[jj] = j_atom2[k];
        }
     if(ii==oxygens_id_local[nh][k][1] && t==hydrogens_id[nh]){
          jlist[jj] = j_atom1[k];
        }
}

This approach looks ugly and dangerous to me but it seems to work. This

this would only work for as simple a case as
you have now and in general is incorrect.
what if the system goes back to the original
topology?

suggest me that the operation neighbor->build does not rebuild the
topology correctly
when the configurations are swapped. Any suggestion on how rebuild
correctly the neighbor list in the swapped configuration?

have a look at fix bond_swap.cpp this shows how to swap two bonds.

axel.