New question #631016 on Yade:
https://answers.launchpad.net/yade/+question/631016
I want to make two-ball clumps from spheres. Would appreciate some help.
Thanks.
#rMean is particles size
utils.readParamsFromTable(rMean=.053,rRelFuzz=.3,maxLoad=1e6,minLoad=1e4)
from yade.params.table import *
from yade import export,pack, plot
#cylinder
O.materials.append(FrictMat(density=2800,young=3e10,poisson=.4,frictionAngle=radians(22),label='walls'))
O.bodies.append(utils.geom.facetCylinder((.5,.5,.75),0.6,1.5,orientation=Quaternion((1,0,0),0),segmentsNumber=16,wallMask=4,material='walls'))
O.bodies.append(utils.geom.facetCylinder((.5,.5,.75),0.6,1.5,orientation=Quaternion((1,0,0),0),segmentsNumber=16,wallMask=2,material='walls'))
#spheres
O.materials.append(FrictMat(density=2650,young=5e8,poisson=.4,frictionAngle=radians(33),label='spheres'))
sp=pack.SpherePack()
sp.makeCloud((0.15,0.15,0.05),(0.9,0.9,3.6),rMean=rMean,rRelFuzz=.0,seed=1)
sp.toSimulation(material='spheres')
O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_L3Geom(),Ig2_Facet_Sphere_L3Geom(),Ig2_Wall_Sphere_L3Geom()],
[Ip2_FrictMat_FrictMat_FrictPhys()],
[Law2_L3Geom_FrictPhys_ElPerfPl()]
),
GravityEngine(gravity=(0,0,-9.81)),
NewtonIntegrator(damping=0.2),
# the label creates an automatic variable referring to this engine
# we use it below to change its attributes from the functions called
PyRunner(command='checkUnbalanced()',realPeriod=2,label='checker'),
]
O.dt=.5*utils.PWaveTimeStep()
def checkUnbalanced():
# at the very start, unbalanced force can be low as there is only few
contacts, but it does not mean the packing is stable
if O.iter<1: return
# the rest will be run only if unbalanced is < .1 (stabilized packing)
if utils.unbalancedForce()>.1: return
# add plate at the position on the top of the packing
# the maximum finds the z-coordinate of the top of the topmost particle
#O.bodies.append(utils.wall(max([b.state.pos[2]+b.shape.radius for b in
O.bodies if isinstance(b.shape,Sphere)]),axis=2,sense=-1))
#m=max([b.state.pos[2]+b.shape.radius for b in O.bodies if
isinstance(b.shape,Sphere)]
b = len(O.bodies)
#print "number of bodies",b
#penetrometer
O.materials.append(FrictMat(density=2800,young=3e10,poisson=.4,frictionAngle=radians(0),label='shaft'))
O.bodies.append(utils.geom.facetCylinder((.5,.5,2.5),0.1,2.0,orientation=Quaternion((1,0,0),0),segmentsNumber=6,wallMask=4,material='shaft'))
O.bodies.append(utils.geom.facetCylinder((.5,.5,2.5),0.1,2.0,orientation=Quaternion((1,0,0),0),segmentsNumber=6,wallMask=1,material='shaft'))
O.bodies.append(utils.geom.facetCone((.5,.5,1.4),0.1,0.0,0.2,orientation=Quaternion((1,0,0),0),segmentsNumber=6,wallMask=4,material='walls'))
#b = len(O.bodies)
#print "number of bodies2",b
#O.pause()
global
plate1,plate2,plate3,plate4,plate5,plate6,plate7,plate8,plate9,plate10,plate11,plate12,plate13,plate14,plate15,plate16,plate17,plate18,plate19,plate20,plate21,plate22,plate23,plate24
plate1=O.bodies[b] # the last body is the plate
plate2=O.bodies[b+1]
plate3=O.bodies[b+2]
plate4=O.bodies[b+3]
plate5=O.bodies[b+4]
plate6=O.bodies[b+5]
plate7=O.bodies[b+6]
plate8=O.bodies[b+7]
plate9=O.bodies[b+8]
plate10=O.bodies[b+9]
plate11=O.bodies[b+10]
plate12=O.bodies[b+11]
plate13=O.bodies[b+12]
plate14=O.bodies[b+13]
plate15=O.bodies[b+14]
plate16=O.bodies[b+15]
plate17=O.bodies[b+16]
plate18=O.bodies[b+17]
plate19=O.bodies[b+18]
plate20=O.bodies[b+19]
plate21=O.bodies[b+20]
plate22=O.bodies[b+21]
plate23=O.bodies[b+22]
plate24=O.bodies[b+23]
# Wall objects are "fixed" by default, i.e. not subject to forces
# prescribing a velocity will therefore make it move at constant velocity
(downwards)
plate1.state.vel=(0,0,-.2)
plate2.state.vel=(0,0,-.2)
plate3.state.vel=(0,0,-.2)
plate4.state.vel=(0,0,-.2)
plate5.state.vel=(0,0,-.2)
plate6.state.vel=(0,0,-.2)
plate7.state.vel=(0,0,-.2)
plate8.state.vel=(0,0,-.2)
plate9.state.vel=(0,0,-.2)
plate10.state.vel=(0,0,-.2)
plate11.state.vel=(0,0,-.2)
plate12.state.vel=(0,0,-.2)
plate13.state.vel=(0,0,-.2)
plate14.state.vel=(0,0,-.2)
plate15.state.vel=(0,0,-.2)
plate16.state.vel=(0,0,-.2)
plate17.state.vel=(0,0,-.2)
plate18.state.vel=(0,0,-.2)
plate19.state.vel=(0,0,-.2)
plate20.state.vel=(0,0,-.2)
plate21.state.vel=(0,0,-.2)
plate22.state.vel=(0,0,-.2)
plate23.state.vel=(0,0,-.2)
plate24.state.vel=(0,0,-.2)
# start plotting the data now, it was not interesting before
O.engines=O.engines+[PyRunner(command='addPlotData()',iterPeriod=200)]
#next time, do not call this function anymore, but the next one
(unloadPlate) instead
#checker.command='unloadPlate()'
checker.command='stopUnloading()'
#def
