particle crossing STL walls when 'fix freeze' is used

Hello Nico

Thanks for your answer. here is my run script-

echo both
log ../DEM/log.liggghts
thermo_log ../DEM/post/thermo.txt

atom_style granular
atom_modify map array
communicate single vel yes

boundary p f p
newton off

units si
processors 14 4 2

# read the restart file
read_restart ../DEM/post/restart/liggghts.restart

neighbor 0.0005 bin
neigh_modify delay 0

# Material properties required for granular pair styles

fix m1 all property/global youngsModulus peratomtype 5.e6 5.e6
fix m2 all property/global poissonsRatio peratomtype 0.45 0.45
fix m3 all property/global coefficientRestitution peratomtypepair 2 0.3 0.3 0.3 0.3
fix m4 all property/global coefficientFriction peratomtypepair 2 0.5 0.5 0.5 0

# pair style
pair_style gran model hertz tangential history # Hertzian without cohesion
pair_coeff * *

# timestep, gravity
timestep 0.0000001
fix ts_check all check/timestep/gran 1000 0.1 0.1
fix gravi all gravity 9.81 vector 0.0 -1.0 0.0

# Defining the walls that can be in contact with the particles; those come from the STL folder
fix inlet all mesh/surface file ../DEM/STL/inlet.stl type 1 curvature_tolerant yes
fix outlet all mesh/surface file ../DEM/STL/outlet.stl type 1 curvature_tolerant yes
fix top all mesh/surface file ../DEM/STL/top.stl type 1 curvature_tolerant yes
fix channelBottom_1 all mesh/surface/surface file ../DEM/STL/channelBottom_1.stl type 1 curvature_tolerant yes
fix channelBottom_2 all mesh/surface/surface file ../DEM/STL/channelBottom_2.stl type 1 curvature_tolerant yes
fix bedBottom all mesh/surface/surface file ../DEM/STL/bedBottom.stl type 1 curvature_tolerant yes
fix front all mesh/surface file ../DEM/STL/front.stl type 1 curvature_tolerant yes
fix back all mesh/surface file ../DEM/STL/back.stl type 1 curvature_tolerant yes

# The wall contacts are only allowed at the box walls
fix wall all wall/gran model hertz tangential history mesh n_meshes 6 meshes top channelBottom_1 channelBottom_2 bedBottom front back

#creating group for bed paricle (type 1)
group bed type 1

# particle distributions and insertion, IFS1 (bridging case)
fix pts1 all particletemplate/sphere 100631591 atom_type 2 density constant 2650 radius constant 0.00078
fix pts2 all particletemplate/sphere 100631807 atom_type 2 density constant 2650 radius constant 0.00049
fix pts3 all particletemplate/sphere 100661111 atom_type 2 density constant 2650 radius constant 0.00042
fix pts4 all particletemplate/sphere 100661399 atom_type 2 density constant 2650 radius constant 0.00036
fix pts5 all particletemplate/sphere 100793137 atom_type 2 density constant 2650 radius constant 0.00030
fix pts6 all particletemplate/sphere 100794431 atom_type 2 density constant 2650 radius constant 0.00025
fix pts7 all particletemplate/sphere 100796713 atom_type 2 density constant 2650 radius constant 0.00021
fix pts8 all particletemplate/sphere 100798699 atom_type 2 density constant 2650 radius constant 0.00018
fix pts9 all particletemplate/sphere 100799183 atom_type 2 density constant 2650 radius constant 0.00015
fix pts10 all particletemplate/sphere 100802419 atom_type 2 density constant 2650 radius constant 0.00011
fix pts11 all particletemplate/sphere 100805377 atom_type 2 density constant 2650 radius constant 0.00006 volume_limit 1e-14
fix pts12 all particletemplate/sphere 100822807 atom_type 2 density constant 2650 radius constant 0.00004 volume_limit 1e-14
fix pdd1 all particledistribution/discrete 100831279 12 pts1 8.8 pts2 4.6 pts3 4.7 pts4 8.2 pts5 9.2 pts6 11.9 pts7 16 pts8 15 pts9 15.2 pts10 5 pts11 1.2 pts12 0.2

#particle insertion (particles at the rate of 0.02 kg/s for 1 sec=0.02 kg particle should be inserted)-for 0.15 m wide flume
fix ins all insert/stream seed 101235583 distributiontemplate pdd1 mass 0.02 massrate 0.02 vel constant 0.675 0.0 0.0 overlapcheck yes insertion_face inlet extrude_length 0.1

# cfd coupling
fix cfd all couple/cfd couple_every 100 mpi
fix cfd2 all couple/cfd/force/implicit
#fix cfd2 all couple/cfd/force/accumulator RongDrag 10 1.5e-3
#fix cfd2 all couple/cfd/force/implicit/accumulated #CrankNicolson 0.5

# apply nve integration to all particles that are inserted as single particles
fix integr all nve/sphere

#freezing the bed particle
fix 2 bed freeze

# center of mass
compute centerOfMass all com

# compute explicit dragforce
compute explDrag all reduce update_on_run_end yes sum f_dragforce[1] f_dragforce[2] f_dragforce[3]

# sum of explicit and implicit drag force given from CFD to DEM
variable totalDragX equal f_cfd2[1]
variable totalDragY equal f_cfd2[2]
variable totalDragZ equal f_cfd2[3]

# explicit drag force given from CFD to DEM
variable explicitDragX equal c_explDrag[1]
variable explicitDragY equal c_explDrag[2]
variable explicitDragZ equal c_explDrag[3]

variable time equal step*dt
fix extra all print 10 "${time} ${explicitDragX} ${explicitDragY} ${explicitDragZ} ${totalDragX} ${totalDragY} ${totalDragZ}" file ../DEM/post/forces.txt title "# time expDrag(X Y Z) totDrag(X Y Z)"

# screen output
compute rke all erotate/sphere
thermo_style custom step atoms ke c_rke vol c_centerOfMass[3] c_explDrag[1] c_explDrag[2] c_explDrag[3] f_cfd2[1] f_cfd2[2] f_cfd2[3]
thermo 10
thermo_modify lost ignore norm no
compute_modify thermo_temp dynamic yes

#dump dmp all custom 5000 ../DEM/post/dump*.liggghts_run id type x y z vx vy vz fx fy fz f_dragforce[1] f_dragforce[2] f_dragforce[3] radius

dump dmp all custom/vtk 500000 ../DEM/post/dump_liggghts_run.*.vtk id type x y z vx vy vz fx fy fz f_dragforce[1] f_dragforce[2] f_dragforce[3] radius

reset_timestep 0
run 0

Another note: I also tried by putting 'fix freeze' command before the 'fix nve/sphare' command. But same results are produced where some large particles (type 1) are thrown out of stl boundary.

Best Regards
Atul

I am still struggling to figure out, what is causing the particles to cross STL walls when 'fix freeze' command is used. Any help will be highly appriciated.

Hello Mschramm

when looking into the log file, I am sure that there are 11258 atoms in group bed. The type-1 atoms are created in initilzation script (in.liggghts_init) and restart file from initialization run is used for cdfem simulation in run script (in.liggghts_run). So all the particles (type-1) from previous simulation should already initialized in run script. As you can see from above attached run script that the type-1 particles are gruopped togetehr as bed and then bed is freezed.

As it is clear that group bed consists of 11258 number of type-1 atoms and indeed there are particle in the created group, this problem of particle crossing STL walls is related top something else.

do you have any other comments/suggestions to solve the issue?

Best Regards
Atul

After struggling for more than a month, I managed to find the cause of it. I am not sure if it because of bug in LIGGGHTS code but here is what I have done and that solved the issue.

Before freezing the bed particles, I force the bed particles to have zero velocity and then freeze them. So in the above script, one line is added before freeze command:

velocity bed set 0.0 0.0 0.0 units box #additional line
fix 2 bed freeze

I hope this would help other facing similar situation and they dont have to spend extra time.

Best Regards
Atul Jaiswal