Find which tetrahedron in a Triangulation contains a Given Facet

3 views (last 30 days)
I have a Delaunay tetrahedralization, defined by a Coordinate Vector "X" (of size n x 3) and tetrahedra IDs "TRI"(of size m x 4). I want to know which tetrahedron in the triangulation contains a given triangular facet, based on IDs of the 3 vertices. MWE:
TRI = [1 2 3 4;5 6 7 8; 9 10 11 12; 13 14 15 16];
Facet = [5 7 8];
and I am obviously looking for
ans = 2
but since the number of columns of TRI and Facet are not the same, I can't use the ismember function. Any other ideas? Thanks.

Accepted Answer

John D'Errico
John D'Errico on 1 Mar 2017
Edited: John D'Errico on 1 Mar 2017
Um, can you say easy, peasy?
1. Sort the facet indices in your target facet.
2. Sort the facet indices in your triangulation.
3. Form a list of all facets. Each tetrahedron has 4 facets, so there will be 4 facets for each simplex in the complex.
FACET = sort(FACET);
TRI = sort(TRI,2);
facetlist = [TRI(:,[1 2 3]);TRI(:,[1 2 4]);TRI(:,[1 3 4]);TRI(:,[2 3 4])];
Now a simple application of ismember will find the facet of interest. Note that the facet may appear once or twice. This depends upon whether the simplex shares that facet with another simplex in the tessellation.
How do you identify which simplex the facets you have just identified came from? One simple solution would be ind2sub, but mod will also suffice.
For example:
tri = delaunayn(rand(7,3))
tri =
3 1 6 5
2 6 1 5
7 3 4 5
7 4 2 5
7 1 3 5
7 2 1 5
tri = sort(tri,2)
tri =
1 3 5 6
1 2 5 6
3 4 5 7
2 4 5 7
1 3 5 7
1 2 5 7
facet = [1 5 7];
facetlist = [tri(:,[1 2 3]);tri(:,[1 2 4]);tri(:,[1 3 4]);tri(:,[2 3 4])];
[~,ind] = ismember(facetlist,facet,'rows');
ind = find(ind)
ind =
17
18
Now, which simplex did each of those facets arise from?
ntri = size(tri,1);
mod(ind-1,ntri) + 1
ans =
5
6
That particular facet was shared between simplex 5 and simplex 6 in the list. If you look back at the original triangulation, it did appear in each of those places. As I said, if a facet is on the boundary of the tessellation, then it will appear only once.
Now you should see why the sorts were necessary, as otherwise the ismember search may well fail.
As I said, easy. A few simple lines of code.

More Answers (0)

Categories

Find more on Delaunay Triangulation in Help Center and File Exchange

Products

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!