###################################Header for General commands#######################################################

atom_style	granular		#Simulation of particles

boundary	f f f			#fixed boundaries -> particles will be deleted if leaving the simulation box

units		si

communicate	single vel yes	#default

newton		off				#default


###################################System variables##################################################################

#Definition of boundaries
variable xmin equal -2.55
variable xmax equal 2.55

variable ymin equal -2.55
variable ymax equal 2.55

variable zmin equal -0.05
variable zmax equal 5.0

#Definition of the timestep

variable dt equal 1e-4			#timestep = 0.0001 second; Each iteration step represents 0.0001 seconds.


###################################Specific variables for current simulation#########################################

variable natoms equal 2 	#1 -> particle #2-> substrate, recoater, hopper

####variable for material properties####f

####Young Modulus####
variable	youngmodulus1	equal 	1e8		#N/mm�
variable	youngmodulus2	equal 	1e8		#N/mm�

####Poission ratio####
variable 	poission1	equal	0.3
variable 	poission2	equal	0.3

####variable for contact properties####

####coefficient of restitution####
variable	CoR11	equal	0.6
variable	CoR12	equal	0.6	
variable	CoR21	equal	0.6
variable	CoR22	equal	0.6

####sliding friction coefficient####
variable	sf11	equal	0.3
variable	sf12	equal	0.36
variable	sf21	equal	0.36
variable	sf22	equal	0.0

####rolling friction coefficient####
variable	rf11	equal	0.8
variable	rf12	equal	0.8
variable	rf21	equal	0.8
variable	rf22	equal	0.8



####variable for particle#####

#Number of particle radius

variable	nradii	equal	8

variable	radius1	equal	0.0075	#m
variable	radius2	equal	0.008	#m
variable	radius3	equal	0.009	#m
variable	radius4	equal	0.010	#m
variable	radius5	equal	0.011	#m
variable	radius6	equal	0.014	#m
variable	radius7	equal	0.018	#m
variable	radius8	equal	0.0225	#m

variable	frac1	equal	0.07	#7%
variable	frac2	equal	0.13
variable	frac3	equal	0.18
variable	frac4	equal	0.25
variable	frac5	equal	0.17
variable	frac6	equal	0.11
variable	frac7	equal	0.06
variable	frac8	equal	0.03

variable	density	equal	1000	#kg/m�





####filling parameters####

variable	filltime	equal	4		#seconds

variable	fillmass	equal	5	#kg

variable	fillmassrate	equal	${fillmass}/${filltime}	#kg/s

variable	fillsteps	equal	${filltime}/${dt}		#Transform time to iteration steps




####settle time####

variable	settletime	equal	3 #second

variable	settlesteps	equal	${settletime}/${dt}		#Transform time to iteration steps




#####recoater blade 1 #####

variable	recoatervel1	equal	-0.4	#m/s

variable	recoatertime1	equal	0.5		#seconds

variable	recoatersteps1	equal	${recoatertime1}/${dt}	#transform time to iteration steps





####recoater blade 2 #####

variable	recoatervel2	equal	0.1		#m/s

variable	recoatertime2	equal	14	#seconds

variable	recoatersteps2	equal	${recoatertime2}/${dt}	#transform time to iteration steps






####discharge time 2 #####

variable	dischargetime2	equal	1	#seconds

variable	dischargesteps2	equal	${dischargetime2}/${dt}	#Transform time to iteration steps



###################################Definition of simulationbox#######################################################

region		reg block ${xmin} ${xmax} ${ymin} ${ymax} ${zmin} ${zmax} units box

create_box	2 reg

neighbor	0.004 bin		#default

neigh_modify	delay 0		#default


#################################Definition of the contact models####################################################

pair_style gran model hertz tangential history rolling_friction epsd2	#contact model

pair_coeff	* *				#default

timestep ${dt}

fix		integrator all nve/sphere			#default

fix		gravi all gravity 9.81 vector 0.0 0.0 -1.0		#gravity of 9.81 m/s� in negative z direction




###################################Definition of Material properties#################################################


fix 		m1 all property/global youngsModulus peratomtype ${youngmodulus1} ${youngmodulus2}

fix 		m2 all property/global poissonsRatio peratomtype ${poission1} ${poission2}

fix 		m3 all property/global coefficientRestitution peratomtypepair ${natoms} 0.6 0.6 0.6 0.6

fix 		m4 all property/global coefficientFriction peratomtypepair ${natoms} ${sf11} ${sf12} ${sf21} ${sf22}

fix 		m5 all property/global coefficientRollingFriction peratomtypepair ${natoms} ${rf11} ${rf12} ${rf21} ${rf22}











###################################Generation and Loading of the Geometry .stl#######################################

fix			substrate	all	mesh/surface file substrate.stl		type 2  scale 0.001		#load mesh from STL file. Type 2 for geometry. Scale down to transform mm to meters

fix			recoater	all	mesh/surface file recoater.stl		type 2  scale 0.001



fix 	walls 	all wall/gran model hertz tangential history rolling_friction epsd2 mesh n_meshes 2 meshes substrate recoater






###################################Generation and Insertion of the particles#########################################




fix	pts1 all particletemplate/sphere 10487 atom_type 1 density constant ${density} radius constant ${radius1}
fix 	pts2 all particletemplate/sphere 11887 atom_type 1 density constant ${density} radius constant ${radius2}
fix 	pts3 all particletemplate/sphere 11897 atom_type 1 density constant ${density} radius constant ${radius3}
fix 	pts4 all particletemplate/sphere 11903 atom_type 1 density constant ${density} radius constant ${radius4}
fix 	pts5 all particletemplate/sphere 11909 atom_type 1 density constant ${density} radius constant ${radius5}
fix 	pts6 all particletemplate/sphere 11923 atom_type 1 density constant ${density} radius constant ${radius6}
fix 	pts7 all particletemplate/sphere 12097 atom_type 1 density constant ${density} radius constant ${radius7}
fix	pts8 all particletemplate/sphere 12197 atom_type 1 density constant ${density} radius constant ${radius8}
fix 	pdd1 all particledistribution/discrete 32452867 ${nradii} pts1 ${frac1} pts2 ${frac2} pts3 ${frac3} pts4 ${frac4} pts5 ${frac5} pts6 ${frac6} pts7 ${frac7} pts8 ${frac8}
fix 	ins_mesh all mesh/surface/planar file Insertionface.stl type 1 scale 0.001
fix 	ins all insert/stream seed 86028157 distributiontemplate pdd1 &
	mass ${fillmass} massrate ${fillmassrate} overlapcheck yes all_in yes vel constant 0 0 -0.75 &
	insertion_face ins_mesh extrude_length 0.25



###################################Dumping of the data for post-processing to visualize############################

shell 	mkdir post

#Definition of the dumptime

variable dumptime equal 0.04	# Every 0.04 s 1 image 

variable dumpstep equal ${dumptime}/${dt}	#Transform to iteration steps



dump    dmpparticle     all custom/vtk ${dumpstep}  post/particles_*.vtk	id type x y z vx vy vz fx fy fz radius mass
dump    dmpsubstrate	all mesh/stl ${dumpstep}    post/Substrate*.stl		substrate
dump    dmprecoater     all mesh/stl ${dumpstep}    post/Recoater*.stl		recoater





####################################RUN the simulation filling###########################################################

run	 ${fillsteps}

unfix	ins

##################################RUN the simulation settling#############################################################

run	 ${settlesteps}