ERROR: Lost atoms: original 1092 current 1091 (src/thermo.cpp:488)

Hello, may I ask that we use MS to adsorb water molecules in CSH model. The total energy size is -9172.284187 kcal/mol, and the force size is 0.3684 kcal/mol/A. Convert the car file of MS into a data file, and then minimize the energy, the force reaches the power of 17, resulting in an error

units real
dimension 3
boundary p p p
atom_style full

pair_style lj/cut/coul/long 12.5
bond_style harmonic
angle_style harmonic
#kspace_style pppm 1.0e-4
kspace_style ewald 1.0e-4

read_data Layer.data

set type 1 charge -1.05 # ob
set type 2 charge 1.36 # cao
set type 3 charge 1.05 # cah
set type 4 charge 2.1 # Si
set type 5 charge -0.82 # o*
set type 6 charge 0.41 # h*

#group H-O type 3 6
group water type 5 6

neigh_modify delay 0

thermo 100
min_style cg
minimize 1.0e-16 1.0e-16 10000 10000

timestep 0.01
thermo 50
thermo_style custom step temp pe etotal press density lx ly lz

dump 1 all custom 2500 relax0.lammpstrj id type x y z q

fix 1 water nvt temp 298 298 $(100*dt)
fix wshake water shake 0.0001 50 0 b 1 a 1
run 1000
unfix 1

The first thing I would do is to check if the atom configuration in the lammps data file is reasonable. If forces are extremely large at the first step, it’s most likely because some atoms are overlapping on each other.

That sounds a lot like you have incorrect force field parameters or have entered them in different units than what is shown in the “units” command or something similar.

I have adjusted the precision of energy minimization from minimize 1.0e-16 1.0e-16 10000 10000 to minimize 1.0e-4 1.0e-14 10000 10000. The energy minimization can pass normally, but it will still lose atoms to the relaxation stage. I wonder why

Thanks for your reply. Accordingly, I have adjusted the precision of energy minimization from minimize 1.0e-16 1.0e-16 10000 10000 to minimize 1.0e-4 1.0e-14 10000 10000 which can pass normally. But when you get to the relaxation stage, you’re still losing atoms, so why is that

You lose atoms, if they move to fast due to large forces. You have large forces, if your force field parameters or simulation settings are bad. There is no enough information here to say which is the case.

What is very odd is the fact that you set the charges with the “set” command. They should have been set in the data file already. So the place you need to look for possible causes of your problems may be in how you created the data file.

LAMMPS data file. msi2lmp v3.9.10 / 10 Mar 2023 / CGCMM for Layer

1097 atoms
358 bonds
179 angles
0 dihedrals
0 impropers

7 atom types
1 bond types
1 angle types

-0.096661948    22.220738052 xlo xhi
 0.054017672    22.224013272 ylo yhi
-0.652364746    22.117635254 zlo zhi
 0.013968515     0.000000000     0.000000000 xy xz yz

Masses

1 15.999400 # ob
2 40.080000 # cao
3 15.999400 # obos
4 40.080000 # cah
5 28.085500 # st
6 15.999400 # o*
7 1.007970 # h*

Pair Coeffs # lj/cut/coul/long

1 0.1554164124 3.1655200879 # ob
2 0.0000050301 5.5666638636 # cao
3 0.1554164124 3.1655200879 # obos
4 0.0000050219 5.5624343529 # cah
5 0.0000018402 3.3019566252 # st
6 0.1554164124 3.1655200879 # o*
7 0.0000000000 0.0000000000 # h*

Bond Coeffs # harmonic

1 553.9350 1.0000 # o*-h*

Angle Coeffs # harmonic

1 45.7530 109.4700 # h*-o*-h*

Atoms # full

The parameters are taken according to the clayff force field, which is no problem when running some types of models, but when running some models, it will report an error

So the parameters are not suitable for that purpose. If you know that your system will show large forces, you can force frequent output to a dump file where you include forces and then you can easily tell which atoms have the unexpectedly strong interaction.
You can then examine your geometry and draw your conclusions.

For example, since you hydrogen atoms have not Lennard Jones interactions, they can come infinitly close to opposite charges. If you have a molecule/geometry where some anion is present or other atom with significant negative partial charge is close, those may cause the kind of collapse you see. The fact that you are using a timestep that is about 100-200 times smaller than what would be typical is an indication that something is not right.

If “MS” stands for “molecular sieve” then there is an obvious problem with your simulation: you are trying to build a molecular sieve out of isolated LJ+charge atoms without any bonds, angles or dihedrals between them. Those particles will collapse into a chaotic mess and your simulation will crash, as you have seen for yourself (and could see for yourself if you tried to visualise the system).

By the way, it is very bad practice to use acronyms or initialisms without explaining what they mean. I don’t know if “MS” might stand for multiple sclerosis, manuscript, methylsulphide, or Microsoft instead. Not everyone has the energy or inclination to try and make helpful random guesses.

1 Like
  • Thank you for your wonderful answers
  • Thank you for your answers