Geometrical intersection between 2 objects

Asked by Kneib François

Hi all,

I have a cylinder and a sphere, and I have to determine if these two objects are penetrated. Geometricaly, if there is an non-void intersection between the two bodies.
I can't do that using O.interactions because simulation is paused, spheres appears and then I have to delete the ones who are penetrated with cylinders ...
I tried using the predicates but I don't think there is a function who determines if the two predicates are penetrated or not.
Is there something in yade who can do this, or must I calculate manually with mathematical stuffs ?

Thanks

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Kneib François
Solved:
Last query:
Last reply:
Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#1

There is no function that you can use straight away, but the code in
Cylinder.cpp lines 56 to 81 is doing exactly that.

I tried to see how it could be isolated in a single function, it looks
like this (below). The only problem is this function can't be wrapped
with python yet (I think) because the arguments are shared pointers.
You'd need to replace the arguments with id1/id2, then define the
cm1/cm2/state1/state2 in the function's body. Finally, you declare the
function in the hpp with ".def()" to use it with python. It should work.

Bruno

bool Ig2_Sphere_ChainedCylinder_CylScGeom::checkOverlap(const
shared_ptr<Shape>& cm1,
                            const shared_ptr<Shape>& cm2,
                            const State& state1, const State& state2)
{
    const State* sphereSt=YADE_CAST<const State*>(&state1);
    const ChainedState* cylinderSt=YADE_CAST<const ChainedState*>(&state2);
    ChainedCylinder *cylinder=YADE_CAST<ChainedCylinder*>(cm2.get());
    Sphere *sphere=YADE_CAST<Sphere*>(cm1.get());
    assert(sphereSt && cylinderSt && cylinder && sphere);
    bool isLast =
(cylinderSt->chains[cylinderSt->chainNumber].size()==(cylinderSt->rank+1));

    shared_ptr<const ChainedState> statePrev;
    if (cylinderSt->rank>0)
        statePrev = YADE_PTR_CAST<const ChainedState>
(Body::byId(cylinderSt->chains[cylinderSt->chainNumber][cylinderSt->rank-1],scene)->state);

    //FIXME : definition of segment in next line breaks periodicity
    shared_ptr<Body> cylinderNext;
    Vector3r segment, branch, direction;
    Real length, dist;
    branch = sphereSt->pos-cylinderSt->pos-shift2;
    if (!isLast) {
        cylinderNext =
Body::byId(cylinderSt->chains[cylinderSt->chainNumber][cylinderSt->rank+1],scene);
        segment = cylinderNext->state->pos-cylinderSt->pos;
        if (segment.dot(branch)>(segment.dot(segment)) return false;
        length = segment.norm();
        direction = segment/length;
        dist = direction.dot(branch);
        if (dist<-interactionDetectionFactor*cylinder->radius &&
            branch.squaredNorm() >
pow(interactionDetectionFactor*(sphere->radius+cylinder->radius),2))
return false;
    } else {//handle the last node with length=0
        segment = Vector3r(0,0,0);
        length = 0;
        direction = Vector3r(0,1,0);
        dist = 0;
        if (branch.squaredNorm() >
interactionDetectionFactor*(sphere->radius+cylinder->radius)) return false;
    }

    //Check sphere-cylinder distance
    Vector3r projectedP = cylinderSt->pos+shift2 + direction*dist;
    branch = projectedP-sphereSt->pos;
    if (branch.squaredNorm()>(pow(sphere->radius+cylinder->radius,2)))
return false;
    return true;
}

Revision history for this message
Kneib François (francois-kneib) said :
#2

OK thanks, I understood but finally I made that in python using the scalar product.
If someone is interested, here is the code :

# encoding: utf-8
from pylab import *

[...]

#Returns the distance between a point and a line in space.
def distPointToLine(A,B,u): #A:the point, #B:a point on the line, #u:vector of the line
 BA=[A[0]-B[0],A[1]-B[1],A[2]-B[2]]
 return norm(cross(BA,u))/norm(u)

#Returns true if the sphere's projection on the cylinder is in the cylinder's segment. Else return false.
def isProjectionInsideSegment(u,v,rs,rc):#
 ratio=inner(u,v)/(norm(u)**2)
 rs=rs/norm(u)
 rc=rc/norm(u)
 if (-rs-rc<ratio<1.0+rs+rc):return True
 else: return False

#Erases all the sphere who are penetrated in any cylinder.
def eraseIntersectedSpheres():
 for i in O.bodies:
  if(type(i.shape).__name__ == 'ChainedCylinder'):
   for j in O.bodies:
    if(type(j.shape).__name__ == 'Sphere'):
     sphereCenter=[j.state.pos[0],j.state.pos[1],j.state.pos[2]]
     sphereRadius=j.shape.radius
     cylinderRadius=i.shape.radius
     cylinderStart=[i.state.pos[0],i.state.pos[1],i.state.pos[2]]
     cylinderSegment=[i.shape.segment[0],i.shape.segment[1],i.shape.segment[2]]
     cylinderSphereSegment=[sphereCenter[0]-cylinderStart[0],sphereCenter[1]-cylinderStart[1],sphereCenter[2]-cylinderStart[2]]
     if (distPointToLine(sphereCenter,cylinderStart,cylinderSegment)<(sphereRadius+cylinderRadius) and isProjectionInsideSegment(cylinderSegment,cylinderSphereSegment,sphereRadius,cylinderRadius)):
      O.bodies.erase(j.id)

[...]

Revision history for this message
Kneib François (francois-kneib) said :
#3

By the way, the main function to be used is the last one : eraseIntersectedSpheres()

Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#4

I have the feeling it doesn't handle intersections between spheres and the nodes connecting cylinders... does it?
You can touch a node between two cylinders while the projections on both cylinders axis are outside the cylinders.
That would explain why your code is simpler.

Revision history for this message
Kneib François (francois-kneib) said :
#5

Yes I think I see the problem, and it's the same if the sphere is just on the extremity of the chained cylinder.
To fix it, I wrote that :
 if (-rs-rc<ratio<1.0+rs+rc):return True
to take care of the sphere and cylinder's radius.
In fact the scalar product have to be inside a larger interval than just 0<ratio<1. I normalized rs and rc :
rs=rs/norm(u)
rc=rc/norm(u)
So I can directly enlarge the interval.
Is this what you were thinking about ?

Revision history for this message
Yade Guide (yade-guide) said :
#6
Revision history for this message
Yade Guide (yade-guide) said :
#7

Hey there! As an automated bot, I've scanned through some relevant threads and compiled a summary for you. Feel free to explore further by clicking on the links below.

Title: "I try to pack the sphere in the cylinder by gravity deposition. But it stop when it touches the funnel."
The script was modified for a smaller sphere, which behaves similarly to a larger sphere when it hits a funnel in the simulation. The unbalanced force condition (unbalancedForce() < 0.01) in the checkUnbalanced() function likely caused the simulation to pause. Printing the value before O.pause() may help resolve this issue.
https://answers.launchpad.net/yade/+question/706797

Title: "O.interactions[0,1] is unable to detect the interactions between sphere and wall, throwing an error stating "No such interaction""
The user is seeking assistance with an issue related to two bodies in a simulation where O.forces.f(0) and O.interactions[0,1] should be the same except that the first will not give an error when there is no interaction. Suggestions include using O.interactions.has(0,1) for clarity.
https://answers.launchpad.net/yade/+question/702412

Title: "extMaterial in cylinderConnection gridpfacet.py not used"
Rohit K. John noticed an inconsistency in the 'cylinderConnection()' function where nodes were assigned 'intMaterial' while connections used 'extMaterial'. Bruno suggested it was likely a typo and could be removed, Rohit questioned whether two materials should be used due to different interactions between nodes and connections.
https://answers.launchpad.net/yade/+question/697249

Revision history for this message
Yade Guide (yade-guide) said :
#8

Hey! As an automated bot, I've scanned through some relevant threads and summarized them below. Feel free to delve deeper into these topics by clicking on the links attached.

Title: "Can not see existed particles in "show 3D" "
Jan provided advice on calculating AABB dimensions for a cylinder shape. Despite using AABB parameters, Jim faced issues with particles exploding during simulation. Jan suggested checking particle visibility before starting the simulation. Leonard thanked Jan and found two solutions: increasing wall's stiffness or reducing O.dt.
https://answers.launchpad.net/yade/+question/683821

Title: "sumForces"
Riccardo is using Yade to simulate a hollow cylinder moving into an aggregate of particles, recording forces with utils.sumForces. He asks if facets not in contact with particles contribute to the forces.
https://answers.launchpad.net/yade/+question/224276