How to determine whether a point belongs to a mesh cell

1

Hi,

This might be a very simple question but can’t find an easy solution to it. Suppose I have a point which lies strictly inside of a mesh cell.
1- What is the fastest way to extract the ID and vertices of that cell?
2- What if that point is exactly a mesh point? Which cell ID will I get as the answer?

Here is a minimal code with a mesh and the x,y,z coordinates of a point which I am sure lies strictly inside of a cell. How can I answer the above questions for this code?

from dolfin import *

meshlevel =10
degree = 1
dim = 2
mesh = UnitDiscMesh.create(MPI.comm_world,meshlevel, degree, dim)

Px = -0.519
Py = 0.33
Pz = 0.0
2

The following code illustrates all your questions:

from dolfin import *

meshlevel =10
degree = 1
dim = 2
mesh = UnitDiscMesh.create(MPI.comm_world,meshlevel, degree, dim)

Px = -0.519
Py = 0.33
Pz = 0.0

p = Point(Px, Py, Pz)
tree = BoundingBoxTree()
tree.build(mesh)
first = tree.compute_first_entity_collision(p)
print(first)
c_to_v = mesh.topology()(mesh.topology().dim(), 0)
print(c_to_v(first), mesh.coordinates()[c_to_v(first)])


p0 = Point(mesh.coordinates()[0])


all_cells = tree.compute_entity_collisions(p0)
print(mesh.coordinates()[0])
for cell in all_cells:
    print(cell, c_to_v(cell), mesh.coordinates()[c_to_v(cell)])

You can choose to either get the first entity collision (as traversed through the boundingboxtree) or get all of them and then choose one of them.

3

Wow. This is exactly what I was looking for. So the key is BoundingBoxTree. Thank you so much.