In bench6.cpp add indices to the vertices, and report the edges on the convex hull

Triangulation/benchmark/Triangulation/bench6d.cpp

@@ -9,7 +9,7 @@
 #include <vector>
 
 typedef CGAL::Epick_d< CGAL::Dimension_tag<6> >  K;
-typedef CGAL::Triangulation_vertex<K>                                                     Vertex;
+typedef CGAL::Triangulation_vertex<K, std::ptrdiff_t>                                     Vertex;
 typedef CGAL::Triangulation_ds_full_cell<void,CGAL::TDS_full_cell_mirror_storage_policy>  DS_full_cell;
 typedef CGAL::Triangulation_full_cell<K,CGAL::No_full_cell_data, DS_full_cell>            Full_cell;
 typedef CGAL::Triangulation_data_structure<CGAL::Dimension_tag<6>,
@@ -18,12 +18,12 @@ typedef CGAL::Triangulation_data_structure<CGAL::Dimension_tag<6>,
 typedef CGAL::Triangulation<K,TDS>                                                        Triangulation;
 
 
-int main()
+int main(int argc, char* argv[])
 {
     const int D = 6;   // we work in Euclidean 6-space
 
     std::vector<Triangulation::Point> points;
-    std::ifstream in("points_6.txt");
+    std::ifstream in(argv[1]);
     Triangulation::Point p;
     int d;
     in >> d;
@@ -41,20 +41,42 @@ int main()
     assert(t.empty());
 
 
-    t.insert(points.begin(), points.end());  // compute triangulation
+    t.insert_and_index(points.begin(), points.end());  // compute triangulation
 
     std::cout << "Time taken to insert "<< points.size() << " points: " << timer.time() << " seconds." << std::endl;
     timer.reset();
     assert( t.is_valid() );
-    // - - - - - - - - - - - - - - - - - - - - - - - - STEP 2
-    typedef Triangulation::Face Face;
-    typedef std::vector<Face> Faces;
-    Faces edges;
-    std::back_insert_iterator<Faces> out(edges);
-    t.tds().incident_faces(t.infinite_vertex(), 1, out);
-    // collect faces of dimension 1 (edges) incident to the infinite vertex
-    std::cout << "There are " << edges.size()
-              << " vertices on the convex hull." << std::endl;
+
+    std::vector<Triangulation::Full_cell_handle> infinite_cells;
+    t.tds().incident_full_cells(t.infinite_vertex(), std::back_inserter(infinite_cells));
+
+    std::set<std::pair<std::ptrdiff_t, std::ptrdiff_t>> edges;
+
+    for(auto ch : infinite_cells) {
+        std::vector<Triangulation::Vertex_handle> vertices;
+        std::cout << "cell" << std::endl;
+        vertices.reserve(D);
+        for(int i = 0; i < D + 1; ++i) {
+            if(ch->vertex(i) != t.infinite_vertex()) {
+                vertices.push_back(ch->vertex(i));
+                std::cout << ch->vertex(i)->data() << " ";
+            }
+        }
+        std::cout << std::endl;
+        for (int i = 0; i < D-1; i++) {
+          for (int j = 0; j < D-1; j++) {
+
+            if((i != j) && (vertices[i]->data() < vertices[j]->data())) {
+                edges.insert(std::make_pair(vertices[i]->data(), vertices[j]->data()));
+                }
+          }
+        }
+    }
+    for(auto pair : edges) {
+        std::cout << "edge  between vertices " << pair.first << " and " << pair.second << std::endl;
+    }
+
     std::cout << "Time taken to find edges: " << timer.time() << " seconds." << std::endl;
+
      return 0;
 }

Triangulation/include/CGAL/Triangulation.h

@@ -22,6 +22,8 @@
 #include <CGAL/Triangulation_vertex.h>
 #include <CGAL/Iterator_project.h>
 #include <CGAL/spatial_sort.h>
+#include <CGAL/Spatial_sort_traits_adapter_d.h>
+#include <CGAL/property_map.h>  // for CGAL::Identity_property_map
 #include <CGAL/Dimension.h>
 #include <CGAL/iterator.h>
 #include <CGAL/Default.h>
@@ -683,6 +685,36 @@ public:
         }
         return number_of_vertices() - n;
     }
+
+
+    template< typename RandomAccessIterator >
+    size_type insert_and_index(RandomAccessIterator start, RandomAccessIterator end)
+    {
+        size_type n = number_of_vertices();
+        std::vector<std::ptrdiff_t> indices(boost::counting_iterator<std::ptrdiff_t>(0),
+                                            boost::counting_iterator<std::ptrdiff_t>(end - start));
+
+        using Point_property_map = boost::iterator_property_map<RandomAccessIterator,
+                                                                CGAL::Identity_property_map<std::ptrdiff_t>>;
+        using Search_traits_d = CGAL::Spatial_sort_traits_adapter_d<Geom_traits, Point_property_map>;
+
+        CGAL::spatial_sort(indices.begin(), indices.end(), Search_traits_d(start));
+        std::vector<Point> points(start, end);
+        spatial_sort(points.begin(), points.end(), geom_traits());
+        Full_cell_handle hint = Full_cell_handle();
+
+      for (auto index : indices) {
+        typename Triangulation::Vertex_handle pos = insert(*(start+index), hint);
+        if (pos != nullptr) {
+          // Save index value as data to retrieve it after insertion
+          pos->data() = n+index;
+          hint = pos->full_cell();
+        }
+      }
+      return number_of_vertices() - n;
+    }
+
+
     Vertex_handle insert(const Point &, Locate_type, const Face &, const Facet &, Full_cell_handle);
     Vertex_handle insert(const Point &, Full_cell_handle start = Full_cell_handle());
     Vertex_handle insert(const Point &, Vertex_handle);

Triangulation/include/CGAL/Triangulation_data_structure.h

@@ -33,6 +33,8 @@
 #include <stack>
 #include <set>
 
+#include <boost/iterator/counting_iterator.hpp>
+
 namespace CGAL {
 
 template<   class Dimen,