how to find the facet-sphere contact

Asked by nobody on 2020-04-17

I want to get the contact information about the shaft and spheres.

below is the code:
########################################
O.reset()
facets = O.bodies.append(geom.facetCylinder((0.2,0.2,0.3),0.2,0.6,wallMask=6,segmentsNumber=6))
print('the ids of the facets:',facets)

x0=0.2;y0=0.2;z0=0.2
cylinderIDS= O.bodies.append(geom.facetCylinder((x0,y0,z0),
 radius=0.08,height=0.3,orientation=Quaternion((1, 0, 0), 0),wallMask=7,segmentsNumber=6, angleRange=None))

print('here is the ids for the cylinder:',cylinderIDS)

facets2 = O.bodies.append(geom.facetCylinder((0.2,0.2,0.9),0.2,0.6,wallMask=4,segmentsNumber=6))

from yade import pack, plot,export
sp=pack.SpherePack()
sp.makeCloud((0.07,0.07,0.6),(0.33,0.33,1.2),rMean=.01)
sp.toSimulation()
print "num_spheres= ", len(O.bodies)

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
 InteractionLoop(
  # handle sphere+sphere and facet+sphere collisions
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
 NewtonIntegrator(gravity=(0,0,-98.1),damping=0.4),
 # call the checkUnbalanced function (defined below) every 2 seconds
 PyRunner(command='checkUnbalanced()',realPeriod=2),
 # call the addPlotData function every 200 steps
 #PyRunner(command='addPlotData()',iterPeriod=100)
]
O.dt=.5*PWaveTimeStep()

O.trackEnergy=True

def checkUnbalanced():
 if unbalancedForce()<.001:
  O.pause()
O.run(10000000,True)
######################### get the shaft_spheres contact ################
for i in O.interactions:
 if not i.isReal: continue
 if isinstance(O.bodies[i.id1].shape,Facet) and isinstance(O.bodies[i.id2].shape,Sphere):
  #for i.id1 in cylinderIDS:
  O.bodies[i.id2].shape.color = (0,0,0)

####################
I want to get the contact information for the inner cylinder (here I call it cylinderIDS).
when the particle falls down , I want to capture the contact force and the contact normal.

 #########: if isinstance(O.bodies[i.id1].shape,Facet) and isinstance(O.bodies[i.id2].shape,Sphere): ###### this line we can get the facet----sphere contact information. but how can I get the specific facet----sphere contact infor.
#### we know the facet IDS contains in "cylinderIDS"..
###### this line here::: for i.id1 in cylinderIDS:-----------I want to put some conditions here to confine this is the inner cylingder, but I don't know how to do it.

################################
besides, I have some another question about the i.id1 and i.id2.
If the contact is between the facet and sphere, I want to know which ids belong to facet group? which ids belong to sphere group? because here is a sentence methioned that "facet are only stored in id1 "

references:
[1]https://github.com/yade/trunk/blob/master/examples/jointedCohesiveFrictionalPM/identifBis.py

thanks!
Yong

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Jan Stránský
Solved:
2020-04-20
Last query:
2020-04-20
Last reply:
2020-04-17
Jérôme Duriez (jduriez) said : #1

Hi,

I probably wrote that sentence that "facets are only stored in i.id1", however I can not remember now what mechanism would insure it...

It's easy to check though if it's true (and adapt your scripts if it's not) with "isinstance(O.bodies[i.id1].shape,Facet)" used in the above link.

Best Jan Stránský (honzik) said : #2

Hello,

> but how can I get the specific facet----sphere contact infor.

it you know "i" is the contact you want, then
i.phys.normalForce
i.geom.normal
etc...

>for i.id1 in cylinderIDS:-----------I want to put some conditions here to confine this is the inner cylingder, but I don't know how to do it.

if i.id1 in cylinderIDS: # if instead of for :-)

if cylinderIDS is "many" IDS, you can make it a set to look up efficiently:
cylinderIDS = set(cylinderIDS)

> besides, I have some another question about the i.id1 and i.id2. If the contact is between the facet and sphere, I want to know which ids belong to facet group? which ids belong to sphere group? because here is a sentence methioned that "facet are only stored in id1 "

I am not sure if such order is strict / will be strict forever. In **my opinion** it is not a bad idea not to depend on specific order and test both id1->Facet id2->Sphere and id1->Sphere id2->Facet combinations:
###
for i in O.interactions:
   if not i.isReal: continue
   b1,b2= O.bodies[i.id1], O.bodies[i.id2]
   if isinstance(b1.shape,Facet) and isinstance(b2.shape,Sphere):
      facet = b1
      sphere = b2
   elif isinstance(b1.shape,Sphere) and isinstance(b2.shape,Facet):
      facet = b2
      sphere = b1
   else: # not facet-sphere contact
      continue
   if facet.id1 in cylinderIDS:
      ...
###

cheers
Jan

nobody (nobody01) said : #3

Thanks Jan Stránský, that solved my question.

nobody (nobody01) said : #4

Hi Jan,

I test the ids of the facets and spheres, the results show that i.id1 are set to facet.
because if I set: if isinstance(O.bodies[i.id1].shape,Facet) and isinstance(O.bodies[i.id2].shape,Sphere):

##### I can get the facet ids.

but if I set another command: if isinstance(O.bodies[i.id1].shape,Sphere) and isinstance(O.bodies[i.id2].shape,Facet):

#### the list from i.ids is empty.

best,
Yong