← Back to team overview

yade-users team mailing list archive

  1. yade-users team
  2. Mailing list archive
  3. Message #23092

Re: [Question #690777]: Polyhedrons escaping the shear box walls during deposition

  Thread Previous • Date Previous • • Thread Next 
Question #690777 on Yade changed:
https://answers.launchpad.net/yade/+question/690777

    Status: Needs information => Open

Ali HAIDAR gave more information on the question:
Hello,
Sorry about the link.
The script is the following:
But it needs the polyhedrons files I am using to run.
Thank you

from yade import plot,polyhedra_utils
import pickle

with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion1", "rb") as fp:
  v1 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion2", "rb") as fp:
  v2 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion3", "rb") as fp:
  v3 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion4", "rb") as fp:
  v4 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion5", "rb") as fp:
  v5 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion6", "rb") as fp:
  v6 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion7", "rb") as fp:
  v7 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion8", "rb") as fp:
  v8 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion9", "rb") as fp:
  v9 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion10", "rb") as fp:
  v10 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion11", "rb") as fp:
  v11 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion12", "rb") as fp:
  v12 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion13", "rb") as fp:
  v13 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion14", "rb") as fp:
  v14 = pickle.load(fp)
with open("/home/alihaidar/Desktop/data/Fusion/pickle_reducedvertices/tolerance=20mm/20-preducedfusion15", "rb") as fp:
  v15 = pickle.load(fp)

m = PolyhedraMat()
m.density = 2607 #kg/m^3 
m.frictionAngle = 0.3 #rad
#m.young = 1E10 #Pa
#m.poisson = 20000/1E6
m2= PolyhedraMat()
m2.frictionAngle = 0.3 #rad
m2.young = 1E100
x=.3
y=.3
z=.2
a=b=0.04
i=0
while a<=(x-0.05):
  while b<=(y-0.05):
    O.bodies.append([
                     polyhedra_utils.polyhedra(material=m,v=v1),
                     polyhedra_utils.polyhedra(material=m,v=v2),
                     polyhedra_utils.polyhedra(material=m,v=v3),
                     polyhedra_utils.polyhedra(material=m,v=v4),
                     polyhedra_utils.polyhedra(material=m,v=v5),
                     polyhedra_utils.polyhedra(material=m,v=v6),
                     polyhedra_utils.polyhedra(material=m,v=v7),
                     polyhedra_utils.polyhedra(material=m,v=v8),
                     polyhedra_utils.polyhedra(material=m,v=v9),
                     polyhedra_utils.polyhedra(material=m,v=v10),
                     polyhedra_utils.polyhedra(material=m,v=v11),
                     polyhedra_utils.polyhedra(material=m,v=v12),
                     polyhedra_utils.polyhedra(material=m,v=v13),
                     polyhedra_utils.polyhedra(material=m,v=v14),
                     polyhedra_utils.polyhedra(material=m,v=v15),
		    ]
		   )
    if i==0 or i==5 or i==13 or i==10:
      O.bodies[-1].state.pos=(a,b,.04-z)
      O.bodies[-2].state.pos=(a,b,.12-z)
      O.bodies[-3].state.pos=(a,b,.2-z)
      O.bodies[-4].state.pos=(a,b,.28-z)
      O.bodies[-5].state.pos=(a,b,.36-z)
      O.bodies[-6].state.pos=(a,b,.44-z)
      O.bodies[-7].state.pos=(a,b,.52-z)
      O.bodies[-8].state.pos=(a,b,.6-z)
      O.bodies[-9].state.pos=(a,b,.68-z)
      O.bodies[-10].state.pos=(a,b,.76-z)
      O.bodies[-11].state.pos=(a,b,.84-z)
      O.bodies[-12].state.pos=(a,b,.92-z)
      O.bodies[-13].state.pos=(a,b,1-z)
      O.bodies[-14].state.pos=(a,b,1.08-z)
      O.bodies[-15].state.pos=(a,b,1.16-z)
    if i==1 or i==6 or i==12 or i==8:
      O.bodies[-15].state.pos=(a,b,.04-z)
      O.bodies[-14].state.pos=(a,b,.12-z)
      O.bodies[-13].state.pos=(a,b,.2-z)
      O.bodies[-12].state.pos=(a,b,.28-z)
      O.bodies[-11].state.pos=(a,b,.36-z)
      O.bodies[-10].state.pos=(a,b,.44-z)
      O.bodies[-9].state.pos=(a,b,.52-z)
      O.bodies[-8].state.pos=(a,b,.6-z)
      O.bodies[-7].state.pos=(a,b,.68-z)
      O.bodies[-6].state.pos=(a,b,.76-z)
      O.bodies[-5].state.pos=(a,b,.84-z)
      O.bodies[-4].state.pos=(a,b,.92-z)
      O.bodies[-3].state.pos=(a,b,1-z)
      O.bodies[-2].state.pos=(a,b,1.08-z)
      O.bodies[-1].state.pos=(a,b,1.16-z)
    if i==4 or i==15 or i==11 or i==9:
      O.bodies[-12].state.pos=(a,b,.04-z)
      O.bodies[-3].state.pos=(a,b,.12-z)
      O.bodies[-9].state.pos=(a,b,.2-z)
      O.bodies[-5].state.pos=(a,b,.28-z)
      O.bodies[-10].state.pos=(a,b,.36-z)
      O.bodies[-15].state.pos=(a,b,.44-z)
      O.bodies[-4].state.pos=(a,b,.52-z)
      O.bodies[-1].state.pos=(a,b,.6-z)
      O.bodies[-6].state.pos=(a,b,.68-z)
      O.bodies[-8].state.pos=(a,b,.76-z)
      O.bodies[-2].state.pos=(a,b,.84-z)
      O.bodies[-7].state.pos=(a,b,.92-z)
      O.bodies[-11].state.pos=(a,b,1-z)
      O.bodies[-14].state.pos=(a,b,1.08-z)
      O.bodies[-13].state.pos=(a,b,1.16-z)
    if i==3 or i==7 or i==2 or i==14:
      O.bodies[-13].state.pos=(a,b,.04-z)
      O.bodies[-14].state.pos=(a,b,.12-z)
      O.bodies[-11].state.pos=(a,b,.2-z)
      O.bodies[-7].state.pos=(a,b,.28-z)
      O.bodies[-2].state.pos=(a,b,.36-z)
      O.bodies[-8].state.pos=(a,b,.44-z)
      O.bodies[-6].state.pos=(a,b,.52-z)
      O.bodies[-1].state.pos=(a,b,.6-z)
      O.bodies[-4].state.pos=(a,b,.68-z)
      O.bodies[-15].state.pos=(a,b,.76-z)
      O.bodies[-10].state.pos=(a,b,.84-z)
      O.bodies[-5].state.pos=(a,b,.92-z)
      O.bodies[-9].state.pos=(a,b,1-z)
      O.bodies[-3].state.pos=(a,b,1.08-z)
      O.bodies[-12].state.pos=(a,b,1.16-z)
    b=b+0.07
    i=i+1
  a=a+0.07
  b=0.04


f1=utils.facet([(0,0,0),(0,0,z),(x,0,0)],material=m2)
f2=utils.facet([(x,0,z),(0,0,z),(x,0,0)],material=m2)
f3=utils.facet([(x,0,z),(x,y,z),(x,0,0)],material=m2)
f4=utils.facet([(x,y,0),(x,y,z),(x,0,0)],material=m2)
f5=utils.facet([(x,y,0),(x,y,z),(0,y,z)],material=m2)
f6=utils.facet([(x,y,0),(0,y,0),(0,y,z)],material=m2)
f7=utils.facet([(0,0,0),(0,y,0),(0,y,z)],material=m2)
f8=utils.facet([(0,0,0),(0,0,z),(0,y,z)],material=m2)
f9=utils.facet([(0,0,0),(0,-.5*y,0),(x,-.5*y,0)],material=m2)
f10=utils.facet([(0,0,0),(x,-.5*y,0),(x,0,0)],material=m2)
f11=utils.facet([(0,0,z),(0,0,5*z),(x,0,z)],material=m2)
f12=utils.facet([(x,0,5*z),(0,0,5*z),(x,0,z)],material=m2)
f13=utils.facet([(x,0,5*z),(x,y,5*z),(x,0,z)],material=m2)
f14=utils.facet([(x,y,z),(x,y,5*z),(x,0,z)],material=m2)
f15=utils.facet([(x,y,z),(x,y,5*z),(0,y,5*z)],material=m2)
f16=utils.facet([(x,y,z),(0,y,z),(0,y,5*z)],material=m2)
f17=utils.facet([(0,0,z),(0,y,z),(0,y,5*z)],material=m2)
f18=utils.facet([(0,0,z),(0,0,5*z),(0,y,5*z)],material=m2)
O.bodies.append([f1,f2,f3,f4,f5,f6,f7,f8,f9,f10,f11,f12,f13,f14,f15,f16,f17,f18])
f19=utils.facet([(0,y,0),(0,1.5*y,0),(x,y,0)],fixed=True,material=m2)
f20=utils.facet([(x,y,0),(x,1.5*y,0),(0,1.5*y,0)],fixed=True,material=m2)
f=0
t=0
n=0
A=.3*.3
V=0.1
F=-13500
O.bodies.append([
		utils.facet([(0,0,0),(0,0,-z),(x,0,0)],material=m2,fixed=True),
		utils.facet([(x,0,-z),(0,0,-z),(x,0,0)],material=m2,fixed=True),
		utils.facet([(x,0,-z),(x,y,-z),(x,0,0)],material=m2,fixed=True),
		utils.facet([(x,y,0),(x,y,-z),(x,0,0)],material=m2,fixed=True),
		utils.facet([(x,y,0),(x,y,-z),(0,y,-z)],material=m2,fixed=True),
		utils.facet([(x,y,0),(0,y,0),(0,y,-z)],material=m2,fixed=True),
		utils.facet([(0,0,0),(0,y,0),(0,y,-z)],material=m2,fixed=True),
		utils.facet([(0,0,0),(0,0,-z),(0,y,-z)],material=m2,fixed=True),
		utils.facet([(x,0,-z),(0,0,-z),(0,y,-z)],material=m2,fixed=True),
		utils.facet([(x,0,-z),(x,y,-z),(0,y,-z)],material=m2,fixed=True),f19,f20
	       ])

O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Polyhedra_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()],verletDist=0.1),
	InteractionLoop(
			[Ig2_Polyhedra_Polyhedra_PolyhedraGeom(),Ig2_Facet_Polyhedra_PolyhedraGeom(),Ig2_Wall_Polyhedra_PolyhedraGeom()],
			[Ip2_PolyhedraMat_PolyhedraMat_PolyhedraPhys()],
			[Law2_PolyhedraGeom_PolyhedraPhys_Volumetric()]
		       ),
	NewtonIntegrator(gravity=(0,0,-9.81),damping=0.5),
	PyRunner(command='checkUnbalanced()',realPeriod=2,label='checker'),	
	  ]	

O.dt=.05*polyhedra_utils.PWaveTimeStep()


def checkUnbalanced():

		if O.iter<100000: return 	
		if unbalancedForce()>.2: return 
		m.frictionAngle=.5 
		m2.frictionAngle=.5
		for i in O.interactions: i.phys.tangensOfFrictionAngle=tan(.5)
		O.bodies.append(wall(max([b.state.pos[2]+0.05 for b in O.bodies if isinstance(b.shape,Polyhedra)]),axis=2,sense=-1,material=m2))
		global plate        
		plate=O.bodies[-1] 
		plate.state.blockedDOFs='xyXYZ'
		plate.state.mass=.1
		O.forces.setPermF(plate.id,(0,0,F))
		checker.command='pos()'
def pos():
		global l
		l=plate.state.pos[2]
		checker.command='shear()'

def shear():
		f1.state.vel=(0,V,0)
		f2.state.vel=(0,V,0)
		f3.state.vel=(0,V,0)
		f4.state.vel=(0,V,0)
		f5.state.vel=(0,V,0)
		f6.state.vel=(0,V,0)
		f7.state.vel=(0,V,0)
		f8.state.vel=(0,V,0)
		f9.state.vel=(0,V,0)
		f10.state.vel=(0,V,0)
		f11.state.vel=(0,V,0)
		f12.state.vel=(0,V,0)
		f13.state.vel=(0,V,0)
		f14.state.vel=(0,V,0)
		f15.state.vel=(0,V,0)
		f16.state.vel=(0,V,0)
		f17.state.vel=(0,V,0)
		f18.state.vel=(0,V,0)	
		O.engines=O.engines+[PyRunner(command='addData()',iterPeriod=100)]
def addData():
		
		f=-sum(O.forces.f(b.id)[1] for b in [f1,f2,f3,f4,f5,f6,f7,f8,plate,f19,f20])
		A=(.3*(.3-V*O.dt*O.iter))
		t=f/A/1000
		print('Shear force',f,'Contact Area',A,'shear stress',t,'KPa')
		n=(O.forces.f(plate.id)[2]+F)/(.3*.3)/1000
		print('normal stress',n,'KPa')
		d=f1.state.pos[1]
		d2=(plate.state.pos[2]-l)
		plot.addData(Shear_displacement=d,Normal_disp=d2,Shear_Stress=t,Normal_Stress=n,i=O.iter)

plot.plots={'Shear_displacement':('Shear_Stress'),'Shear_displacement
':('Normal_Stress'),'Shear_displacement  ':'Normal_disp'}

plot.plot()
yade.qt.View()	
O.run()
O.saveTmp()

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.
  Thread Previous • Date Previous • • Thread Next