Lammps: Lost atoms with fix wall/lj93

Dear Users,

I want to simulate a liquid crystal composed of ellipsoidal particles that interact through the Gay-Berne potential. The system has periodic boundary conditions except on the axis-z, the limits in this dimension are walls that interact with the particles through lj93 potential. The problem is that the system lose atoms (ERROR: Lost atoms: original 256 current 238 (…/thermo.cpp441)) Can anyone sort out what is my mistake? I attach the input.

input ------------------------------------------

units lj
atom_style ellipsoid
dimension 3
boundary p p f

lattice fcc 1.88
region simbox2 block 0 4 0 4 8 12
region simbox3 block 0 20 0 20 0 20
create_box 1 simbox2
create_atoms 1 box
fix zwalls all wall/lj93 zlo 2.0 1.0 1.0 2.5 zhi 18.0 1.0 1.0 2.5

#set group all type/fraction 2 0.1 73263
#set type 1 mass 1.0
set type 1 shape 1 1 0.5
set group all quat/random 27193

compute rot all temp/asphere
group spheroid type 1
variable dof equal count(spheroid)
compute_modify rot extra ${dof}
compute kerot all erotate/asphere
velocity all create 1.607 348943 loop geom

pair_style gayberne 1.0 1.0 2.0 4.0
pair_coeff * * 1.0 1.0 1 1 1 1 1 1

neighbor 0.3 bin
neigh_modify every 10 delay 0 check no

restart 2500 restart # punto de salvado

thermo_style custom step c_rot temp ke c_kerot pe etotal press vol
thermo 10
compute orient all property/atom quati quatj quatk quatw
compute shape all property/atom shapex shapey shapez

dump 1 all custom 10 dump.lammpatrj id type x y z vx vy vz &
c_orient[1] c_orient[2] c_orient[3] c_orient[4] &
c_shape[1] c_shape[2] c_shape[3]

timestep 0.0002

fix 1 all nvt/asphere temp 0.5 0.5 0.01
compute_modify 1_temp extra ${dof}
#fix 1 all nve/asphere
run 1000

write_restart restart.end


Lost atoms are usually a consequence of high potential energy, typically due to random placement of (point) particles or due to random orientations.

The way to deal with it is to “drain” the excess potential energy from the system (that usually gets converted into kinetic energy). fix nvt/asphere is not quite as good at that as fix langevin or other dissipative thermostats. the key is to have a sufficiently short simulation time step and thermostat relaxation time. after the system can maintain a stable kinetic energy, you can choose less restrictive settings and/or switch back to nose-hoover thermostats.

Please also note that the fix wall/lj93 interaction only “knows” point particles and thus will only interact with your ellipsoids based on their center. The contributed fix wall/ees style is the only wall style in the LAMMPS distribution that I know to be aware of ellipsoids.