In the context of 2 meshes colliding/intersecting (contact mechanics), what is the best way to get the colliding paired-nodes?

General
1

Hi,

I am dealing with 2 dynamic submeshes contacting one to each other. I would like to get the paired nodes that are colliding.
To do so, I am trying to use BoundingBoxTree and in particular its method compute_collisions(tree) (see code below) but I do not understand how to properly use it.

Here are my questions:

  • What exactly does BoundingBoxTree.compute_collisions(tree) method return? I thought it was two lists: list of DOFs in mesh 1 and list of DOFs in mesh 2 that are colliding between each other, and classified by colliding pairs (i.e. collinding_pair_1 = list1[0], list2[0]). Is it? When I plot the pairs, it does not seem to be so and the pairs are distant?

  • What is the best way to get colliding paired-nodes?

Many thanks for your help,
Anne

dimension = 3
V = fenics.VectorFunctionSpace(mesh, "CG", 1)

mesh1 = fenics.SubMesh(mesh, sudomains, 1)
mesh2 = fenics.SubMesh(mesh, sudomains, 2)

# build mesh1 & mesh2 trees
bbtree1 = fenics.BoundingBoxTree()
bbtree1.build(mesh1)

bbtree2= fenics.BoundingBoxTree()
bbtree2.build(mesh2)

# compute colliding entities between 2 trees
colliding_entities_DOFs_1, colliding_entities_DOFs_2= bbtree1.compute_collisions(bbtree2)

# Get node_index to dofs correspondance map
nodeindex_to_dofs = fenics.vertex_to_dof_map(V)
nodeindex_to_dofs = nodeindex_to_dofs.reshape((-1, dimension)) 

# Convert DOF into node_index in colliding_entities_DOFs_1 
colliding_nodes_indices_1 = np.zeros( (len(colliding_entities_DOFs_1)), dtype=int ) 
for i, dof in enumerate(np.array(colliding_entities_DOFs_1)):
     if len(np.where(nodeindex_to_dofs[:,0] == dof)[0]) != 0: # find node_index value for which DOFs contain dof (get array with node index value)
          colliding_nodes_indices_1[i] = np.where(nodeindex_to_dofs[:,0] == dof)[0][0] # get node_index value

     elif len(np.where(nodeindex_to_dofs[:,1] == dof)[0]) != 0: 
          colliding_nodes_indices_1[i] = np.where(nodeindex_to_dofs[:,1] == dof)[0][0] 
     else:
          colliding_nodes_indices_1[i] = np.where(nodeindex_to_dofs[:,2] == dof)[0][0]
 
# Same conversion for colliding_entities_DOFs_2 
colliding_nodes_indices_2 = np.zeros( (len(colliding_entities_DOFs_2)), dtype=int ) 
for i, dof in enumerate(np.array(colliding_entities_DOFs_2)):
     if len(np.where(nodeindex_to_dofs[:,0] == dof)[0]) != 0: 
          colliding_nodes_indices_2[i] = np.where(nodeindex_to_dofs[:,0] == dof)[0][0] 

     elif len(np.where(nodeindex_to_dofs[:,2] == dof)[0]) != 0: 
          colliding_nodes_indices_2[i] = np.where(nodeindex_to_dofs[:,1] == dof)[0][0] 
     else:
          colliding_nodes_indices_2[i] = np.where(nodeindex_to_dofs[:,2] == dof)[0][0]

# Get the colliding-nodes pairs indices
pairs_colliding_nodes = np.zeros((len(colliding_entities_DOFs_1), 2), dtype=int)
pairs_colliding_nodes[:,0], pairs_colliding_nodes[:,1] = colliding_nodes_indices_1[:], colliding_nodes_indices_2[:]
2

I would suggest reading: Bitbucket
as it explains it in quite some detail.

3

Many thanks for this help @dokken. Reading that detail makes me think of using compute_entity_collisions(bbtree) rather than compute_collisions(bbtree).
Am I correct?

Many thanks again,
Anne

4

Yes, if you want to know what cells/facets collide/intersect with eachother, this is what you should use.

5

Ok, thank you.
Actually, I would like to know nodes that collide with eachother. → Which method may I use?

Another additional question to better understand the output result of methods:
Using I get lists containing dofs, like:

colliding_entities_DOFs_1, colliding_entities_DOFs_2 =
[13372, 13372, 13372, 13372, 13372, 13372], [60062, 21788, 59070, 59069, 7852, 48755]

→ Does 1st dof ‘13372’ refer to a cell/facet or a node in case of compute_entity_collisions(bbtree) ? And in case of compute_collisions(bbtree)?

Many thanks again,
Anne

6

Nodes are points, they rarely collide (only when they are at the same coordinate).
You can get the vertices of each cell that is colliding:

from IPython import embed
import fenics


class Domain1(fenics.SubDomain):
    def inside(self, x, on_boundary):
        return x[0] <= 0.5


mesh = fenics.UnitSquareMesh(10, 10)


subdomains = fenics.MeshFunction("size_t", mesh, mesh.topology().dim(), 1)
Domain1().mark(subdomains, 2)
mesh1 = fenics.MeshView.create(subdomains, 1)
mesh2 = fenics.MeshView.create(subdomains, 2)

# build mesh1 & mesh2 trees
bbtree1 = fenics.BoundingBoxTree()
bbtree1.build(mesh1, mesh1.topology().dim())

bbtree2 = fenics.BoundingBoxTree()
bbtree2.build(mesh2, mesh2.topology().dim())

cells_1, cells_2 = bbtree1.compute_entity_collisions(bbtree2)

c_to_v1 = mesh1.topology()(mesh1.topology().dim(), 0)
c_to_v2 = mesh2.topology()(mesh2.topology().dim(), 0)

for i, (cell1, cell2) in enumerate(zip(cells_1, cells_2)):
    vertices1 = c_to_v1(cell1)
    vertices2 = c_to_v2(cell2)
    print("\n"+10*"-"+"\n")
    print(mesh1.coordinates()[vertices1])
    print(mesh2.coordinates()[vertices2])

For these vertices you could do matching vertex to vertex and use the vertex to dof map

7

Many thanks @dokken!