Out of range while using compute_closest_entity

1

Consider a circular mesh with associated function space:

mesh_sphere = Mesh()
with XDMFFile(mesh_path + "circle_mesh.xdmf") as infile:
    infile.read(mesh_sphere)

L1_sph = FiniteElement("Lagrange", mesh_sphere.ufl_cell(), 1)
V_sph = FunctionSpace(mesh_sphere, L1_sph)

Now, launching

tree_sph = mesh_sphere.bounding_box_tree()
tree_sph.compute_closest_entity(Point(.5,.86602541))[0]

returns 4294967295. You can find the mesh here: Circle mesh - Google Drive.

What can I do to overcome this? My intent is to find the closest cell to that point, and for other points, even very far away from the circle, this worked well.

Also for the point e.g. Point(.5,.86602542)

2

This number

is np.uintc(-1)

which I would interpret as the algorithm returning its input:
https://bitbucket.org/fenics-project/dolfin/src/74a6ac2095059ea6582fe4ea11686c1c723e14f5/dolfin/geometry/GenericBoundingBoxTree.cpp?at=master#lines-282
which means that the point found in Bitbucket
is as close as any point in the mesh (i.e. a vertex/mipoint is the closest entity).

This seems like a bug to me, and you could consider adding an issue on fenics-project / DOLFIN / issues — Bitbucket

However, as DOLFINx is now the recommened version to use, ref: The new DOLFINx solver is now recommended over DOLFIN
it is unlikely that it will be fixed in the old version.

For DOLFINx, the syntax would be:

from IPython import embed
import dolfinx
import numpy as np
import ufl
from mpi4py import MPI

mesh_path = "./"

with dolfinx.io.XDMFFile(MPI.COMM_WORLD, mesh_path + "circle_mesh.xdmf", "r") as infile:
    mesh_sphere = infile.read_mesh(name="Grid")

L1_sph = ufl.FiniteElement("Lagrange", mesh_sphere.ufl_cell(), 1)
V_sph = dolfinx.fem.FunctionSpace(mesh_sphere, L1_sph)

dim = mesh_sphere.topology.dim
tree_sph = dolfinx.geometry.BoundingBoxTree(mesh_sphere, dim)

num_entities_local = mesh_sphere.topology.index_map(
    dim).size_local + mesh_sphere.topology.index_map(dim).num_ghosts
entities = np.arange(num_entities_local, dtype=np.int32)
midpoint_tree = dolfinx.geometry.create_midpoint_tree(
    mesh_sphere, dim, entities)

points = np.array([[.5, .86602541, 0]])

closest_entities = dolfinx.geometry.compute_closest_entity(
    tree_sph, midpoint_tree, mesh_sphere, points)
squared_distance = dolfinx.geometry.squared_distance(
    mesh_sphere, dim, closest_entities, points)
print(closest_entities, squared_distance)


coords = mesh_sphere.geometry.x
x_dofmap = mesh_sphere.geometry.dofmap
coordinate_dofs = coords[x_dofmap.links(closest_entities[0])]
distance_vec = dolfinx.geometry.compute_distance_gjk(
    coordinate_dofs, points[0])
print(f"Point {points[0]} Closest point on mesh {points[0] + distance_vec}")

which yields you the following image:

distance
distance1461×1037 15.1 KB

2 Likes