Bug report on fix evaporate while working on large system

Dear Lammps developers and users,

I’m going to model salt (CaSO4) precipitation in aqueous solution using LAMMPS, stable version 12 Aug 2018.
To do so, firstly a large system consisting of 66000 water molecules dissolving 10 CaSO4 was generated. Then, to model super-saturation, fix evaporate was employed to reduce the number of water molecules, that follows, and consequently increase the concentration of dissolved CaSO4:

fix 3 owater evaporate 100 9 surface 38277 molecules yes

here, owater is group_ID of the water’s oxygen atoms and “surface” is the whole simulation box. In this way, 3 water molecules would be deleted at every 100-time steps.

Due to the restriction imposed by mailing list on the size of attachments, I cannot send associating input and data file, but, it will be OK to send anyone if want.

Unfortunately, the simulation goes wrong! the number of removed atoms is more than the desired one, i.e., 9. Also, the deleting manner is apparently erratic, deleted number of molecules varies at every time step.

I tried that scenario for a smaller system with 6000 waters and 25 CaSO4 and surprisingly the system works well, deleting atoms as invoked by preceding fix evaporate.

I suppose this issue is a bug with the large system.

Any hint is welcome.

Many thanks in advance

Dear Lammps developers and users,

I'm going to model salt (CaSO4) precipitation in aqueous solution using LAMMPS, stable version 12 Aug 2018.
To do so, firstly a large system consisting of 66000 water molecules dissolving 10 CaSO4 was generated. Then, to model super-saturation, *fix evaporate* was employed to reduce the number of water molecules, that follows, and consequently increase the concentration of dissolved CaSO4:

fix 3 owater evaporate 100 9 surface 38277 molecules yes

here, owater is group_ID of the water’s oxygen atoms and "surface" is the whole simulation box. In this way, 3 water molecules would be deleted at every 100-time steps.

Due to the restriction imposed by mailing list on the size of attachments, I cannot send associating input and data file, but, it will be OK to send anyone if want.

Unfortunately, the simulation goes wrong! the number of removed atoms is more than the desired one, i.e., 9. Also, the deleting manner is apparently erratic, deleted number of molecules varies at every time step.

I tried that scenario for a smaller system with 6000 waters and 25 CaSO4 and surprisingly the system works well, deleting atoms as invoked by preceding fix evaporate.

I suppose this issue is a bug with the large system.

there is nothing here that suggests a bug. bug that happen only with
larg(er) system sizes are *extremely* rare, and fix evaporate has been
around a very long time, so if there was a bug that obvious, people
would have noticed it.

it is *much* more likely, that there is something wrong or
inconsistent with your system or your input.

can you please insert the following commands into your input after the
fix evaporate command, right before the run command and report the
resulting log file output?

compute allmol all chunk/atom molecule nchunk every compress yes
compute watmol owater chunk/atom molecule nchunk every compress yes
compute mols all reduce max c_allmol c_watmol

thermo_style custom step pe temp atoms c_mols[*]
thermo 100
thermo_modify lost warn

run 1000

axel.

Dear Axel
Thanks for your kind reply
I followed your commands and the resulting log file is attached.
Would you please hint me on this issue?

log.lammps (9.05 KB)

Dear Axel
Thanks for your kind reply
I followed your commands and the resulting log file is attached.
Would you please hint me on this issue?

the problem is, that you have in your data file multiple molecules with the same molecule id.
also, the molecule ids of the Ca and SO4 ions overlap with those of water.

Step PotEng Temp Atoms c_mols[1] c_mols[2]
100 -315821.81 240.70073 180060 9999 9999

have a look at the last two columns. the value for c_mols[1] should be larger than the one for c_mols[2] (by 20, if you have 10 Ca cations and 10 SO4 anions) and it c_mols[2] should be 60000 corresponding the the number of oxygen atoms.

200 -326904.82 286.01374 180042 9998 9998
300 -330464.59 314.92403 180024 9997 9997
400 -331686.81 335.063 180006 9996 9996
500 -334151.02 353.67335 179988 9995 9995
600 -336489.4 358.42022 179970 9994 9994
700 -340515.98 356.7316 179952 9993 9993
800 -345255.97 356.84949 179934 9992 9992
900 -349808.68 356.34511 179916 9991 9991
1000 -354401.71 356.24121 179898 9990 9990
1100 -358610.55 355.26899 179880 9989 9989

as you can see in every invocation of fix evaporate, one “molecule” is removed. this happens, because fix evaporate removes one atom, then all atoms with the same molecule id, then find that this is equal or more than the required 9 atoms and then stop removing until 100 more steps have passed.

the number 9999 is quite special, as that is the largest residue id supported by the .pdb file format.
so let me speculate, that you created your system in .pdb format and then converted it to a data file using the residue id as molecule id. most .pdb based tools make atoms unique by assigning different segment ids, but LAMMPS does not support that.

for your smaller system, this will have (mostly worked), since this molecule id wraparound does not apply to all atoms, but occasionally might have hit molecules with shared molecule ids (assuming you prepared the input the same way).

in short, this is certainly not due to a bug in LAMMPS, but due to an incorrectly prepared data file.

axel.