Find simplex that contains point in triangulation

Suppose I have a triangulation on a domain in \mathbb{R}^n, where the element simplexes are given by their vertices.

What is the recommended approach for finding which simplex contains a given point? Naively, I would just test each with something like

function is_point_in_simplex(x, zs) # unoptimized, no checks, etc
    w = vcat(reduce(hcat, zs), ones(1, length(zs))) \ vcat(x, 1)
    all(w .≥ 0)

S = [[0, 0], [0, 1], [1, 0]]
is_point_in_simplex([0.3, 0.3], S) # true
is_point_in_simplex([1, 1], S) # false

but I am not sure if this is the right approach.

Binning strategies like tree-based domain decomposition could reduce it further, but I don’t think there’s a much better approach (at least from a finite element toolbox).

This isn’t done often in FEM because relevant points in a simplex are usually mapped from a reference simplex explicitly.


The package VoronoiDelaunay.jl has a method locate whose source implementation seems quite similar to your example code, which seems like a reasonable indication that your approach is workable.

EDIT: looking at the code for a minute, their findindex function is a custom function (at a glance I assumed it was something in Base). It is slightly more clever and looks to use a tree-like method to search more efficiently. So a bit more thoughtful than a brute-force search.

1 Like

There is also the Jump and Walk Algorithm which does not need any additional data structures beyond the triangulation itself.