← Back to team overview

yade-users team mailing list archive

Re: [Question #215627]: Translation and rotation engines

 

Question #215627 on Yade changed:
https://answers.launchpad.net/yade/+question/215627

m.asd posted a new comment:
dear chareyre,
I had changed the script and nothing changed, please let me know if there is any mistake in the script, as I cannot find and fix it.


overburden=8
from yade import pack,plot,qt
idmat=O.materials.append(FrictMat(density=2700,young=1e9,poisson=0.5,frictionAngle=0))
O.bodies.append(utils.wall(position=0,axis=2,sense=1))
O.bodies.append(utils.wall(position=0,axis=1,sense=1))
O.bodies.append(utils.wall(position=0.01,axis=1,sense=-1))
O.bodies.append(utils.wall(position=0,axis=0,sense=1))
O.bodies.append(utils.wall(position=0.01,axis=0,sense=-1))
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(0.01,0.01,0.02),rMean=.0008,rRelFuzz=0)
sp.toSimulation(color=(0,0,1)) # blue

O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Wall_Aabb()],verletDist=0),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
		[Ip2_FrictMat_FrictMat_FrictPhys()],
		[Law2_ScGeom_FrictPhys_CundallStrack()]
	),
	GravityEngine(gravity=(0,0,-9.81)),
	NewtonIntegrator(damping=0.6),
	qt.SnapshotEngine(fileBase='3d-',iterPeriod=233000,label='snapshot'),
	PyRunner(command='sphereData()',iterPeriod=2000000),
	PyRunner(command='interactionData()',iterPeriod=2000000),
	PyRunner(command='unbalancedData()',iterPeriod=2000000),
	PyRunner(command='plateData()',iterPeriod=2000000),
	PyRunner(command='side1Data()',iterPeriod=2000000),
	PyRunner(command='side2Data()',iterPeriod=2000000),
	PyRunner(command='checkUnbalanced()',realPeriod=1,label='checker'),
]

O.dt=1e-7

def checkUnbalanced():
	if O.iter<5000000: return
	if utils.unbalancedForce()>.05: return
	posi=O.bodies[5].state.pos[2]+O.bodies[5].shape.radius
	m=5
	while m<len(O.bodies):
		b=O.bodies[m]
		try:
			if isinstance(b.shape,Sphere):
				posi2=b.state.pos[2]+b.shape.radius
				if posi2>=posi: posi=posi2
				m=m+1
		except: m=m+1
	O.bodies.append(utils.wall(position=posi,axis=2,sense=-1)) 
	global plate
	plate=O.bodies[-1]
	plate.state.vel=(0,0,-.002)
	checker.command='shearing()' 
	O.materials[0].frictionAngle=.2 
	for i in O.interactions: i.phys.tangensOfFrictionAngle=tan(.2)

def shearing():
	if abs(O.forces.f(plate.id)[2])<overburden: return
	plate.state.vel*=0
	O.bodies[3].state.angVel=(0,0.1,0)
	O.bodies[4].state.angVel=(0,0.1,0)
	checker.command='stopShearing()'

def stopShearing():
	if O.bodies[3].state.rot<(0,0.001,0): return
	O.pause(); O.wait()
	plot.saveGnuplot('shearplot')
	plot.saveDataTxt('shearunbal.txt.bz2',vars=('itperi','unbalanced'))
	plot.saveDataTxt('shearsphere.txt.bz2',vars=('itperi','bodyid','radi','posiX','posiY','posiZ','oriX','oriY','oriZ','oripi'))
	plot.saveDataTxt('shearinteraction.txt.bz2',vars=('itperi','bodynum1','bodynum2','normstiff','shearstiff','normforceX','normforceY','normforceZ','shearforceX','shearforceY','shearforceZ'))
	plot.saveDataTxt('shearwall1.txt.bz2',vars=('itperi','rotor1X','rotor1Y','rotor1Z','rotor1pi','F1x','F1z'))
	plot.saveDataTxt('shearwall2.txt.bz2',vars=('itperi','rotor2X','rotor2Y','rotor2Z','rotor2pi','F2x','F2z'))
	plot.saveDataTxt('shearplate.txt.bz2',vars=('itperi','Fz','w'))
	print 'Finished'

def unbalancedData():
	plot.addData(itperi=O.iter,unbalanced=utils.unbalancedForce())
	
def sphereData():
	i=5
	while i<len(O.bodies):
		if isinstance(O.bodies[-1].shape,Wall):
			plot.addData(); return
		b=O.bodies[i]
		plot.addData(itperi=O.iter,bodyid=b.id,posiX=b.state.pos[0],posiY=b.state.pos[1],posiZ=b.state.pos[2],radi=b.shape.radius,oriX=b.state.ori[0],oriY=b.state.ori[1],oriZ=b.state.ori[2],oripi=b.state.ori[3])
		i=i+1

def interactionData():
	i=5  #id num of bodies for interaction detection starting from spheres
	while i<len(O.bodies):
		id1=i
		j=0
		while j<len(O.bodies):
			id2=j
			try:
				interact=O.interactions[id1,id2]
				plot.addData(itperi=O.iter,bodynum1=id1,bodynum2=id2,normstiff=interact.phys.kn,shearstiff=interact.phys.ks,normforceX=interact.phys.normalForce[0],normforceY=interact.phys.normalForce[1],normforceZ=interact.phys.normalForce[2],shearforceX=interact.phys.shearForce[0],shearforceY=interact.phys.shearForce[1],shearforceZ=interact.phys.shearForce[2])
				j=j+1
			except:
				j=j+1
		i=i+1

def plateData():
	if not isinstance(O.bodies[-1].shape,Wall):
		plot.addData(); return
	Fz=O.forces.f(plate.id)[2]
	w=plate.state.refPos[2]-plate.state.pos[2]
	plot.addData(Fz=Fz,w=w)

def side1Data():
	if not isinstance(O.bodies[3].shape,Wall):
		plot.addData(); return
	rotor1X=O.bodies[3].state.ori[0]
	rotor1Y=O.bodies[3].state.ori[1]
	rotor1Z=O.bodies[3].state.ori[2]
	rotor1pi=O.bodies[3].state.ori[3]
	F1x=O.forces.f(O.bodies[3].id)[0]
	F1z=O.forces.f(O.bodies[3].id)[2]
	plot.addData(itperi=O.iter,F1x=F1x,F1z=F1z,rotor1X=rotor1X,rotor1Y=rotor1Y,rotor1Z=rotor1Z,rotor1pi=rotor1pi)

def side2Data():
	if not isinstance(O.bodies[4].shape,Wall):
		plot.addData(); return
	rotor2X=O.bodies[4].state.ori[0]
	rotor2Y=O.bodies[4].state.ori[1]
	rotor2Z=O.bodies[4].state.ori[2]
	rotor2pi=O.bodies[4].state.ori[3]
	F2x=O.forces.f(O.bodies[4].id)[0]
	F2z=O.forces.f(O.bodies[4].id)[2]
	plot.addData(itperi=O.iter,F2x=F2x,F2z=F2z,rotor2X=rotor2X,rotor2Y=rotor2Y,rotor2Z=rotor2Z,rotor2pi=rotor2pi)

plot.plots={'itperi':('unbalanced',)}
plot.plot() 

Gl1_Sphere.stripes=True
qt.View()

rr=yade.qt.Renderer()
rr.shape=True

-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.