yade-users team mailing list archive

[Clinton Schramm <question294531@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Answered => Open

Clinton Schramm is still having a problem:
Thank you again Jan. I have done what you said, but still I can cannot
see the spheres changing color.

If I put the pyrunner('Changecolor') in the main engine, it does change
the spheres, how ever I need the spheres to change colour only after the
cylinder is given it's velocity.

To work around this I tried calling a pyrunner after this piece of code
was started. It is below, and sadly doesn't work. I think the procedure
of code implementation is correct though?

Another method I thought was running the code in the main engine, and
then setting the spheres to the base color just before the cylinder was
initiated but this did not seem to work either.

readParamsFromTable(rMean=.03,rRelFuzz=.003,maxLoad=1e5,minLoad=1e4)
# make rMean, rRelFuzz, maxLoad accessible directly as variables later
from yade.params.table import *
numspherespack1=1000
numspherespack2 = 173
Zextentbox = 0.600 #was .200
width = 0.154
length = 0.441
cylinderheight = 0.130
cylinderradius = 0.05
nosides = 10

#Set materials for spheres, and box
O.materials.append(FrictMat(young=78.71e6,poisson=0.27,frictionAngle=.5934,density=2600,label='spheres')) #Need to find youngs modulus
O.materials.append(FrictMat(young=30e9,poisson=0.2,frictionAngle=0,label='walls'))

# create box with free top, and ceate loose packing inside the box
from yade import pack, plot
O.bodies.append(geom.facetBox((0,0,0),(width/2,length/2,Zextentbox),wallMask=31,material='walls'))   #(centre),(x,y,z) Old Box with correct dimensions   
TBM=O.bodies.append(geom.facetCylinder(center=(-0.077,.1555,-.100),radius=cylinderradius,height=cylinderheight,orientation=Quaternion((1,0,0),3.1415926/2),segmentsNumber=nosides,wallMask=7,angleRange=None,closeGap=False,wire=False,highlight=False))

sp=pack.SpherePack()
sp.makeCloud((-0.077,0.180,-0.200),(0.077,-0.200,-0.155),rMean=rMean,rRelFuzz=rRelFuzz,num=numspherespack2) #Cloud for under tunnel
#May need to make cloud larger tos full up the model
sp.makeCloud((-0.077,0.180,-0.05),(0.077,-0.180,1),rMean=rMean,rRelFuzz=rRelFuzz,num=numspherespack1) #Move cloud in +z 400
sp.toSimulation(material='spheres',color=(0,1,1)) 

O.engines=[
   ForceResetter(),
   # sphere, facet, wall
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
   InteractionLoop(
      # the loading plate is a wall, we need to handle sphere+sphere, sphere+facet, sphere+wall
      [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
      [Ip2_FrictMat_FrictMat_FrictPhys()],
      [Law2_ScGeom_FrictPhys_CundallStrack()]
   ),	
   NewtonIntegrator(gravity=(0,0,-9.81),damping=0.5),
   SnapshotEngine(iterPeriod=300,fileBase='/home/clinton/Desktop/Video/1/Centri-',label='cent'),
   # the label creates an automatic variable referring to this engines
   PyRunner(command='checkUnbalanced1()',iterPeriod = 200,label='checker'),
   #PyRunner(command='Changecolor()',iterPeriod = 200),
]
O.dt=.5*PWaveTimeStep()

# the following checkUnbalanced, unloadPlate and stopUnloading functions are all called by the 'checker'
# (the last engine) one after another; this sequence defines progression of different stages of the
# simulation, as each of the functions, when the condition is satisfied, updates 'checker' to call
# the next function when it is run from within the simulation next time

# check whether the gravity deposition has already finished
# if so, add wall on the top of the packing and start the oedometric test
def checkUnbalanced1():	
   # 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<5000: return 
   # the rest will be run only if unbalanced is < .1 (stabilized packing)
   if 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(wall(max([b.state.pos[2]+b.shape.radius for b in O.bodies if isinstance(b.shape,Sphere)]),axis=2,sense=-1,material='walls'))
   global plate        # without this line, the plate variable would only exist inside this function
   plate=O.bodies[-1]  # the last particles is the plate
   # Wall objects are "fixed" by default, i.e. not subject to forces
   # prescribing a velocity will therefore make it move at constant velocity (downwards)
   plate.state.vel=(0,0,-.7)
   # 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='unloadPlate1()'

def unloadPlate1():
   # if the force on plate exceeds maximum load, start unloading
   if abs(O.forces.f(plate.id)[2])>maxLoad :
      plate.state.vel*=-3
      # next time, do not call this function anymore, but the next one (stopUnloading) instead
      checker.command='stopUnloading()'
      
def stopUnloading():

   if abs(O.forces.f(plate.id)[2])<minLoad:
	plate.state.vel*=0
	areacylinder = (.5*nosides*cylinderradius**2*(sin(degrees(360)/10)))/2
	Volcylinder = areacylinder*cylinderheight
       	height = max([b.state.pos[2]+b.shape.radius for b in O.bodies if isinstance(b.shape,Sphere)])+Zextentbox
	print 'Height of heighest sphere is: ',height
	v1=Volcylinder+sum([4./3.*pi*(b.shape.radius)**3. for b in O.bodies if isinstance(b.shape,Sphere)])
	v2=width*length*height
	Voidratio=(v2-v1)/v2*100.
	print'The Void Ratio:',Voidratio,'%','This is the second one'
	
    	checker.command='checkUnbalancedtunnelmovement()'
	
def checkUnbalancedtunnelmovement():
	O.engines=O.engines+[TranslationEngine(translationAxis=[0,1,0],velocity=0.3,ids=TBM)]
	# start plotting the data now, it was not interesting before
   	#O.engines=O.engines+[PyRunner(command='addPlotDataForceTunnel()',iterPeriod=200)]
	# we use it below to changefor s in O.bodies if isinstance(s.shape,Sphere):
	global CylinderforceX,  CylinderforceY, CylinderforceZ
	O.engines=O.engines+[PyRunner(command='Changecolor()',iterPeriod=200)]
	
   	CylinderforceX = utils.sumForces([10,14,16,22,26,30,34,36,42,46],(1,0,0))
	CylinderforceY = utils.sumForces([10,14,16,22,26,30,34,36,42,46],(0,1,0))
	CylinderforceZ = utils.sumForces([10,14,16,22,26,30,34,36,42,46],(0,0,1))
 	
	Cylindernormalforce = sumFacetNormalForces([10,14,16,22,26,30,34,36,42,46],1)
	#print 'The sum of normal forces are: ',Cylindernormalforce
   
	plot.addData(CylinderforceX=CylinderforceX,tx=(O.bodies[20].state.pos[1]-O.bodies[20].state.refPos[1]),CylinderforceY=CylinderforceY,ty=(O.bodies[20].state.pos[1]-O.bodies[20].state.refPos[1]),CylinderforceZ=CylinderforceZ,tz=(O.bodies[20].state.pos[1]-O.bodies[20].state.refPos[1]),Cylindernormalforce=Cylindernormalforce,t=(O.bodies[20].state.pos[1]-O.bodies[20].state.refPos[1]),i=O.iter)
   	# start plotting the data now, it was not interesting before
   	#O.engines=O.engines+[PyRunner(command='addPlotDataForceTunnel()',iterPeriod=200)]
   	# next time, do not call this function anymore, but the next one (unloadPlate) instead
   	#O.pause()

#plot.saveDataTxt(O.tags['d.id']+'.txt')

def Changecolor():   
   for s in O.bodies:
	s.shape.color=scalarOnColorScale(s.state.displ().norm(),0,2)
	#print 0.1*s.state.displ()
		# for changing color based on the movement of sheres. its attributes from the functions called
 	#print 'Change'  


def addPlotData():

   if not isinstance(O.bodies[-1].shape,Wall):
      plot.addData(); return
   Fz=O.forces.f(plate.id)[2]
   plot.addData(Fz=Fz,w=plate.state.pos[2]-plate.state.refPos[2],unbalanced=unbalancedForce(),i=O.iter)
   
#Plotting the force on the tunnel face
def addPlotDataForceTunnel():
  print 'test => ',CylinderforceX
  plot.addData(CylinderforceX,tx=(O.bodies[20].state.pos[1]-O.bodies[20].state.refPos[1]),i=O.iter)


# besides unbalanced force evolution, also plot the displacement-force diagram
plot.plots={'tx':('CylinderforceX',),'ty':('CylinderforceY',),'tz':('CylinderforceZ',),'t':('Cylindernormalforce',)}
#plot.plots={'i':('unbalanced',),'w':('Fz',),'tx':('Facex',)}

plot.plot()


#O.saveTmp()
from yade import qt
v=qt.View()
v.eyePosition=(-1,0.3,0.25);v.upVector=(0,0,1); v.lookAt=(0,0,0); v.axes=False; v.sceneRadius=1.9

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.