yade-users team mailing list archive

[procrastinator <question684615@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Answered => Open

procrastinator is still having a problem:
Hi Jan,

I think I figured out this problem, this code doesn't report any error.
But I found another problem.

as you can see from my code, that I want to divide the particles into
three layers and get the per particle stress for each layer separately.
when I print the number of ids in the ball_list3 (which means the
particle numbers in the top layer).  the print number is not the same as
the particle numbers show in the show3D window view.

here is the MWE.

from yade import pack, plot, export
O.bodies.append(geom.facetBox((.5,.5,.5),(.5,.5,.5),wallMask=31))
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(1,1,1),rMean=.08,rRelFuzz=.5)
sp.toSimulation()
O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
		[Ip2_FrictMat_FrictMat_FrictPhys()],
		[Law2_ScGeom_FrictPhys_CundallStrack()]
	),
	PyRunner(command='checkUnbalanced()',realPeriod=2),
	PyRunner(command='addPlotData()',iterPeriod=100),
	PyRunner(command='subbox()',iterPeriod=100),
	PyRunner(command='stress_export()',iterPeriod=100),
	NewtonIntegrator(gravity=(0,0,-9.81),damping=0.4),
]
O.dt=.5*PWaveTimeStep()
print(len(O.bodies))
def checkUnbalanced():
	if unbalancedForce()<.05:
		O.pause()
		#plot.saveDataTxt('bbb.txt')
def addPlotData():
	plot.addData(i=O.iter,unbalanced=unbalancedForce())
#plot.plots={'i':('unbalanced')}
#plot.plot()
######################## sub-region and coloring#############################
ball_list1 =[]
ball_list2 =[]
ball_list3 =[]
def subbox():
	for b in O.bodies:
		if 0< b.state.pos[2] <= 0.15:
			if isinstance(b.shape,Sphere):
				b.shape.color =Vector3(255,255,255) # white
				m = b.id
				ball_list1.append(m)
		elif 0.15 < b.state.pos[2] <= 0.3:
			if isinstance(b.shape,Sphere):
				b.shape.color =Vector3(0,100,000) # green
				n = b.id
				ball_list2.append(n)			
		elif    b.state.pos[2] > 0.3:
			if isinstance(b.shape,Sphere):
				b.shape.color =Vector3(00,00,255) # blue
				q = b.id
				ball_list3.append(q)
######################## sub-region and coloring#############################

######################################## get the per-perticle stess#####################
TW=TesselationWrapper()
TW.setState()
TW.computeVolumes()
vtk = export.VTKExporter("stress_field")
def stress_export():
	s=bodyStressTensors()
	for i in ball_list3:
  		b = O.bodies[i]
		if isinstance(b.shape,Sphere):
			b.mystress = s[b.id]*4.*pi/3.*b.shape.radius**3/TW.volume(b.id)	
	vtk.exportSpheres(ids=ball_list3,what=[('pos','b.state.pos'),('stress','b.mystress')],useRef=False)
######################################## get the per-perticle stess#####################


print("this is the number in region 3:",len(ball_list3))  #### this number is not the same as in the 3D window view, at the same time this is not the same as in the Paraview plot.


I think maybe there is something wrong with the    ball_list. append, but right now I haven't figured it out,  can you help me to run the code and give me some suggestion.

best,
yong

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