Concept

TriangulationDataStructure_3

Definition

3D-triangulation data structures are meant to maintain the combinatorial information for 3D-geometric triangulations.

In Cgal, a triangulation data structure is a container of cells (3-faces) and vertices (0-faces). Following the standard vocabulary of simplicial complexes, an i-face fi and a j-face fj (0 j < i 3) are said to be incident in the triangulation if fj is a (sub)face of fi, and two i-faces (0 i 3) are said to be adjacent if they share a common incident (sub)face.

Each cell gives access to its four incident vertices and to its four adjacent cells. Each vertex gives direct access to one of its incident cells, which is sufficient to retrieve all the incident cells when needed.

The four vertices of a cell are indexed with 0, 1, 2 and 3. The neighbors of a cell are also indexed with 0, 1, 2, 3 in such a way that the neighbor indexed by i is opposite to the vertex with the same index (see Figure 40.1).

Edges (1-faces) and facets (2-faces) are not explicitly represented: a facet is given by a cell and an index (the facet i of a cell c is the facet of c that is opposite to the vertex of index i) and an edge is given by a cell and two indices (the edge (i,j) of a cell c is the edge whose endpoints are the vertices of indices i and j of c).

As Cgal explicitly deals with all degenerate cases, a 3D-triangulation data structure in Cgal can handle the cases when the dimension of the triangulation is lower than 3 (see Section 40.1).

Thus, a 3D-triangulation data structure can store a triangulation of a topological sphere Sd of d+1, for any d {-1,0,1,2,3}.




The second template parameter of the basic triangulation class (see Chapter 39 ) Triangulation_3 is a triangulation data structure class. (See Chapter 40.)

To ensure all the flexibility of the class Triangulation_3, a model of a triangulation data structure must be templated by the base vertex and the base cell classes (see 40.1): TriangulationDataStructure_3<TriangulationVertexBase_3,TriangulationCellBase_3>. The optional functionalities related to geometry are compulsory for this use as a template parameter of Triangulation_3.




A class that satisfies the requirements for a triangulation data structure class must provide the following types and operations.

Types

TriangulationDataStructure_3::Vertex
Vertex type

TriangulationDataStructure_3::Cell
Cell type


TriangulationDataStructure_3::size_type
Size type (unsigned integral type)

TriangulationDataStructure_3::difference_type
Difference type (signed integral type)

Vertices and cells are usually manipulated via handles, which support the two dereference operators operator* and operator->.

TriangulationDataStructure_3::Vertex_handle
TriangulationDataStructure_3::Cell_handle

Requirements for Vertex and Cell are described in TriangulationDataStructure_3::Vertex and TriangulationDataStructure_3::Cell .

template <typename Vb2>
TriangulationDataStructure_3:: struct Rebind_vertex
This nested template class allows to get the type of a triangulation data structure that only changes the vertex type. It has to define a type Other which is a rebound triangulation data structure, that is, the one whose TriangulationDSVertexBase_3 will be Vb2.

template <typename Cb2>
TriangulationDataStructure_3:: struct Rebind_cell
This nested template class allows to get the type of a triangulation data structure that only changes the cell type. It has to define a type Other which is a rebound triangulation data structure, that is, the one whose TriangulationDSCellBase_3 will be Cb2.

typedef Triple<Cell_handle, int, int>
Edge; (c,i,j) is the edge of cell c whose vertices indices are i and j. (See Section 40.1.)
typedef std::pair<Cell_handle, int>
Facet; (c,i) is the facet of c opposite to the vertex of index i. (See Section 40.1.)

The following iterators allow one to visit all the vertices, edges, facets and cells of the triangulation data structure. They are all bidirectional, non-mutable iterators.

TriangulationDataStructure_3::Cell_iterator
TriangulationDataStructure_3::Facet_iterator
TriangulationDataStructure_3::Edge_iterator
TriangulationDataStructure_3::Vertex_iterator

The following circulators allow us to visit all the cells and facets incident to a given edge. They are bidirectional and non-mutable.

TriangulationDataStructure_3::Facet_circulator
TriangulationDataStructure_3::Cell_circulator

Iterators and circulators are convertible to the corresponding handles, thus the user can pass them directly as arguments to the functions.

Creation

TriangulationDataStructure_3 tds;
Default constructor.


TriangulationDataStructure_3 tds ( tds1);
Copy constructor. All vertices and cells are duplicated.

TriangulationDataStructure_3& tds = tds1 Assignment operator. All vertices and cells are duplicated, and the former data structure of tds is deleted.

Vertex_handle tds.copy_tds ( tds1, Vertex_handle v = Vertex_handle())
tds1 is copied into tds. If v != Vertex_handle(), the vertex of tds corresponding to v is returned, otherwise Vertex_handle() is returned.
Precondition: The optional argument v is a vertex of tds1.

void tds.swap ( & tds1) Swaps tds and tds1. There is no copy of cells and vertices, thus this method runs in constant time. This method should be preferred to tds=tds1 or tds(tds1) when tds1 is deleted after that.
void tds.clear () Deletes all cells and vertices. tds is reset as a triangulation data structure constructed by the default constructor.

Operations

Access Functions

int tds.dimension () const The dimension of the triangulated topological sphere.
size_type tds.number_of_vertices () const The number of vertices. Note that the triangulation data structure has one more vertex than an associated geometric triangulation, if there is one, since the infinite vertex is a standard vertex and is thus also counted.
size_type tds.number_of_cells () const The number of cells. Returns 0 if tds.dimension()<3.

Non constant-time access functions

size_type tds.number_of_facets () const The number of facets. Returns 0 if tds.dimension()<2.
size_type tds.number_of_edges () const The number of edges. Returns 0 if tds.dimension()<1.

Setting

void tds.set_dimension ( int n) Sets the dimension to n.

Queries

bool tds.is_vertex ( Vertex_handle v) const
Tests whether v is a vertex of tds.

bool tds.is_edge ( Cell_handle c, int i, int j) const
Tests whether (c,i,j) is an edge of tds. Answers false when dimension() <1 .
Precondition: i,j {0,1,2,3}
bool tds.is_edge ( Vertex_handle u, Vertex_handle v, Cell_handle & c, int & i, int & j) const
Tests whether (u,v) is an edge of tds. If the edge is found, it computes a cell c having this edge and the indices i and j of the vertices u and v, in this order.
bool tds.is_edge ( Vertex_handle u, Vertex_handle v) const
Tests whether (u,v) is an edge of tds.

bool tds.is_facet ( Cell_handle c, int i) const
Tests whether (c,i) is a facet of tds. Answers false when dimension() <2 .
Precondition: i {0,1,2,3}
bool
tds.is_facet ( Vertex_handle u,
Vertex_handle v,
Vertex_handle w,
Cell_handle & c,
int & i,
int & j,
int & k)
const
Tests whether (u,v,w) is a facet of tds. If the facet is found, it computes a cell c having this facet and the indices i, j and k of the vertices u, v and w, in this order.

bool tds.is_cell ( Cell_handle c) const
Tests whether c is a cell of tds. Answers false when dimension() <3 .

bool
tds.is_cell ( Vertex_handle u,
Vertex_handle v,
Vertex_handle w,
Vertex_handle t,
Cell_handle & c,
int & i,
int & j,
int & k,
int & l)
const
Tests whether (u,v,w,t) is a cell of tds. If the cell c is found, it computes the indices i, j, k and l of the vertices u, v, w and t in c, in this order.

There is a method has_vertex in the cell class. The analogous methods for facets are defined here.

bool tds.has_vertex ( Facet f, Vertex_handle v, int & j) const
If v is a vertex of f, then j is the index of v in the cell f.first, and the method returns true.
Precondition: tds.dimension()=3
bool tds.has_vertex ( Cell_handle c, int i, Vertex_handle v, int & j) const
Same for facet (c,i). Computes the index j of v in c.
bool tds.has_vertex ( Facet f, Vertex_handle v) const
bool tds.has_vertex ( Cell_handle c, int i, Vertex_handle v) const
Same as the first two methods, but these two methods do not return the index of the vertex.

The following three methods test whether two facets have the same vertices.

bool tds.are_equal ( Facet f, Facet g) const
bool tds.are_equal ( Cell_handle c, int i, Cell_handle n, int j) const
bool tds.are_equal ( Facet f, Cell_handle n, int j) const
For these three methods:
Precondition: tds.dimension()=3.

Flips

Two kinds of flips exist for a three-dimensional triangulation. They are reciprocal. To be flipped, an edge must be incident to three tetrahedra. During the flip, these three tetrahedra disappear and two tetrahedra appear. Figure 40.7(left) shows the edge that is flipped as bold dashed, and one of its three incident facets is shaded. On the right, the facet shared by the two new tetrahedra is shaded.

Figure 40.7:  Flips.

Flips

The following methods guarantee the validity of the resulting 3D combinatorial triangulation. Moreover the flip operations do not invalidate the vertex handles, and only invalidate the cell handles of the affected cells.

Flips for a 2d triangulation are not implemented yet

bool tds.flip ( Edge e)
bool tds.flip ( Cell_handle c, int i, int j)
Before flipping, these methods check that edge e=(c,i,j) is flippable (which is quite expensive). They return false or true according to this test.

void tds.flip_flippable ( Edge e)
void tds.flip_flippable ( Cell_handle c, int i, int j)
Should be preferred to the previous methods when the edge is known to be flippable.
Precondition: The edge is flippable.

bool tds.flip ( Facet f)
bool tds.flip ( Cell_handle c, int i) Before flipping, these methods check that facet f=(c,i) is flippable (which is quite expensive). They return false or true according to this test.

void tds.flip_flippable ( Facet f)
void tds.flip_flippable ( Cell_handle c, int i)
Should be preferred to the previous methods when the facet is known to be flippable.
Precondition: The facet is flippable.

Insertions

The following modifier member functions guarantee the combinatorial validity of the resulting triangulation.

Vertex_handle tds.insert_in_cell ( Cell_handle c)
Creates a new vertex, inserts it in cell c and returns its handle. The cell c is split into four new cells, each of these cells being formed by the new vertex and a facet of c.
Precondition: tds.dimension() = 3 and c is a cell of tds.

Vertex_handle tds.insert_in_facet ( Facet f) Creates a new vertex, inserts it in facet f and returns its handle. In dimension 3, the two incident cells are split into 3 new cells; in dimension 2, the facet is split into 3 facets.
Precondition: tds.dimension() 2 and f is a facet of tds.

Vertex_handle tds.insert_in_facet ( Cell_handle c, int i)
Creates a new vertex, inserts it in facet i of c and returns its handle.
Precondition: tds.dimension() 2, i {0,1,2,3} in dimension 3, i=3 in dimension 2 and (c,i) is a facet of tds.

Vertex_handle tds.insert_in_edge ( Edge e) Creates a new vertex, inserts it in edge e and returns its handle. In dimension 3, all the incident cells are split into 2 new cells; in dimension 2, the 2 incident facets are split into 2 new facets; in dimension 1, the edge is split into 2 new edges.
Precondition: tds.dimension() 1 and e is an edge of tds.

Vertex_handle tds.insert_in_edge ( Cell_handle c, int i, int j)
Creates a new vertex, inserts it in edge (i,j) of c and returns its handle.
Precondition: tds.dimension() 1. i j, i,j {0,1,2,3} in dimension 3, i,j {0,1,2} in dimension 2, i,j {0,1} in dimension 1 and (c,i,j) is an edge of tds.

Vertex_handle tds.insert_increase_dimension ( Vertex_handle star = Vertex_handle())
Transforms a triangulation of the sphere Sd of d+1 into the triangulation of the sphere Sd+1 of d+2 by adding a new vertex v: v is linked to all the vertices to triangulate one of the two half-spheres of dimension (d+1). Vertex star is used to triangulate the second half-sphere (when there is an associated geometric triangulation, star is in fact the vertex associated with its infinite vertex). See Figure 40.8.
The numbering of the cells is such that, if f was a face of maximal dimension in the initial triangulation, then (f,v) (in this order) is the corresponding face in the new triangulation. This method can be used to insert the first two vertices in an empty triangulation.
A handle to v is returned.
Precondition: tds.dimension() = d < 3. When tds.number_of_vertices() >0, star Vertex_handle() and star is a vertex of tds.

Figure 40.8:  insert_increase_dimension (1-dimensional case).

insert_increase_dimension} (1-dimensional case)

template <class CellIt>
Vertex_handle tds.insert_in_hole ( CellIt cell_begin, CellIt cell_end, Cell_handle begin, int i)
Creates a new vertex by starring a hole. It takes an iterator range [cell_begin; cell_end[ of Cell_handles which specifies a set of connected cells (resp. facets in dimension 2) describing a hole. (begin, i) is a facet (resp. an edge) on the boundary of the hole, that is, begin belongs to the set of cells (resp. facets) previously described, and begin->neighbor(i) does not. Then this function deletes all the cells (resp. facets) describing the hole, creates a new vertex v, and for each facet (resp. edge) on the boundary of the hole, creates a new cell (resp. facet) with v as vertex. v is returned.
Precondition: tds.dimension() 2, the set of cells (resp. facets) is connected, and its boundary is connected.

template <class CellIt>
Vertex_handle
tds.insert_in_hole ( CellIt cell_begin,
CellIt cell_end,
Cell_handle begin,
int i,
Vertex_handle newv)
Same as above, except that newv will be used as the new vertex, which must have been allocated previously with e.g. create_vertex.

Removal

void tds.remove_decrease_dimension ( Vertex_handle v, Vertex_handle w = v)
This operation is the reciprocal of insert_increase_dimension(). It transforms a triangulation of the sphere Sd of d+1 into the triangulation of the sphere Sd-1 of d by removing the vertex v. Delete the cells incident to w, keep the others.
Precondition: tds.dimension() = d -1. tds.degree(v) = degree(w) = tds.number_of_vertices() -1.

Cell_handle tds.remove_from_maximal_dimension_simplex ( Vertex_handle v)
Removes v. The incident simplices of maximal dimension incident to v are replaced by a single simplex of the same dimension. This operation is exactly the reciprocal to tds.insert_in_cell(v) in dimension 3, tds.insert_in_facet(v) in dimension 2, and tds.insert_in_edge(v) in dimension 1.
Precondition: tds.degree(v) = tds.dimension()+1.

Dimension Manipulation

The following operation, decrease_dimension, is necessary when the displacement of a vertex decreases the dimension of the triangulation.

void tds.decrease_dimension ( Cell_handle c, int i)
The link of a vertex v is formed by the facets disjoint from v that are included in the cells incident to v. When the link of v = c->vertex(i) contains all the other vertices, decrease_dimension crushes the triangulation of the sphere Sd of d+1 onto the triangulation of the sphere Sd-1 of d formed by the link of v augmented with the vertex v itself, for d==2,3; this one is placed on the facet (c, i) (see Fig. 40.9).
Precondition:  The dimension must be 2 or 3. The degree of v must be equal to the total number of vertices of the triangulation data structure minus 1.

Figure 40.9:  From an Sd data structure to an Sd-1 data structure (top: d==2, bottom: d==3).

Lowering dimension from 3D to 2D

Other modifiers

The following modifiers can affect the validity of the triangulation data structure.

void tds.reorient () Changes the orientation of all cells of the triangulation data structure.
Precondition: tds.dimension() 1.

Vertex_handle tds.create_vertex ( Vertex v = Vertex())
Adds a copy of the vertex v to the triangulation data structure.

Vertex_handle tds.create_vertex ( Vertex_handle v)
Creates a vertex which is a copy of the one pointed to by v and adds it to the triangulation data structure.

Cell_handle tds.create_cell ( Cell c = Cell())
Adds a copy of the cell c to the triangulation data structure.

Cell_handle tds.create_cell ( Cell_handle c) Creates a cell which is a copy of the one pointed to by c and adds it to the triangulation data structure.

Cell_handle tds.create_cell ( Vertex_handle v0, Vertex_handle v1, Vertex_handle v2, Vertex_handle v3)
Creates a cell and adds it into the triangulation data structure. Initializes the vertices of the cell, its neighbor handles being initialized with the default constructed handle.

Cell_handle
tds.create_cell ( Vertex_handle v0,
Vertex_handle v1,
Vertex_handle v2,
Vertex_handle v3,
Cell_handle n0,
Cell_handle n1,
Cell_handle n2,
Cell_handle n3)
Creates a cell, initializes its vertices and neighbors, and adds it into the triangulation data structure.

void tds.delete_vertex ( Vertex_handle v)
Removes the vertex from the triangulation data structure.
Precondition: The vertex is a vertex of tds.

void tds.delete_cell ( Cell_handle c) Removes the cell from the triangulation data structure.
Precondition: The cell is a cell of tds.

template <class VertexIt>
void tds.delete_vertices ( VertexIt first, VertexIt last)
Calls delete_vertex over an iterator range of value type Vertex_handle.

template <class CellIt>
void tds.delete_cells ( CellIt first, CellIt last)
Calls delete_cell over an iterator range of value type Cell_handle.

Traversing the triangulation

Cell_iterator tds.cells_begin () const Returns cells_end() when tds.dimension() <3.
Cell_iterator tds.cells_end () const
Cell_iterator tds.raw_cells_begin () const Low-level access to the cells, does not return cells_end() when tds.dimension() <3.
Cell_iterator tds.raw_cells_end () const
Facet_iterator tds.facets_begin () const Returns facets_end() when tds.dimension() <2.
Facet_iterator tds.facets_end () const
Edge_iterator tds.edges_begin () const Returns edges_end() when tds.dimension() <1.
Edge_iterator tds.edges_end () const
Vertex_iterator tds.vertices_begin () const
Vertex_iterator tds.vertices_end () const

Cell and facet circulators

Cell_circulator tds.incident_cells ( Edge e) const
Starts at an arbitrary cell incident to e.
Precondition: tds.dimension() =3
Cell_circulator tds.incident_cells ( Cell_handle c, int i, int j) const
As above for edge (i,j) of c.
Cell_circulator tds.incident_cells ( Edge e, Cell_handle start) const
Starts at cell start.
Precondition: tds.dimension() =3 and start is incident to e.
Cell_circulator tds.incident_cells ( Cell_handle c, int i, int j, Cell_handle start) const
As above for edge (i,j) of c.

The following circulators on facets are defined only in dimension 3, though facets are defined also in dimension 2: there are only two facets sharing an edge in dimension 2.

Facet_circulator tds.incident_facets ( Edge e) const
Starts at an arbitrary facet incident to e.
Precondition: tds.dimension() =3
Facet_circulator tds.incident_facets ( Cell_handle c, int i, int j) const
As above for edge (i,j) of c.
Facet_circulator tds.incident_facets ( Edge e, Facet start) const
Starts at facet start.
Precondition: start is incident to e.
Facet_circulator tds.incident_facets ( Edge e, Cell_handle start, int f) const
Starts at facet of index f in start.
Facet_circulator tds.incident_facets ( Cell_handle c, int i, int j, Facet start) const
As above for edge (i,j) of c.
Facet_circulator tds.incident_facets ( Cell_handle c, int i, int j, Cell_handle start, int f) const
As above for edge (i,j) of c and facet (start,f).

Traversal of the incident cells, facets and edges, and the adjacent vertices of a given vertex

template <class OutputIterator>
OutputIterator tds.incident_cells ( Vertex_handle v, OutputIterator cells) const
Copies the Cell_handles of all cells incident to v to the output iterator cells. Returns the resulting output iterator.
Precondition: tds.dimension() =3, v Vertex_handle(), tds.is_vertex(v).

template <class OutputIterator>
OutputIterator tds.incident_facets ( Vertex_handle v, OutputIterator facets) const
Copies the Facets incident to v to the output iterator facets. Returns the resulting output iterator.
Precondition: tds.dimension() >1, v Vertex_handle(), tds.is_vertex(v).

template <class OutputIterator>
OutputIterator tds.incident_edges ( Vertex_handle v, OutputIterator edges) const
Copies all Edges incident to v to the output iterator edges. Returns the resulting output iterator.
Precondition: tds.dimension() >0, v Vertex_handle(), tds.is_vertex(v).

template <class OutputIterator>
OutputIterator tds.adjacent_vertices ( Vertex_handle v, OutputIterator vertices) const
Copies the Vertex_handles of all vertices adjacent to v to the output iterator vertices. If tds.dimension() <0, then do nothing. Returns the resulting output iterator.
Precondition: v Vertex_handle(), tds.is_vertex(v).

size_type tds.degree ( Vertex_handle v) const
Returns the degree of a vertex, that is, the number of incident vertices.
Precondition: v Vertex_handle(), tds.is_vertex(v).

Traversal between adjacent cells

int tds.mirror_index ( Cell_handle c, int i) const
Returns the index of c in its ith neighbor.
Precondition: i {0, 1, 2, 3}.
Vertex_handle tds.mirror_vertex ( Cell_handle c, int i) const
Returns the vertex of the ith neighbor of c that is opposite to c.
Precondition: i {0, 1, 2, 3}.
Facet tds.mirror_facet ( Facet f) const Returns the same facet seen from the other adjacent cell.

Checking

bool tds.is_valid ( bool verbose = false) const
Checks the combinatorial validity of the triangulation by checking the local validity of all its cells and vertices (see functions below). (See Section 40.1.) Moreover, the Euler relation is tested.
When verbose is set to true, messages are printed to give a precise indication on the kind of invalidity encountered.

bool tds.is_valid ( Vertex_handle v, bool verbose = false) const
Checks the local validity of the adjacency relations of the triangulation. It also calls the is_valid member function of the vertex. When verbose is set to true, messages are printed to give a precise indication on the kind of invalidity encountered.

bool tds.is_valid ( Cell_handle c, bool verbose = false) const
Checks the local validity of the adjacency relations of the triangulation. It also calls the is_valid member function of the cell. When verbose is set to true, messages are printed to give a precise indication on the kind of invalidity encountered.

I/O

istream& istream& is >> & tds Reads a combinatorial triangulation from is and assigns it to tds

ostream& ostream& os << tds Writes tds into the stream os

The information stored in the iostream is: the dimension, the number of vertices, the number of cells, the indices of the vertices of each cell, then the indices of the neighbors of each cell, where the index corresponds to the preceding list of cells. When dimension < 3, the same information is stored for faces of maximal dimension instead of cells.

Has Models

CGAL::Triangulation_data_structure_3

See Also

TriangulationDataStructure_3::Vertex
TriangulationDataStructure_3::Cell