# direct shear v1, setup geometry in .svg file # Simulation procedure: To simulation our experiment procedure, we first randomly and separately pour particles above the shear cell, and let the kinetic energy goes to zero. This by experience 10 secs is enough. Then, we shear the packing by displacing the top tube at speed 0.1 mm/s, upto 1/3 (maxTvl) of the tube diameter. Then the simulation is relaxed for an additional sec and stopped. 1 sec = 1e6 simulation step ############################################### ### Tunable Physical parameters # cell radius variable cellR equal 480 # number of particle variable nPrt equal 250 # max travel variable maxTvl equal 300.0 # simulation seed variable seedR equal 12 # particle and particle-wall interaction variable kn equal 1.2e6 variable kt equal 1.3e6 variable damp equal 1.5e4 variable mu equal 0.1 variable muW equal 0.36 ################################################ ### setup parameters # box x reach variable boxL equal ${cellR}+${maxTvl} # define and realize a box region in 3d, the same unit as particles region reg block -${cellR} ${boxL} -${cellR} ${cellR} 0 2535 units box # shear time (unit timestep) variable shearN equal ${maxTvl}*1000000 # total simulation step variable totalRun equal 10000000+${shearN}+1000000 # define cylinder motion position variable t1 equal step>=10000000 variable t2 equal step>(10000000+v_shearN) variable r equal (step-10000000)/v_shearN variable dispX equal (v_r*v_t1+(1-v_r)*v_t2)*v_maxTvl # bottom tube: cylindrical container with top open region wall1 cylinder z 0 0 ${cellR} 0 830 open 2 units box # top tube: cylindrical container which will be movable region wall4 cylinder z 0 0 ${cellR} 835 1835 open 1 open 2 units box & move v_dispX NULL NULL # pourin region: insure region size + molecule size does not overlap wall variable pourR equal ${cellR}-150 region pourinregion cylinder z 0 0 ${pourR} 1835 2535 units box ############################################### ### set up simulation # finite size particle: use sphere atom_style sphere # map for molecular systems atom_modify map array # fixed, fixed, fixed-and-shrink-wrapped with a minimum value boundary f f fm # Newton's 3rd law is not communicated between processors newton on # ghost particles communication between processors of lammps # 1. velocity is communicated for compute granular pairwise interaction # 2. look into neigbouring proccessors for ghost atoms with distance cutoff comm_modify vel yes cutoff 2.0 # initiate simulation box create_box 1 reg # modify neighbor list rebuild skin distance neighbor 0.2 bin # modify neighbor list rebuild time step neigh_modify delay 0 # specify granular interactions pair_style gran/hertz/history ${kn} ${kt} ${damp} ${damp} ${mu} 1 pair_coeff * * # add time step: 1e-6 timestep 0.000001 # define molecule IDs for fix rigid, extended to ghost atoms as well fix prop all property/atom mol ghost yes # define hexapods molecule hexapod ../simulationTemplate/star30v3.hexapod fix hexarigid all rigid/small molecule mol hexapod # exclude inter molecular interaction neigh_modify exclude molecule/intra all ############################################### ### add boundary conditions and pour particles # add gravity: g = 9.8e4 fix 2 all gravity 9.8e4 spherical 0.0 -180.0 # add bottom tube fix w1 all wall/gran/region hertz/history ${kn} ${kt} ${damp} ${damp} ${muW} 1 region wall1 # add top tube fix w4 all wall/gran/region hertz/history ${kn} ${kt} ${damp} ${damp} ${muW} 1 region wall4 # Insert finite-size particles or molecules fix ins all pour ${nPrt} 0 ${seedR} id next vol 0.4 20 region pourinregion mol hexapod rigid hexarigid ###################################### ### output # compute inter-atomic forces compute 3 all property/local patom1 patom2 compute 2 all pair/local dist fx fy fz p1 p2 p3 # compute translational rigid body ke compute keT all ke/rigid hexarigid # compute rigid body net torque compute 4 all rigid/local hexarigid mol tqx tqy tqz # specify thermodynamic output thermo_style custom step atoms c_keT v_dispX fmax # do not normalize extensive quantities such as kinetic energy thermo_modify lost warn norm no # print thermodynamic information every 1000000 step thermo 1000000 # output images variable colors string "red green blue yellow white purple pink orange lime gray" variable mol2 atom mol%10 dump 2 all image 1000000 images/image.*.jpg v_mol2 type zoom 2.0 adiam 30.0 dump_modify 2 pad 5 amap 0 10 sa 1 10 ${colors} # output atomic coordination and forces dump id all custom 1000000 atomInfo.dump id mol x y z fx fy fz dump_modify id format float "%.15g" # output interparticle forces dump 3 all local 1000000 contactForce.dump c_3[*] c_2[*] dump_modify 3 format line "%0.0f %0.0f %.15g %.15g %.15g %.15g %.15g %.15g %.15g" # output rigid body forces dump 4 all local 1000000 bodyTorque.dump c_4[*] dump_modify 4 format line "%0.0f %.15g %.15g %.15g" # additional simulation run ${totalRun}