lost atoms error when using 'fix reax/c/species ' command

Dear Lammps-User,
My name is wenlong Chen. I came across a ‘lost atoms’ error when I tried to simulate the mixing of water and silica.

This happens when the ‘fix reax/c/species’ command is on. Its weird since the fix reax/c/species command does not change the update of the position of the atoms but only calculate the species based on the obtained bond order. Below is my input script as well as the log file. Any help/discussion will be very much appreciated.

Input script:

Initialization

units real
dimension 3
boundary p p p
atom_style charge

lattice diamond 1.0 origin 1e-9 1e-9 1e-9 spacing 3.0 3.0 3.0

Atom Definition

read_data data.water extra/atom/types 2

region silicabox block EDGE EDGE 10.0 22.0 EDGE EDGE units box

molecule sili silica.mol toff 2

mass 3 28.0855
mass 4 15.9994

delete_atoms region silicabox
create_atoms 0 random 80 96896 silicabox mol sili 45892 units box

pair_style reax/c lmp_control
pair_coeff * * ffield.reax.AgSiCHON H O Si O

neighbor 2.0 bin
neigh_modify every 10 delay 0 check no

variable dt equal 0.2
variable Nrun equal 100000

timestep ${dt}

Settings

fix 1 all nvt temp 300.0 300.0 100.0
fix 2 all qeq/reax 1 0.0 10.0 1e-6 reax/c
#fix 3 all temp/berendsen 300.0 300.0 100

fix 4 all reax/c/bonds 1000 bonds.reaxc
fix 5 all reax/c/species 100 10 1000 species.out element H O Si O

dump 1 all atom 500 silica.lmptrj

Output

thermo_style one
thermo 5000
#thermo_modify lost ignore

Run the simulation

run ${Nrun}

#write_data data.silicawatEq
#write_restart SilicaWatMix.res

log.lammps

LAMMPS (3 Mar 2020)
OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (…/comm.cpp:94)
using 1 OpenMP thread(s) per MPI task

Initialization

units real
dimension 3
boundary p p p
atom_style charge

lattice diamond 1.0 origin 1e-9 1e-9 1e-9 spacing 3.0 3.0 3.0
Lattice spacing in x,y,z = 3 3 3

Atom Definition

read_data data.water extra/atom/types 2
orthogonal box = (-0.799456 -1.16546 0.204719) to (33.1695 32.8035 10.0693)
2 by 3 by 1 MPI processor grid
reading atoms …
1020 atoms
reading velocities …
1020 velocities
read_data CPU = 0.00765974 secs

region silicabox block EDGE EDGE 10.0 22.0 EDGE EDGE units box

molecule sili silica.mol toff 2
Read molecule sili:
3 atoms with max type 4
0 bonds with max type 0
0 angles with max type 0
0 dihedrals with max type 0
0 impropers with max type 0

mass 3 28.0855
mass 4 15.9994

delete_atoms region silicabox
Deleted 365 atoms, new total = 655
create_atoms 0 random 80 96896 silicabox mol sili 45892 units box
Created 240 atoms
create_atoms CPU = 0.000167551 secs

pair_style reax/c lmp_control
pair_coeff * * ffield.reax.AgSiCHON H O Si O

neighbor 2.0 bin
neigh_modify every 10 delay 0 check no

variable dt equal 0.2
variable Nrun equal 100000

timestep ${dt}
timestep 0.2

Settings

fix 1 all nvt temp 300.0 300.0 100.0
fix 2 all qeq/reax 1 0.0 10.0 1e-6 reax/c
#fix 3 all temp/berendsen 300.0 300.0 100

fix 4 all reax/c/bonds 1000 bonds.reaxc
fix 5 all reax/c/species 100 10 1000 species.out element H O Si O
WARNING: Resetting reneighboring criteria for fix reax/c/species (…/fix_reaxc_species.cpp:97)

dump 1 all atom 500 silica.lmptrj

Output

thermo_style one
thermo 5000
#thermo_modify lost ignore

Run the simulation

run {Nrun} run 100000 Neighbor list info ... update every 1000 steps, delay 0 steps, check no max neighbors/atom: 2000, page size: 100000 master list distance cutoff = 12 ghost atom cutoff = 12 binsize = 6, bins = 6 6 2 2 neighbor lists, perpetual/occasional/extra = 2 0 0 (1) pair reax/c, perpetual attributes: half, newton off, ghost pair build: half/bin/newtoff/ghost stencil: half/ghost/bin/3d/newtoff bin: standard (2) fix qeq/reax, perpetual, copy from (1) attributes: half, newton off, ghost pair build: copy stencil: none bin: none Per MPI rank memory allocation (min/avg/max) = 88.52 | 89.57 | 90.54 Mbytes Step Temp E_pair E_mol TotEng Press 0 214.91895 -57882.045 0 -57309.32 -6367858.3 **ERROR: Lost atoms: original 895 current 886 (../thermo.cpp:438)** Last command: run {Nrun}

Kind regards
Wenlong Chen

Dear Lammps-User,
My name is wenlong Chen. I came across a ‘lost atoms’ error when I tried to simulate the mixing of water and silica.

This happens when the ‘fix reax/c/species’ command is on. Its weird since the fix reax/c/species command does not change the update of the position of the atoms but only calculate the species based on the obtained bond order. Below is my input script as well as the log file. Any help/discussion will be very much appreciated.

fix reax/c/species does not change the positions, but if you average over data, it requires that the neighbor list is not changed/updated during this period. Thus it changes the neighbor list update settings. You can see a corresponding warning in your log and also will see that the neighborlist update frequency changed. In your case to a rather unreasonable value and that is why you see the lost atoms.

Axel.