How to get global vertex index of a point starting from coordinates

mesh
1

Hi everybody, I’m using the c++ interface, and I’m trying to get the global index of a vertex starting from the corresponing coordinate inside the mesh.
I found a way, but it is very unefficient, for completeness I report it down below.

auto Mesh = std::make_shared<Mesh>(BoxMesh::create({p0,p1}, {{nx,ny,nz}}, CellType::Type::tetrahedron)); 
dolfin::Point pp(x[0], x[1], x[2]);
std::shared_ptr<dolfin::BoundingBoxTree> tree = Mesh->bounding_box_tree();
auto iTet = tree->compute_first_entity_collision(pp);
	
ufc::cell ufc_cell; 
dolfin::Cell cell(*Mesh, iTet); 
cell.get_cell_topology(ufc_cell); 
std::vector<size_t> vertex_idxs = ufc_cell.entity_indices[0]; 
	
int check=0;
int index=0;
for(size_t i=0; i<4; i++)
	{
		dolfin::MeshGeometry m_geom = Th->geometry(); 
		dolfin::Point pp2 = m_geom.point(vertex_idxs[i]); //trovo coord dato indice
		if (!check && pp2.x() == pp.x() && pp2.y() == pp.y() && pp2.z() == pp.z())
                    {
                    index = vertex_idxs[i];
			        check=1;
                    }
	}

Does somebody know a more compact and efficient way to do the same thing?
I need to know because I’m doing this operation a very large number of times and the computation takes very much.
Thanks to everybody is going to help me!

2

This is an inherently-costly operation for a general unstructured mesh. I would first consider: Is it really necessary? If you have the exact spatial coordinates of a mesh vertex, then presumably they were originally obtained from the Mesh object and would have been accessed via an index. You might be able to re-structure the rest of your code to retain that index alongside the coordinates.

If, instead, you are relying on the fact that BoxMesh is uniform and structured to get the spatial coordinates of vertices, you may be able to reverse-engineer the indexing scheme for the vertices by looking at small meshes (or by inspecting the code for BoxMesh::create), to get an O(1) mapping between coordinates and vertex indices. However, that is somewhat of a hack.

1 Like
4

The mesh I am using is structured and uniform.
I need to do that because I’m assigning some grey values of an image to the rhs of my variational form and I have those values collected into a vector ordered like the global indeces of my vertex (to summarize: in each vertex I need to put the corresponding grey value).
So all the proceedure is to be done inside the eval method of the dolfin Expression that corresponds to my rhs.
I don’t know if this is clear enough.

An alternative way in which I did that before was creating a map to recover the position inside the vector starting from the coordinate value, but sice the coordinate values are real numbers and not integer I think that there might be some roundness problem in that way.
That’s why I was searching for an alternative using the dolfin functions.

5

Would it work to instead represent the grey values as a Function, f_grey, in the space V_grey = FunctionSpace(mesh,"CG",1)? If you already have grey values ordered by vertex index, you could use vertex_to_dof_map(V_grey) (or its C++ API equivalent) to fill in the entries of f_grey.vector(). Then f_grey could be evaluated in the usual way during assembly of the right-hand side vector.

1 Like
6

It worked!
I just added one special class to pass from f_grey Function to an f Expression, otherwise I can’t use it in the variational form
Thank you so much!