include data.heead modify_timing off echo both variable basepath string "../../.." variable neighborDistance equal 2*${radiusP} ################################################ ### compute the Rayleigh timestep if it has not been set, yet. ##################################################################################### if "${SIMTS} == 0" then "variable SIMTS equal ${percentRayleigh}*(PI*${radiusP}*sqrt(${densityP}/(${youngsModulusP}/(2*(1+${poissonsRatioP}))))/((0.1631*${poissonsRatioP})+0.8766))" ###################################################################################### ### convert simulation durations into timesteps ###################################################################################### variable dump_time equal round(${DumpInterv}/${SIMTS}) ###################################################################################### ### This model simulates a slip-stick plate with constant input of particles. The ### particles are tracked by a massflow surface beneath the end of the plat. Massflow ### data is written into the responding folders. After the simulation a python script ### is called that calculates the mass at different time steps and the least sqaure ### difference to the experimental results at those time steps. ###################################################################################### ### initialization ###################################################################################### atom_style granular atom_modify map array boundary f f f newton off communicate single vel yes processors 8 1 1 units si region reg block -0.85 0.65 -0.08 0.08 -0.18 0.06 create_box 2 reg neighbor ${neighborDistance} bin neigh_modify delay 0 ###################################################################################### ### assign material properties to bulk solid and walls ###################################################################################### #fix fixBal all balance 2000 x 10 1.05 #Material properties required for new pair styles fix m1 all property/global youngsModulus peratomtype ${youngsModulusP} ${youngsModulusW} fix m2 all property/global poissonsRatio peratomtype ${poissonsRatioP} ${poissonsRatioW} fix m3 all property/global coefficientRestitution peratomtypepair 2 & ${coefRestitutionPP} ${coefRestitutionPW} & ${coefRestitutionPW} ${coefRestitutionWW} & fix m4 all property/global coefficientFriction peratomtypepair 2 & ${coefFrictionPP} ${coefFrictionPW} & ${coefFrictionPW} ${coefFrictionWW} & fix m5 all property/global coefficientRollingFriction peratomtypepair 2 & ${coefRollingFrictionPP} ${coefRollingFrictionPW} & ${coefRollingFrictionPW} ${coefRollingFrictionWW} & #New pair style pair_style gran model hertz tangential history rolling_friction epsd2 #Hertzian without cohesion pair_coeff * * timestep ${SIMTS} fix gravi all gravity 9.81 vector 0 0 -1 shell mkdir post massflow #granular walls #SlipstickplateMesh fix plateMesh all mesh/surface file meshes/Chute_modi.stl type 2 scale 0.001 fix move all move/mesh mesh plateMesh viblin axis -1. 0. 0 order 6 amplitude 6.3969629646e-03 4.5640435724e-03 2.4099411107e-03 1.4776328568e-03 8.8529316960e-05 4.9801942458e-04 phase -8.9355625141e-01 -4.0490753469e+00 1.5271174328e+00 -2.7216769702e+00 2.4800649251e+00 -1.5843410458e+00 period 2.4721455701e+02 2.1360336651e-01 1.3681424351e+02 1.0680168694e-01 5.8203380003e+01 7.1201139663e-02 fix conti all mesh/surface file meshes/Container_modi.stl type 2 scale 0.001 fix geometry all wall/gran model hertz tangential history rolling_friction epsd2 mesh n_meshes 2 meshes plateMesh conti #distributions for insertion fix pts1 all particletemplate/sphere 32452843 atom_type 1 density constant ${densityP} radius constant ${radiusP} #region for insertion #region insBox block 0.495 0.501645 -0.0033225 0.0033225 0.034 0.084 fix ins_mesh all mesh/surface file meshes/Insertionface_15x15mm.stl type 2 scale 0.001 #distribution fix insDistr all particledistribution/discrete 32452867 1 pts1 1 #particle insertion #fix partInsertion all insert/stream seed 49979687 distributiontemplate insDistr mass 0.550 vel constant 0 0 -0.0705 massrate 0.011093 insert_every 5000 overlapcheck yes insertion_face ins_mesh fix partInsertion all insert/stream seed 49979687 distributiontemplate insDistr mass 0.550 vel constant 0 0 -0.0705 massrate 0.011093 overlapcheck yes insertion_face ins_mesh extrude_length 0.02 #apply nve integration to all particles fix integr all nve/sphere region delBlock block -0.78 -0.063 -0.076 0.076 -0.112 -0.012 #output settings compute rke all erotate/sphere compute zmag all reduce/region delBlock max z thermo_style custom step time cpu atoms ke #c_rke vol c_zmag variable zMax equal c_zmag variable ThermoOut equal round(0.10/${SIMTS}) thermo ${ThermoOut} thermo_modify lost ignore norm no run 1 variable min_mass equal 0.230 variable mass_lead equal mass(all,delBlock) variable runsteps_mass230 equal round(76/${SIMTS}) #dump settings dump dmp all custom ${dump_time} post/plate*.liggghts id type type x y z ix iy iz vx vy vz fx fy fz omegax omegay omegaz radius dump dmpMesh all mesh/stl ${dump_time} post/Chutemove*.stl plateMesh dump_modify dmp region delBlock run 1 label loop3 print "---> Checking mass of particles ......,mass_check = ${mass_lead}" if "${mass_lead}>${min_mass}" then "jump in.Chute endloop3" print "---> mass_in_container = ${mass_lead}" run 20000 jump in.Chute loop3 label endloop3 variable runsteps equal round(60/${SIMTS}) #length of the experimental data 30 seconds to get 300 grams in container variable settletime equal round(2/${SIMTS} delete_atoms region delBlock # surface for massflowCounter fix inface1 all mesh/surface file meshes/Massflowsurafce_modi.stl type 2 scale 0.001 fix mass1 all massflow/mesh mesh inface1 vec_side 0 0 -1 screen no file massflow/massflow.dat writeTime count once delete_atoms no #fix move all move/mesh mesh plateMesh viblin axis 1. 0. 0 order 6 amplitude 6.3969629646e-03 4.5640435724e-03 2.4099411107e-03 1.4776328568e-03 8.8529316960e-05 4.9801942458e-04 phase -8.9355625141e-01 -4.0490753469e+00 1.5271174328e+00 -2.7216769702e+00 2.4800649251e+00 -1.5843410458e+00 period 2.4721455701e+02 2.1360336651e-01 1.3681424351e+02 1.0680168694e-01 5.8203380003e+01 7.1201139663e-02 #fix mark all property/atom/tracer region_mark delBlock mark_step 1 mark_style dirac #compute nparticles all nparticles/racer/region region_count delBlock tracer mark reset_marker no run 1 variable mass_threshold equal 0.299 variable n_particles equal count(all,delBlock) variable total_mass equal mass(all,delBlock) compute zmag all reduce/region update_on_run_end yes delBlock max z label loop2 print "---> Checking no. of particles ......,np = ${n_particles}" print "---> Checking mass of particles ......,mp = ${total_mass}" if "${total_mass}>${mass_threshold}" then "jump in.Chute endloop2" print "---> mp = ${total_mass}" run 5000 jump in.Chute loop2 label endloop2 unfix move unfix partInsertion run {settletime} variable xmax equal c_zma print "Height = " group grp_atoms region delBlock variable zMax equal bound(grp_atoms,zmax) print "Max height =${zMax}" print "Mass in block = ${total_mass}" shell python OctaveFuns/massflow_simLabExp.py