New question #429604 on Yade:
https://answers.launchpad.net/yade/+question/429604
Hi all,
I am having a problem when I use a cylinder with periodic boundary condition.
I made a simple script to reproduce the problem. The script is a column of
grain falling inside a box under gravity, with a cylinder fixed on the lower
(infinite) plane of the box. When I do that, everything works fine (option
periodic = 0 and cylinderBottom = 1 in the script below).
If now, I remove the sides of the box and replace it with periodic boundary
condition (periodic = 1 and cylinderBottom = 1 in script below), some particles
literally goes inside the cylinder. And this whatever the stiffness of the
interaction, the time step used, the contact law used
(Law2_ScGeom_ViscElPhys_Basic or Law2_ScGeom_FrictPhys_CundallStrack) and the
length of the cylinder (which I made much larger than the periodic cell).
This does not happen if I consider instead of the cylinder a row of fixed
spheres (periodic = 1 and cylinderBottom = 0 in the script below).
There is also a strange behavior when defining the cylinder position and
extents, which might be linked to the present problem:
- It is not possible to define a cylinder which extents corresponds to an
integer of the periodic cell extent, e.g. cylinder(begin =
(L/2.,-100*L,ground+dp/2.),end =(L/2.,100*L,ground+dp/2.),.... It considers
that the length is zero at this moment, probably applying
However, if you interchange begin and end, it works, i.e. cylinder(begin =
(L/2.,100*L,ground+dp/2.),end =(L/2.,-100*L,ground+dp/2.),....
- The 3D view shows very strange position of the cylinder, which do not
correspond to its actual position, when you consider the position of the
cylinder nodes.
Do you have any idea what the problem could be linked to ?
Raphael
Script:
from yade import pack
from yade.gridpfacet import *
dp= 6e-3
periodic = 1
cylinderBottom = 1
ground = 0.2*dp
H = 300*dp
L = 20*dp
O.engines = [
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb(),Bo1_GridConnection_Aabb()],allowBiggerThanPeriod
= True),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom(),Ig2_Sphere_GridConnection_ScGridCoGeom()],
[Ip2_ViscElMat_ViscElMat_ViscElPhys()],
[Law2_ScGeom_ViscElPhys_Basic()]
# [Ip2_FrictMat_FrictMat_FrictPhys()],
# [Law2_ScGeom_FrictPhys_CundallStrack()]
),
PyRunner(command = 'check()', virtPeriod = 0.1),
NewtonIntegrator(damping=0., gravity = (0,0,-9.81))
]
#Material creation
O.materials.append(ViscElMat(en=0.5, et=1., young=5e6, poisson=0.5,
density=2500., frictionAngle=0.4, label='Mat'))
#O.materials.append(FrictMat(young = 5e6, poisson = 0.5,
density=2500,frictionAngle=0.4, label='Mat'))
#Periodic Cell or containing box
if periodic:
O.periodic = True
O.cell.setBox(L,L,H)
O.bodies.append(box(center=
(L/2.,L/2.,ground),extents=(100,100,0),fixed=True,color = (0.,1.,0.),wire =
True,material = 'Mat'))
else:
O.bodies.append(box(center= (L/2.,L/2.,ground),extents=(100,100,0),wire
= True,fixed=True,color = (0.,1.,0.),material = 'Mat'))#Made invisible to see
below
O.bodies.append(box(center= (0,L/2.,H/2.),extents=(0,L/2.,H/2.),wire =
True,fixed=True,color = (0.,1.,0.),material = 'Mat'))#Made invisible to see
O.bodies.append(box(center=
(L,L/2.,H/2.),extents=(0,L/2.,H/2.),fixed=True,color = (0.,1.,0.),material =
'Mat'))
O.bodies.append(box(center= (L/2.,0,H/2.),extents=(L/2.,0,H/2.),wire =
True,fixed=True,color = (0.,1.,0.),material = 'Mat'))#Made invisible to see
O.bodies.append(box(center=
(L/2.,L,H/2.),extents=(L/2.,0,H/2.),fixed=True,color = (0.,1.,0.),material =
'Mat'))
#Bottom fixed cylinder or fixed particles
if cylinderBottom:
n = len(O.bodies)
cylinder(begin = (L/2.,100*L,ground+dp/2.),end
=(L/2.,-100*L,ground+dp/2.),radius = dp/2.,fixed = True,color =
(0,0,1),intMaterial = 'Mat',extMaterial = 'Mat')#Made invisible to see inside
else:
for n in range(0,int(L/dp)+1):
O.bodies.append(sphere(center = (L/2.,n*dp,ground+dp/2.),radius
= dp/2.,fixed = True,color = (0,0,1),wire = True,material = 'Mat'))#Made
invisible to see inside
#Particle cloud for gravity deposition
partCloud = pack.SpherePack()
partCloud.makeCloud(minCorner=(4*L/10.,4*L/10.,ground+dp),maxCorner=(6*L/10.,6*L/10.,H),rRelFuzz=0.,
rMean=dp/2., num = 2000)
partCloud.toSimulation(material='Mat')
O.saveTmp()
O.dt = 1e-5
#Function to check if the center of a particle is contained inside the cylinder
or pseudo cylinder
def check():
for b in O.bodies:
if b.dynamic and b.state.pos[2]<ground+dp:
posCylXZ = Vector3(L/2.,0.,ground+dp/2.)
pRelX = b.state.pos[0]%L - posCylXZ[0]
pRelZ = b.state.pos[2] - posCylXZ[2]
if sqrt(pow(pRelX,2) + pow(pRelZ,2))<dp/2.:#If the
particle center is closer than the radius (i.e. if it is completely inside!)
delta = dp/2.-sqrt(pow(pRelX,2) + pow(pRelZ,2))
print('particle {0} is inside the cylinder from
an amount of {1} diameter'.format(b.id,delta/dp))
O.pause()
--
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
