## DEM Shear Test of wheat straw in LIGGGHTS ### Material Independent Properties variable particle_diameter equal 2.50e-3 # particles diameter in meters variable bond_out_diameter equal 2.50e-3 # fiber outer diameter in meters variable bond_in_diameter equal 0.0 # fiber inner diameter in meters variable bond_length equal 2.50e-3 # distance, in meters, between two particles in bond variable bond_damp_type equal 1 # Use Yu Guo Damping variable bond_damp_val equal 10 # NA variable particle_density equal 2500.0 # kg/m3 variable fiber_contact_youngs_modulus equal 1.0e7 # Pa variable wall_contact_youngs_modulus equal 180.0e9 # Pa variable bond_youngs_modulus equal 1.0e7 # Pa variable particle_poissons_ratio equal 0.3 # NA variable wall_poissons_ratio equal 0.3 # NA variable coef_res_pp equal 0.0600 # particle-particle coefficient of restitution variable coef_res_pw equal 0.0600 # particle-wall coefficient of restitution variable coef_res_ww equal 0.500 # wall-wall coefficient of restitution variable coef_fri_pp equal 0.400 # particle-particle coefficient of friction variable coef_fri_pw equal 0.600 # particle-wall coefficient of friction variable coef_fri_ww equal 0.200 # wall-wall coefficient of friction ### Simulation Independent Parameters variable ke_tol equal 1.0e-4 # Energy that we will run the simulation to obtain variable res_tol equal 1.0e-6 # Change of height tolerance that must be met variable shear_pressure equal 50000 # Pressure in N/m2 variable shear_box_length equal 0.101 # Length of shear box in meters variable mass_shear_top_plate equal 0.401 # Mass of the top shear plate variable restart_time equal 1.0e-2 # how often we write a restart file variable fileprint_time equal 1.0e-2 # how often we print to the file in seconds variable thermo_time equal 1.0e-3 # how often we print to the screen in seconds variable output_time equal 1.0e-3 # how often data is written to csv file in seconds ### Material Dependent Properties variable particle_radius equal 0.5*${particle_diameter} variable bond_shear_modulus equal ${bond_youngs_modulus}/(2.0*(1.0+${particle_poissons_ratio})) variable bond_out_per equal ${bond_out_diameter}/${particle_diameter} variable bond_in_per equal ${bond_in_diameter}/${particle_diameter} ### Calculate dt using the bond model variable r2 equal ${particle_radius}*${particle_radius} variable r3 equal ${r2}*${particle_radius} variable K equal ${bond_youngs_modulus}*PI*${r2}/${bond_length} variable m equal 4.0*PI*${r3}*${particle_density}/3.0 variable w equal sqrt($K/$m) variable dt equal 0.8/((1.0+3.0*${bond_damp_val})*$w) ### Simulation Dependent Parameters variable restart_step equal ceil(${restart_time}/${dt}) variable fileprint_step equal ceil(${fileprint_time}/${dt}) variable thermo_step equal ceil(${thermo_time}/${dt}) variable output_step equal ceil(${output_time}/${dt}) variable shear_force equal ${shear_pressure}*(${shear_box_length}*${shear_box_length})+${mass_shear_top_plate}*9.81 ## Build Simulation # Set Atom Style atom_style hybrid granular bond/gran n_bondtypes 2 bonds_per_atom 50 atom_modify map array hard_particles yes # Fix boundaries boundary f f f newton off communicate single vel yes units si # Create simulation domain and tell liggghts that you have 2 types of particles region domain block -0.125 0.225 -0.125 0.125 -0.01 0.125 units box create_box 2 domain # Use hertz-mindlin and granular fiber bonds pair_style gran model hertz tangential history bond_style gran # Create a folder for restart files shell mkdir restarts # Set bin size and to check for a neighbor update every step neighbor 0.001 bin neigh_modify delay 0 # Set coefficients for contact and bond models pair_coeff * * # none are needed bond_coeff 1 ${bond_out_per} ${bond_in_per} ${bond_youngs_modulus} ${bond_shear_modulus} ${bond_damp_type} ${bond_damp_val} 1 1.0e32 1.0e32 bond_coeff 2 ${bond_out_per} ${bond_in_per} ${bond_youngs_modulus} ${bond_shear_modulus} ${bond_damp_type} ${bond_damp_val} 1 1.0e32 1.0e32 # Set material properties fix m1 all property/global youngsModulus peratomtype ${fiber_contact_youngs_modulus} ${wall_contact_youngs_modulus} fix m2 all property/global poissonsRatio peratomtype ${particle_poissons_ratio} ${wall_poissons_ratio} fix m3 all property/global coefficientRestitution peratomtypepair 2 ${coef_res_pp} ${coef_res_pw} & ${coef_res_pw} ${coef_res_ww} fix m4 all property/global coefficientFriction peratomtypepair 2 ${coef_fri_pp} ${coef_fri_pw} & ${coef_fri_pw} ${coef_fri_ww} # Import the geometry files and rotate/move to correct locations fix BotShear all mesh/surface/stress file STL/ShearBox_TopPlate.STL type 2 & scale 0.001 move 0.0 0.0 0.0 curvature_tolerant yes # Top plate has a mass of 401 grams fix TopShear all mesh/surface/stress/servo file STL/ShearBox_TopPlate.STL type 2 & scale 0.001 rotate axis 1 0 0 angle 180 move 0.0 0.1015 0.080 & curvature_tolerant yes com 0.05 0.05 0.08 ctrlPV force & axis 0.0 0.0 -1.0 target_val ${shear_force} vel_max 5.0 & kp 5.0e-2 ki 0.0 kd 0.0 fix BotHalfShear all mesh/surface/stress file STL/ShearBox_BottomHalf.STL & type 2 scale 0.001 move 0.0 0.0 0.0 curvature_tolerant yes fix TopHalfShear all mesh/surface/stress file STL/ShearBox_TopHalf.STL & type 2 scale 0.001 move 0.0 0.0 0.025 curvature_tolerant yes # Make all of the geometries able to interact with particles fix wall all wall/gran model hertz tangential history mesh n_meshes 4 meshes & BotShear TopShear BotHalfShear TopHalfShear # Add primitive walls for particle insertion fix wall_x1 all wall/gran model hertz tangential history primitive type 2 xplane 0.0 fix wall_x2 all wall/gran model hertz tangential history primitive type 2 xplane 0.1015 fix wall_y1 all wall/gran model hertz tangential history primitive type 2 yplane 0.0 fix wall_y2 all wall/gran model hertz tangential history primitive type 2 yplane 0.1015 # Add gravity fix grav all gravity 9.81 vector 0.0 0.0 -1.0 # Integrate using velocity verlet fix integr all nve/sphere # Create our fibers fix pts1 all particletemplate/multiplespheres 15485863 atom_type 1 & density constant ${particle_density} nspheres 4 ntry 5000 spheres & 0.0000 0.0 0.0 0.00125 & 0.0025 0.0 0.0 0.00125 & 0.0050 0.0 0.0 0.00125 & 0.0075 0.0 0.0 0.00125 & bonded yes # Set the distribution of particles fix pdd1 all particledistribution/discrete 32452843 1 pts1 1.0 # Create a region to insert the particles region fill_box block 0.0 0.1015 0.0 0.1015 0.0065 0.0735 units box # Particle insertion fix ins all insert/pack seed 32452867 distributiontemplate pdd1 & maxattempt 500 insert_every once overlapcheck yes orientation random & all_in yes vel constant 0.0 0.0 -0.0 region fill_box & particles_in_region 5000 ntry_mc 10000 check_dist_from_subdomain_border no # Set the timestep to use timestep ${dt} # Set skin for for bond creation variable bond_skin equal 1.000001*${particle_diameter} # Bond the atoms together to form fibers run 1 fix bondcr all bond/create/gran 1 1 1 ${bond_skin} 1 6 #every itype jtype cutoff btype newperts run 1 fix_modify bondcr every 0 #do not create new bonds after this line run 1 # Write the restart file write_restart restarts/restart0.liggghts