Lost atoms error

Hi all,

I have been trying to generate an a-SiO2 slab using a forcefield I found on a paper using lammps. First step is to melt a crystobalite sample. I am fairly new to lammps and my input file to melting a crystobalite slab is as below.

units           real
atom_style      full

boundary        p p p

pair_style      hybrid/overlay buck 8.0 coul/wolf 0.231714 14.0
read_data       crys.dat

pair_coeff      1  1  buck 64522.092 2.893 0.0  # SI-SI
pair_coeff      2  2  buck 25839.738 4.407 601.88949  # O-O
pair_coeff      1  2  buck 532886.67 5.098 3221.6077  # O-Si

pair_coeff      *  *  coul/wolf

neighbor        2.0 bin
neigh_modify    every 1 delay 0 check yes

group mobile    type 1 2

velocity        mobile create 298 1103 dist gaussian

fix             1 mobile nve
fix             2 mobile press/berendsen iso 1.0 1.0 1000 modulus 360000
fix             3 mobile temp/berendsen 298 4000 100

timestep        1.0

run_style       verlet

thermo_style    custom step cpu temp pe ke press vol density lx ly lz
thermo          100

dump            traj mobile xyz 100 melt.xyz
dump_modify     traj element Si O
restart         500000 restart.1 restart.2

run             4000000type or paste code here

After running this I keep getting the following error.

Step CPU Temp PotEng KinEng Press Volume Density Lx Ly Lz 
       0            0          298 2.664738e+09     510.7618 1.8095785e+09    35332.099   0.54216622      21.0155      21.0155           80 
ERROR: Lost atoms: original 576 current 574 (src/thermo.cpp:427)
Last command: run               4000000

can someone help me to resolve this error?

You potential energy is extremely large. That is usually an indication of a bad initial geometry.

Interesting, the initial geometry is a completely ordered crystalline structure. Thank you for the comment Dr. Akohlmey.

Try dumping the per-atom potential energy. Perhaps you can identify some areas where atoms have high potential energy (e.g. at the box boundaries due to miscalculated box dimensions).

A simple test for close contacts is to insert a delete_atoms overlap 0.2 all all command before the first minimize or run command. If it deletes atoms, the there are some close contacts.

When using a slab, you usually don’t need a barostat and using isotropic change is a bad idea anyway.
Finally, you have to check whether the force field parameters you use are suitable for a slab or only for bulk systems. For surfaces, you often either have models with explicit bonds or some form of polarization or charge equilibration.

1 Like

Thank you very much for the Ideas.