Why do particles end up on or inside fix wall surface? (Seemingly programming them not to)

Hello! I am having trouble understanding why I am receiving the error (ERROR on proc 0: Particle on or inside fix wall surface) in the code I am trying to run:

------------- Initialization

units lj
dimension 2
atom_style sphere
boundary f f p
pair_style lj/cut 33.67384588865670
pair_modify shift yes

------------- System definition

region testbox block -10 10 -10 10 -0.25 0.25
create_box 2 testbox
create_atoms 1 single 5.0 5.0 0
create_atoms 1 single 5.0 -5.0 0
create_atoms 2 single -5.0 -5.0 0
create_atoms 2 single -5.0 5.0 0
group 1 type 1
group 2 type 2
fix wallxhi all wall/lj93 xhi EDGE 1.0 1.0 0.1
fix wallxlo all wall/lj93 xlo EDGE 1.0 1.0 0.1
fix wallylo all wall/lj93 ylo EDGE 1.0 1.0 0.1
fix wallyhi all wall/lj93 yhi EDGE 1.0 1.0 0.1
fix 1 1 brownian 100 44830747 gamma_t 1.0 rng gaussian
fix 2 2 brownian 50000 34490088 gamma_t 1.0 rng gaussian
fix 3 1 viscous 0.0001884955592153876
fix 4 2 viscous 0.0001884955592153876
pair_coeff 1 1 0.00000000000431 30.00000000000000 33.67384588865670
pair_coeff 1 2 0.00000000107932 30.00000000000000 33.67384588865670
pair_coeff 2 2 0.00000000215433 30.00000000000000 33.67384588865670
dump mydmp all xyz 1 dump_(16_26_16).xyz
timestep 5e-06
run 100000

Output:

LAMMPS (29 Sep 2021 - Update 3)
Created orthogonal box = (-10.000000 -10.000000 -0.25000000) to (10.000000 10.000000 0.25000000)
1 by 1 by 1 MPI processor grid
Created 1 atoms
using lattice units in orthogonal box = (-10.000000 -10.000000 -0.25000000) to (10.000000 10.000000 0.25000000)
create_atoms CPU = 0.000 seconds
Created 1 atoms
using lattice units in orthogonal box = (-10.000000 -10.000000 -0.25000000) to (10.000000 10.000000 0.25000000)
create_atoms CPU = 0.000 seconds
Created 1 atoms
using lattice units in orthogonal box = (-10.000000 -10.000000 -0.25000000) to (10.000000 10.000000 0.25000000)
create_atoms CPU = 0.000 seconds
Created 1 atoms
using lattice units in orthogonal box = (-10.000000 -10.000000 -0.25000000) to (10.000000 10.000000 0.25000000)
create_atoms CPU = 0.000 seconds
2 atoms in group 1
2 atoms in group 2
Neighbor list info …
update every 1 steps, delay 10 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 33.973846
ghost atom cutoff = 33.973846
binsize = 16.986923, bins = 2 2 1
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair lj/cut, perpetual
attributes: half, newton on
pair build: half/bin/atomonly/newton
stencil: half/bin/2d
bin: standard
Setting up Verlet run …
Unit style : lj
Current step : 0
Time step : 5e-06
Per MPI rank memory allocation (min/avg/max) = 5.156 | 5.156 | 5.156 Mbytes
Step Temp E_pair E_mol TotEng Press
0 0 0.0023089621 0 0.0023089621 0.00013863795
ERROR on proc 0: Particle on or inside fix wall surface (…/fix_wall_lj93.cpp:95)
Last command: run 100000

You are using a cutoff of 0.1 sigma. You are not giving the wall a chance to repel any particle.

P.S.: why do you use 4 fix wall/lj93 instances? you could just have one for all walls.

Thank you for your response! I did not think you could write one instance for all walls, but thank you, I shall give that a try. As for the cutoff. I had used 1.0, 1.5 with the same result. What value would you recommend I use for the situation presented in my code?

I cannot make a recommendation because the closer I look at your input, the less sense it makes to me. You have reduced units and a 20 sigma by 20 sigma 2-d box but then your pair potential coefficient have a sigma parameter of 30 sigma and an even larger cutoff. since it is all the same LJ sigma, why not use 1.0? And why does it have to be larger than the box? Those particles will obviously repel but then the walls have a tiny sigma and tinier cutoff in comparison.

The other thing I don’t understand is the massively large temperature you chose in fix brownian. Why? This is what is slamming your particles into the walls with so much energy that they reach the boundary.

I have two more questions/remarks on this:

  • the concept of temperature does not make much sense for just two atoms. Why such a small system?
  • what is the purpose of using atom style sphere? I know it is required by fix brownian/sphere, but since your lj/cut pair style computes interactions solely based on distance, there is no meaning to having particles that have the ability of rotate. That rotation is only set to some random numbers.

In addition to Axel’s questions:

  • shouldn’t you use the fix enforce2d to ensure that atoms stay in their plane? (its a genuine question, I’ve never tried to run a 2D simulation without it)