temper and Periodic boundary condition

Dear lammps users and developers,

I’m trying to write a modification of the temper.cpp code where, before attempting a swap between the replicas, I drive the system out of equilibrium by heating (or cooling) the system to the partner’s temperature. After heating or cooling the system to the target temperature if the swap is rejected the system has to go back to the state before the heating/cooling part started. If the swap is accepted the system keep evolving at constant temperature.

I have set two arrays where I copy the position and velocities of all the atoms before starting the heating/cooling part. If the swap is rejected I copy this positions and velocities to the lammps x and v arrays. Also, after copying the coordinates and velocities I’m using neighbor->decide() and neighbor->build() to determine if I need to reconstruct neighbors list and rebuild it if necessary. To deal with the pbc I have tried using domain->pbc() and domain->remap(x[i],image[i]).

The code I’m using is:

if (n_pre_exchange) modify->pre_exchange();
if (triclinic) domain->x2lamda(atom->nlocal);
domain->pbc(); or remap
//for (i = 0; i < nlocal; i++) domain->remap(x[i],image[i]);
if (domain->box_change) {
domain->reset_box();
comm->setup();
if (neighbor->style) neighbor->setup_bins();
}
timer->stamp();
comm->exchange();
if (sortflag && ntimestep >= atom->nextsort) atom->sort();
comm->borders();
if (triclinic) domain->lamda2x(atom->nlocal+atom->nghost);
timer->stamp(TIME_COMM);
if (n_pre_neighbor) modify->pre_neighbor();
neighbor->build();
timer->stamp(TIME_NEIGHBOR);

Nevertheless, the system explodes after a couple of attempts when the swap is rejected. The error message that I get is:
“ERROR on proc 0: Bond atoms 14396 14394 missing on proc 0 at step 1002 (neigh_bond.cpp:55)”

Some extra facts:

- I’ve tried writing a restart file just before the simulation explodes and then reading this file. When I do this the simulations keep on running beyond the point where usually explodes. I believe the problem is either in building the neighbors list or because there is an overlap between system atoms and images/ghost atoms. I have gone over the read_restart.cpp and read_dump.cpp to see the different steps necessary set up domains, pbc and neighbors lists but I haven’t able to get it to work.

- The atoms/bonds that go missing are always water molecules that are at the edge of the box

- Also, if I set my simulations in such a way that the swaps are always accepted I runs just fine.

- Reducing the time-step does not solve the problem

I’m running the latest version of lammps.

Any help is very very appreciated,
Regards,

Ignacia Echeverria

There are many things that can go wrong with trying to save a snapshot
and restore it later. The prd and tad commands do this carefully.
It is especially tricky if you are running with multiple
processors per replica which temper allows you to do.

You will have to debug what is going wrong
and see why you are getting the error message, which is typically
b/c LAMMPS cannot find an atom it needs.

Steve