Variable for the bounds of region style 'block' and for 'final' style of command 'change_box'


In order to save computational resources, I need to resize the simulation box once every N steps for a system with a constant total density. One way of doing this is to run the command ‘change_box’ in a loop of N timesteps. Since the system density is not perfectly homogeneous, one needs to adjust the number of particle deletions or insertions by the commands ‘fix evaporate’ or ’ fix deposit’, respectively, so that the total density remains constant after shrinking the box. These commands need a region to be specified as the deletion/insertion region and a proper style of the region would be ‘block’ if the box dimension is orthogonal.

I think that enabling the equal-style variable for the region bounds of style ‘block’ as well as for the final values of command ‘change_box’ would be one way to achieve the desied constant density. Please see the script below as an example:

1- variable rho equal 10
2- variable r_1 equal 90
3- thermo_modify lost ignore

4- region deletion block -95.0 95.0 -95.0 95.0 -95.0 95.0 side out
5- delete_atoms region deletion compress no
6- change_box all x final -95.0 95.0 y final -95.0 95.0 z final -95.0 95.0
7- variable n_All equal count(all)
8- variable n_desired equal ${rho}*${190}^3
9- region reg1 sphere 0 0 0 ${r_1} side out
10- variable n equal abs(v_n_All-v_n_desired)

11- label loop1
12- variable try1 loop $n
13- if '${n_All}<${n_desired}' then 'fix dep1 all deposit 1 1 1 12345 region reg1' else 'fix evap1 all evaporate 1 1 reg1 54321'
14- run 1
15- next try1
16- jump SELF loop1

In lines 4 and 6, enabling equal-style variables for the new box dimensions allows us to repeat this block of the script in a systematic way such that resizing the box can be performed once a certain condition is met.

I would like to ask the developers to include such changes in the future versions of LAMMPS, in case there is no better way of keeping density constant when resizing the simulation box.

Thank you in advance for any suggestions.

I fail to be convinced by your reasoning.

If you have periodic boundaries then shrinking the box and deleting atoms at the box “surface” will cause unphysical shockwaves due to the sudden change of environment of the neighboring atoms near the boundaries and thus require significant time spent on re-equilibrating the system. So where are the time savings?

In case you have non-periodic boundaries, then density is ill-defined since the boundary is essentially arbitrary (and the density zero since the system is infinite) so why would there be a need to keep it constant within the box?

thank you for your reply, Axel. Generally, I agree with your point.
For the system I am simulating, the most important phenomena of interest occur far from the box boundaries (a microparticle centered at the box). In addition, I am running DPD simulations for which the relaxation is fast enough such that the shock waves would not affect the microparticle. So I think there are cases where the sudden change of box bounds would be useful.

I disagree. What you say is all conjecture and no proof. To merely save some computational effort (and then how much?) does in my opinion not justify the risk of invalidating the correctness of your findings.
This is even more true given that it is currently quite easy (compared to the situation 10 or more years ago) to get access to significant additional computational resources either that the institutional or national level.

I think it would be doing LAMMPS users a disservice to make such operations easier. It is an extreme corner case, anyway.

Thank you for your reply.

Now I think you are right that is an extreme corner. Therefore, I changed the source code of relevant fixes to achieve this shrinking procedure. I also found that the effect of shock waves vanishes quickly far before reaching the centered microparticle and, besides, the final configuration remains unaffected compared to the case without shrinking protocol.

Thank you again.