yade-dev team mailing list archive

[luc scholtes <lscholtes63@xxxxxxxxx>]

Hi all,

Maybe I should have report a bug for that, but I do not really get where is
the issue...

I am trying to identify spheres which could belong to a joint in a given
(bonded) granular mass using python. The idea is to work on interactions in
order to have every sphere couples which interaction is "cut" by a joint. In
addition, joints can be persistent (infinite extents) or not, and can have
as much as orientation as possible.

Until now, I have chosen to work with facets in order to define the joint
plane. This could allow a real flexibility with, for example, an import of
the joint geometry from any kind of mesh. So, given an initial spheres
packing and a joint defined by a facet, I loop on interactions (after 1
iteration) to get the spheres that interact with the facet (see the python
script).  Nothing exciting until now, but already an issue... Indeed, it
seems that some spheres-facet interactions are not listed in the interaction
container (see the snapshot I attached where the "SphereOnPlane" should be
red). I tried with several orientation of the facet, and I always get an
"error" in the estimation except for the case where the facet is parallel to
one of the axis. Moreover, depending on the iteration, spheres are not the
same...

If you have any advices...

  Luc

Attachment: snapshot_SphereOnJoint.png
Description: PNG image

# -*- coding: utf-8 -*-
# encoding: utf-8

O=Omega()

from yade import pack,utils

############################ Definition of the ground (with facetBox)
facetMat=O.materials.append(GranularMat(young=2e7,frictionAngle=radians(30),poisson=0.2,density=3000))
Ground=O.bodies.append(utils.facetBox((0,0,-0.5),(10,10,-0.5),wallMask=16,color=(0,0,1),wire=False,material=facetMat)) 

############################ definition of a joint via a facet
Joint=O.bodies.append(utils.facet([[-5,-5,0],[0,6,5],[6,0,5]],dynamic=False,noBound=False,wire=False,color=(1,1,1)))
## How to get the id of the Joint? 
#O.bodies['Joint'] does not work
## Is the normal to the facet can be define?
## Rk: the facet will not interact if noBound=True -> OK but how can we identify therefore the jointed interaction?

############################ Definition of the sphere assembly
sphereMat=O.materials.append(GranularMat(young=2e7,frictionAngle=radians(30),poisson=0.2,density=3000))
pred=pack.inAlignedBox((-1,-1,0.1),(1,1,6.1))
#O.bodies.append(pack.randomDensePack(pred,radius=0.3,rRelFuzz=.2,memoizeDb='/tmp/triaxPackCache.sqlite',material=sphereMat))
O.bodies.append(pack.regularOrtho(pred,radius=0.3,gap=0,material=sphereMat))
#O.bodies.append(pack.regularHexa(pred,radius=0.3,gap=-0.05,material=sphereMat))

############################# no color argument in randomDensepack -> to assign a color (blue) to the assembly
for i in O.bodies:
    i.shape['diffuseColor']=(0,0,1)

print len(O.bodies)
###########################

O.initializers=[BoundDispatcher([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()])]

O.engines=[

	BexResetter(),
	BoundDispatcher([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
	InsertionSortCollider(),

	InteractionGeometryDispatcher([Ig2_Sphere_Sphere_Dem3DofGeom(),Ig2_Facet_Sphere_Dem3DofGeom()]),
	InteractionPhysicsDispatcher([Ip2_BMP_BMP_CSPhys()]),
	LawDispatcher([Law2_Dem3Dof_CSPhys_CundallStrack()]),

	GravityEngine(gravity=[0,0,-9.81]),
	NewtonIntegrator(damping=0.10)

]

############################  

O.dt=0.3*utils.PWaveTimeStep()
print 'Timestep',O.dt

from yade import qt
v=qt.Controller()
v=qt.View()

O.run(1); O.wait()

############## Determining spheres on joint through the interactions and giving them another color (red) -> does not work!!!

for i in O.interactions:
    if not i.isReal : continue

    if O.bodies[i.id1].shape.name=='Facet' or O.bodies[i.id2].shape.name=='Facet': 

	O.bodies[i.id1].shape['diffuseColor']=(1,0,0)
	O.bodies[i.id2].shape['diffuseColor']=(1,0,0)

## How to access from one grain to its neighbours ??????????? -> O.interactions.withBody(id)

O.save('/media/DATAPART1/Luc/YADE/WorkingArea/bzr/bin/spheresOnJoint.xml');