Periodic boundaries and shake atoms missing

Hello,

I am modeling SPC/E water with periodic boundary conditions, NVT, and I ran into trouble when I tried using change_box to expand the box (still keeping the PBC). I was getting an error like this:

ERROR on proc 323: Shake atoms 24847 24848 24849 missing on proc 323 at step 1250 (…/fix_shake.cpp:535)

This error would repeat a few hundred times in the log. I would not get a lost atoms warning however. After trying to change my neighbor lists/cutoff distances/timestep/trying different # processors, I found by looking at a snapshot that some hydrogen atoms were entirely outside the box boundary. So I tried deleting all atoms around the perimeter of the box and this allowed the change_box command to work without the shake error.

I am wondering if there is a way to prevent atoms from extending slightly past the box boundary so I don’t run into this error, or any other possible solutions to look into.

Thank you,

Eric Bird

Hi Erica,

Without an input script it is hard to diagnose exactly what your problem is; however, did you have the remap option applied to your change_box command? Something else to consider is just building your system at the desired size you would set the change_box command to instead, and see if that also throws the error.

Otherwise, what lammps version are you using, and would you be able to supply a copy of your input (data file + input script) if it is small enough? I’ve had shake atom errors come up with SPC/E water in the past, usually due to a bug in my input somewhere, so it would be helpful to see if that is the case.

Best,
Zeke

Hello Zeke,

Thank you for responding to my question.

Typically I like to equilibrate a liquid system in a cube and then use change_box to create a liquid-gas interface. I tried using the remap command before on change_box and it prevented the shake error, but this is not my desired effect of using change_box. I want all water molecules to retain their absolute position after the change_box. I am unaware of an option to remap atoms yet keep their original positions pre-change_box. Also I would like to use the same data file on multiple types of simulations if possible and just replicate, say 2nm of water, to whatever dimension I want. So yes I could place the water in the full sized box to begin with, but this is a bit inconvenient because I would have to make a new data file for several types of simulations.

I am using lammps/16Mar18. I included my data file and input script in the email. In the included input script, I am trying equilibrate a cube of water (about 2nm side length before replication) and then expand the box in all directions, eventually creating a water droplet. The input script as-is should give the shake error, and I commented out the code I used to delete the boundary molecules, which prevents the error. Pardon any confusion on account of all the variables.

I’m mostly just looking for a cleaner solution than having to delete atoms before I change_box. I asked a colleague of mine and he said he did a similar workaround where he replicated the cube of water an extra 2 times in each direction and deleted atoms in the replicated regions.

Thank you,

Eric Bird

C.Water.Droplet.Test (4.76 KB)

r.Water.Setup.dat (921 KB)

Hi Erica,

I looked through your input scripts and didn’t see something particularly obvious as the cause of the issue, your system is pretty big so I didn’t end up running it. I’m wondering if it is a case of the sudden change in your system size (going straight from water box to vacuum). You might consider using fix deform (https://lammps.sandia.gov/doc/fix_deform.html) to try changing the box size gradually over time instead of all at once. This might help prevent the issue you are causing. I’m also not sure how change_box interacts with the way the shake atoms are handled by the processors based on the note in the documentation for that command (The simulation box size/shape can be changed by arbitrarily large amounts by this command. This is not a problem, except that the mapping of processors to the simulation box is not changed from its initial 3d configuration; see the processors command. Thus, if the box size/shape changes dramatically, the mapping of processors to the simulation box may not end up as optimal as the initial mapping attempted to be. You may wish to re-balance the atoms by using the balance command if that is the case.). ​

I hope one of these two fixes solves your problem; otherwise, hopefully someone else will have something for you to try!

Best,
Zeke

please note that the change_box command (as well as fix deform) will by default rescale all coordinates of the atoms in the group they operate on according to the change in the box dimensions. that is the documented behavior but probably not desired for the scenario described in this discussion thread.
it may be helpful to first define an empty group (e.g. with group none empty) and then operate on that group, so that all atoms stay in place.

axel.