Question #691285 on Yade changed:
https://answers.launchpad.net/yade/+question/691285
NGANDU KALALA Gauthier posted a new comment:
Hi Jan,
Thank you for the response. This is a MWE i propose for a better understanding
of my problem.
As it is not possible to provide the mill mesh file in gmsh format, The mill
is simulated here as a cylinder of the same size, rotating around z-axis. The
balls are added to the mill using sphere pack, (p=pack.SpherePack());
However the behavior is not quite what I expected, especially since the mill is
not equipped with a lifter to make the balls circulate in the mill as it is
physically done. However there are interactions with the mill that normally
have to produce normal and tangential forces, That i can't save actually.
I also noticed that the increased rigidity, some balls pass through the walls
of the mill, instead of remaining inside.*
Moreover, I expected that at the beginning of the simulation, packages of balls
could fall along the y-axis, to disintegrate afterwards. But i notice that it
already disintegrates before falling, and go in the whole direction of mill,
despite the acceleration of gravity is defined according to the y-axis.*
*
Finaly find the same recording problem mentioned above.
This is the MWE.
# Material
rho = 7640.0 # [kg/m^3] density (of the balls)
Eyoung = 900.0e9 # [N/m] equivalent normal stiffness
vp = 0.6 # [N/m] equivalent tangential stiffness
en = 0.3 # [-] normal coefficient of restitution
(sphere/sphere)
mu = 2 # [-] friction coefficient
rMean = 0.8/2.0
length = 0.4
rBall = 0.0075 # [m] ball radius
sphRadFuzz = 0.8
# Simulation
pcv = 70.0 # [%] percentage of critical velocity
tEnd = 10.0 # [s] total simulation time
O.dt = 1.0e-5 # [s] fixed time step
tRec = 1.0e-2 # [s] recording interval
#material of cylinder
cylMat =
O.materials.append(FrictMat(poisson=vp,young=Eyoung,density=rho,frictionAngle=radians(mu)))
#material of spheres
ballMat =
O.materials.append(FrictMat(poisson=vp,young=Eyoung,density=rho,frictionAngle=radians(mu)))
ids=O.bodies.append(geom.facetCylinder((0,0,0),0.8,0.4,segmentsNumber=64,
wallMask=7, angleRange=None,
radiusTopInner=-1,
radiusBottomInner=-1,material=cylMat))
#walls = O.bodies.append([
# utils.wall((0,0,length/2),2,material=Mat1,color=(1,1,1),sense=+1),
# utils.wall((0,0,-length/2),2,material=Mat1,color=(1,1,1),sense=+1)
# ])
# define domains for initial cloud of red and blue spheres
packHt=.8*rMean # size of the area
bboxes=[(Vector3(-.5*length,-.5*packHt,-.5*packHt),Vector3(.5*length,0,.5*packHt)),(Vector3(-.5*length,-.5,-.5*packHt),Vector3(.5*length,.5*packHt,.5*packHt))]
colors=(1,0,0),(0,0,1)
for i in (0,1): # red and blue spheres
sp=pack.SpherePack();
bb=bboxes[i];
vol=(bb[1][0]-bb[0][0])*(bb[1][1]-bb[0][1])*(bb[1][2]-bb[0][2])
sp.makeCloud(bb[0],bb[1],rBall,sphRadFuzz,int(.25*vol/((4./3)*pi*rBall**3)),False)
O.bodies.append([utils.sphere(s[0],s[1],material=ballMat,color=colors[i]) for s
in sp])
O.engines = [
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),
Ig2_Facet_Sphere_ScGeom(),
Ig2_Wall_Sphere_ScGeom()],
[Ip2_FrictMat_FrictMat_MindlinPhys()],
[Law2_ScGeom_MindlinPhys_Mindlin()]
),
#GravityEngine(gravity=(0,-3e4,0)), # gravity artificially high, to make it
faster going ;-)
RotationEngine(ids=ids,angularVelocity=1000,rotateAroundZero=True,rotationAxis=(0,0,1)),
NewtonIntegrator(damping=0,gravity=[0,0,-9.81]),
PyRunner(command='exportForces()',iterPeriod=1)
]
def exportForces():
# Mechanical power
data = [ ] # List of data for each interaction
for i in O.interactions:
# Force by sphere on facet
if isinstance(O.bodies[i.id1].shape,Sphere) \
and isinstance(O.bodies[i.id2].shape,Facet):
fn = i.phys.normalForce
fs = i.phys.shearForce
cp = i.geom.contactPoint
pd = i.geom.penetrationDepth
d =
dict(normalForce=fn,shearForce=fs,contactPoint=cp,penetrationDepth=pd,
id1=i.id1,id2=i.id2)
data.append(d)
print(data)
# Save data
dataSim = "data{:06d}.json".format(O.iter)
with open(DataSim,"w") as f:
json.dump(data,f,indent=2)
O.dt = 1e-6
###
O.run(int(math.ceil(tEnd/O.dt))+1)
Thank's
Gauthier
--
You received this question notification because your team yade-users is
an answer contact for Yade.
_______________________________________________
Mailing list: https://launchpad.net/~yade-users
Post to : [email protected]
Unsubscribe : https://launchpad.net/~yade-users
More help : https://help.launchpad.net/ListHelp
