1 #include <CGAL/Triangulation_data_structure.h>
2 #include <CGAL/Combination_enumerator.h>
3 #include <CGAL/assertions.h>
4
5 #include <iostream>
6 #include <iterator>
7 #include <vector>
8
9 template< typename TDS >
10 void find_face_from_vertices(const TDS & tds,
11 const std::vector<typename TDS::Vertex_handle> & face_vertices,
12 typename TDS::Face & face);
13
14 template< typename TDS >
barycentric_subdivide(TDS & tds,typename TDS::Full_cell_handle fc)15 void barycentric_subdivide(TDS & tds, typename TDS::Full_cell_handle fc)
16 { /* This function builds the barycentric subdivision of a single
17 full cell |fc| from a triangulation data structure |tds|. */
18 typedef typename TDS::Vertex_handle Vertex_handle;
19 typedef typename TDS::Face Face;
20 const int dim = tds.current_dimension();
21
22 // First, read handles to the cell's vertices
23 std::vector<Vertex_handle> vertices;
24 std::vector<Vertex_handle> face_vertices;
25 for( int i = 0; i <= dim; ++i ) vertices.push_back(fc->vertex(i));
26
27 // Then, subdivide the cell |fc| once by inserting one vertex
28 tds.insert_in_full_cell(fc);
29 // From now on, we can't use the variable |fc|...
30
31 // Then, subdivide faces of |fc| in order of decreasing dimension
32 for( int d = dim-1; d > 0; --d )
33 {
34 face_vertices.resize(d+1);
35 // The following class
36 // enumerates all (d+1)-tuple of the set {0, 1, ..., dim}
37 CGAL::Combination_enumerator<unsigned int> combi(d+1, 0, dim);
38 while ( !combi.finished() )
39 {
40 for( int i = 0; i <= d; ++i )
41 face_vertices[i] = vertices[combi[i]];
42 // we need to find a face with face_vertices
43 Face face(dim);
44 find_face_from_vertices(tds, face_vertices, face);
45 tds.insert_in_face(face);
46 ++combi;
47 }
48 }
49 }
50
51 template< typename TDS >
find_face_from_vertices(const TDS & tds,const std::vector<typename TDS::Vertex_handle> & face_vertices,typename TDS::Face & face)52 void find_face_from_vertices( const TDS & tds,
53 const std::vector<typename TDS::Vertex_handle> & face_vertices,
54 typename TDS::Face & face)
55 { /* The main goal of this function is to find a full cell that
56 contains a given set of vertices |face_vertices|. Then, it
57 builds a corresponding |face|. */
58 typedef typename TDS::Vertex_handle Vertex_handle;
59 typedef typename TDS::Full_cell_handle Full_cell_handle;
60 typedef typename TDS::Full_cell::Vertex_handle_iterator Vertex_h_iterator;
61
62 // get the dimension of the face we want to build
63 std::size_t fdim(face_vertices.size() - 1);
64 if( fdim <= 0) exit(-1);
65
66 // find all full cells incident to the first vertex of |face|
67 typedef std::vector<Full_cell_handle> Cells;
68 Cells cells;
69 std::back_insert_iterator<Cells> out(cells);
70 tds.incident_full_cells(face_vertices[0], out);
71 // Iterate over the cells to find one which contains the face_vertices
72 for( typename Cells::iterator cit = cells.begin(); cit != cells.end(); ++cit){
73 // find if the cell *cit contains the Face |face|
74 std::size_t i = 0;
75 for( ; i <= fdim; ++i ) {
76 Vertex_handle face_v(face_vertices[i]);
77 bool found(false);
78 Vertex_h_iterator vit = (*cit)->vertices_begin();
79 for( ; vit != (*cit)->vertices_end(); ++vit ) {
80 if( *vit == face_v ) {
81 found = true;
82 break;
83 }
84 }
85 if( ! found )
86 break;
87 }
88 if( i > fdim ) {// the full cell |*cit| contains |face|
89 face.set_full_cell(*cit);
90 for( std::size_t i = 0; i <= fdim; ++i )
91 {
92 face.set_index(static_cast<int>(i),
93 (*cit)->index(face_vertices[i]));
94 }
95 return;
96 }
97 }
98 std::cerr << "Could not build a face from vertices"<<std::endl;
99 CGAL_assertion(false);
100 }
101
102
main()103 int main()
104 {
105 const int sdim = 5; // dimension of TDS with compile-time dimension
106 typedef CGAL::Triangulation_data_structure<CGAL::Dimension_tag<sdim> >
107 TDS;
108 TDS tds(sdim);
109
110 TDS::Vertex_handle one_vertex = tds.insert_increase_dimension();
111 for( int i = 1; i < sdim+2; ++i )
112 tds.insert_increase_dimension(one_vertex);
113 // we get a triangulation of space of dim sdim homeomorphic to
114 // the boundary of simplex of dimension sdim+1 with sdim+2 vertices
115 CGAL_assertion( sdim == tds.current_dimension() );
116 CGAL_assertion( 2+sdim == tds.number_of_vertices() );
117 CGAL_assertion( 2+sdim == tds.number_of_full_cells() );
118
119 barycentric_subdivide(tds, tds.full_cells_begin());
120
121 // The number of full cells should be twice the factorial of
122 // |tds.current_dimension()+1|. Eg, 1440 for dimension 5.
123 std::cout << "Triangulation has "
124 << tds.number_of_full_cells() << " full cells";
125 CGAL_assertion( tds.is_valid() );
126 std::cout << " and is valid!"<<std::endl;
127 return 0;
128 }
129