Improve its_split for large number of parts

sandboxes/its_neighbor_index/ItsNeighborIndex.cpp

@@ -325,7 +325,7 @@ FaceNeighborIndex its_create_neighbors_index_4(const indexed_triangle_set &its)
 // Create an index of faces belonging to each vertex. The returned vector can
 // be indexed with vertex indices and contains a list of face indices for each
 // vertex.
-static std::vector<std::vector<size_t>> create_vertex_faces_index(const indexed_triangle_set &its)
+std::vector<std::vector<size_t>> create_vertex_faces_index(const indexed_triangle_set &its)
 {
     std::vector<std::vector<size_t>> index;
 

sandboxes/its_neighbor_index/ItsNeighborIndex.hpp

@@ -12,4 +12,6 @@ std::vector<std::array<size_t, 3>> its_create_neighbors_index_6(const indexed_tr
 std::vector<std::array<size_t, 3>> its_create_neighbors_index_7(const indexed_triangle_set &its);
 FaceNeighborIndex its_create_neighbors_index_8(const indexed_triangle_set &its);
 std::vector<Vec3crd> its_create_neighbors_index_9(const indexed_triangle_set &its);
+
+std::vector<std::vector<size_t>> create_vertex_faces_index(const indexed_triangle_set &its);
 }

sandboxes/its_neighbor_index/main.cpp

@@ -2,6 +2,7 @@
 #include <fstream>
 #include <vector>
 #include <tuple>
+#include <random>
 
 #include "ItsNeighborIndex.hpp"
 
@@ -10,13 +11,15 @@
 
 namespace Slic3r {
 
+enum { IndexCreation, Split };
 struct MeasureResult
 {
-    double t_index_create = 0;
+    static constexpr const char * Names[] = {
-    double t_split = 0;
+        "Index creation [s]",
-    double memory = 0;
+        "Split [s]"
+    };
 
-    double full_time() const { return t_index_create + t_split; }
+    double measurements[std::size(Names)] = {0.};
 };
 
 template<class IndexCreatorFn>
@@ -30,27 +33,73 @@ static MeasureResult measure_index(const indexed_triangle_set &its, IndexCreator
         ItsNeighborsWrapper itsn{its, fn(its)};
         b.stop();
 
-        r.t_index_create += b.getElapsedSec();
+        r.measurements[IndexCreation] += b.getElapsedSec();
 
         b.start();
         auto res = its_split(itsn);
         b.stop();
 
-        if (res.size() != 2 || res[0].indices.size() != res[1].indices.size() )
+//        if (res.size() != 2 || res[0].indices.size() != res[1].indices.size() )
-            std::cerr << "Something is wrong, split result invalid" << std::endl;
+//            std::cerr << "Something is wrong, split result invalid" << std::endl;
 
-        r.t_split += b.getElapsedSec();
+        r.measurements[Split] += b.getElapsedSec();
     }
 
-    r.t_index_create /= 10;
+    r.measurements[IndexCreation] /= 10;
-    r.t_split /= 10;
+    r.measurements[Split] /= 10;
 
     return r;
 }
 
+const auto Seed = 0;// std::random_device{}();
+
+static indexed_triangle_set make_sphere_rnd(double radius, double detail)
+{
+    using namespace Slic3r;
+
+    auto sphere = its_make_sphere(radius, detail);
+
+    auto vfidx = create_vertex_faces_index(sphere);
+
+    const size_t vertexnum = sphere.vertices.size();
+    const size_t facenum   = sphere.indices.size();
+
+    std::mt19937 rng{Seed};
+    std::uniform_int_distribution<size_t> distv(sphere.vertices.size() / 2, sphere.vertices.size() - 1);
+    std::uniform_int_distribution<size_t> distf(sphere.indices.size() / 2, sphere.indices.size() - 1) ;
+
+    std::vector<bool> was(vertexnum / 2, false);
+
+    for (size_t i = 0; i < vertexnum / 2; ++i) {
+        size_t image = distv(rng);
+        if (was[image - vertexnum / 2]) continue;
+        was[image - vertexnum / 2] = true;
+
+        std::swap(sphere.vertices[i], sphere.vertices[image]);
+        for (size_t face_id : vfidx[i]) {
+            for (int &vi : sphere.indices[face_id])
+                if (vi == int(i)) vi = image;
+        }
+
+        for (size_t face_id : vfidx[image]) {
+            for (int &vi : sphere.indices[face_id])
+                if (vi == int(image)) vi = i;
+        }
+
+        std::swap(vfidx[i], vfidx[image]);
+    }
+
+    for (size_t i = 0; i < facenum / 2; ++i) {
+        size_t image = distf(rng);
+        std::swap(sphere.indices[i], sphere.indices[image]);
+    }
+
+    return sphere;
+}
+
 static indexed_triangle_set two_spheres(double detail)
 {
-    auto sphere1 = its_make_sphere(10., 2 * PI / detail), sphere2 = sphere1;
+    auto sphere1 = make_sphere_rnd(10., 2 * PI / detail), sphere2 = sphere1;
 
     its_transform(sphere1, identity3f().translate(Vec3f{-5.f, 0.f, 0.f}));
     its_transform(sphere2, identity3f().translate(Vec3f{5.f, 0.f, 0.f}));
@@ -60,22 +109,47 @@ static indexed_triangle_set two_spheres(double detail)
     return sphere1;
 }
 
+static indexed_triangle_set make_spheres(unsigned N, double detail)
+{
+    indexed_triangle_set ret, sphere = make_sphere_rnd(10., 2. * PI / detail);
+
+    for (unsigned i = 0u ; i < N; ++i)
+        its_merge(ret, sphere);
+
+    return ret;
+}
+
 constexpr double sq2 = std::sqrt(2.);
 
 static const std::pair<const std::string, indexed_triangle_set> ToMeasure[] = {
-    {"two_spheres_1x", two_spheres(60.)},
-    {"two_spheres_2x", two_spheres(120.)},
-    {"two_spheres_4x", two_spheres(240.)},
-    {"two_spheres_8x", two_spheres(480.)},
-    {"two_spheres_16x", two_spheres(2 * 480.)},
-    {"two_spheres_32x", two_spheres(2 * 2 * 480.)},
-
 //    {"two_spheres_1x", two_spheres(60.)},
-//    {"two_spheres_2x", two_spheres(sq2 * 60.)},
+//    {"two_spheres_2x", two_spheres(120.)},
-//    {"two_spheres_4x", two_spheres(2 * 60.)},
+//    {"two_spheres_4x", two_spheres(240.)},
-//    {"two_spheres_8x", two_spheres(sq2 * 2. * 60.)},
+//    {"two_spheres_8x", two_spheres(480.)},
-//    {"two_spheres_16x", two_spheres(4. * 60.)},
+//    {"two_spheres_16x", two_spheres(2 * 480.)},
-//    {"two_spheres_32x", two_spheres(sq2 * 4. * 60.)},
+//    {"two_spheres_32x", two_spheres(2 * 2 * 480.)},
+
+    {"two_spheres_1x", two_spheres(60.)},
+    {"two_spheres_2x", two_spheres(sq2 * 60.)},
+    {"two_spheres_4x", two_spheres(2 * 60.)},
+    {"two_spheres_8x", two_spheres(sq2 * 2. * 60.)},
+    {"two_spheres_16x", two_spheres(4. * 60.)},
+    {"two_spheres_32x", two_spheres(sq2 * 4. * 60.)},
+    {"two_spheres_64x", two_spheres(8. * 60.)},
+    {"two_spheres_128x", two_spheres(sq2 * 8. * 60.)},
+    {"two_spheres_256x", two_spheres(16. * 60.)},
+    {"two_spheres_512x", two_spheres(sq2 * 16. * 60.)}
+
+//    {"2_spheres", make_spheres(2, 60.)},
+//    {"4_spheres", make_spheres(4, 60.)},
+//    {"8_spheres", make_spheres(8, 60.)},
+//    {"16_spheres", make_spheres(16,  60.)},
+//    {"32_spheres", make_spheres(32, 60.)},
+//    {"64_spheres", make_spheres(64, 60.)},
+//    {"128_spheres", make_spheres(128, 60.)},
+//    {"256_spheres", make_spheres(256, 60.)},
+//    {"512_spheres", make_spheres(512, 60.)},
+//    {"1024_spheres", make_spheres(1024, 60.)}
 };
 
 static const auto IndexFunctions = std::make_tuple(
@@ -86,27 +160,32 @@ static const auto IndexFunctions = std::make_tuple(
     std::make_pair("vojta's vertex->face", [](const auto &its) { return measure_index(its, its_create_neighbors_index_9); }),
     std::make_pair("tamas's std::sort based", [](const auto &its) { return measure_index(its, its_create_neighbors_index_6); }),
     std::make_pair("tamas's tbb::parallel_sort based", [](const auto &its) { return measure_index(its, its_create_neighbors_index_7); }),
-    std::make_pair("tamas's map based", [](const auto &its) { return measure_index(its, its_create_neighbors_index_8); })/*,
+    std::make_pair("tamas's map based", [](const auto &its) { return measure_index(its, its_create_neighbors_index_8); }),
     std::make_pair("TriangleMesh split", [](const auto &its) {
 
-        MeasureResult ret;
+        MeasureResult r;
         for (int i = 0; i < 10; ++i) {
             TriangleMesh m{its};
             Benchmark b;
+
+            b.start();
+            m.repair(); // FIXME: this does more than just create neighborhood map
+            b.stop();
+            r.measurements[IndexCreation] += b.getElapsedSec();
+
             b.start();
-            m.repair();
             auto res = m.split();
             b.stop();
+            r.measurements[Split] += b.getElapsedSec();
 
-            if (res.size() != 2 || res[0]->size() != res[1]->size())
+//            if (res.size() != 2 || res[0]->size() != res[1]->size())
-                std::cerr << "Something is wrong, split result invalid" << std::endl;
+//                std::cerr << "Something is wrong, split result invalid" << std::endl;
-
-            ret.t_split += b.getElapsedSec();
         }
-        ret.t_split /= 10;
+        r.measurements[IndexCreation] /= 10;
+        r.measurements[Split] /= 10;
 
-        return ret;
+        return r;
-    })*/
+    })
 
 //    std::make_pair("Vojta's vertex->face index", [](const auto &its){
 //        Benchmark b;
@@ -147,6 +226,9 @@ int main(const int argc, const char * argv[])
         auto &m = ToMeasure[i];
         auto &name = m.first;
         auto &mesh = m.second;
+
+//        its_write_obj(mesh, (std::string(name) + ".obj").c_str());
+
         std::cout << "Mesh " << name << " has " << mesh.indices.size() << " faces and " << mesh.vertices.size() << " vertices." << std::endl;
         libnest2d::opt::metaloop::apply([&mesh, i, &results, &funcnames](int N, auto &e) {
             MeasureResult r = e.second(mesh);
@@ -165,6 +247,9 @@ int main(const int argc, const char * argv[])
 
     std::ostream &out = outfile.is_open() ? outfile : std::cout;
 
+    for (size_t m = 0; m < std::size(MeasureResult::Names); ++m) {
+        out << MeasureResult::Names[m] << "\n";
+        out << std::endl;
         out << "model;" ;
         for (const std::string &funcname : funcnames) {
             out << funcname << ";";
@@ -177,10 +262,11 @@ int main(const int argc, const char * argv[])
             const std::string &name = ToMeasure[i].first;
             out << name << ";";
             for (auto &r : result_row)
-            out << r.t_index_create << ";";
+                out << r.measurements[m] << ";";
 
             out << std::endl;
         }
+    }
 
     return 0;
 }

src/libslic3r/MeshSplitImpl.hpp

@@ -29,7 +29,7 @@ template<class NeighborIndex>
 std::vector<size_t> its_find_unvisited_neighbors(
     const indexed_triangle_set &its,
     const NeighborIndex &       neighbor_index,
-    std::vector<bool> &         visited)
+    std::vector<char> &         visited)
 {
     using stack_el = size_t;
 
@@ -92,26 +92,27 @@ void its_split(const Its &m, OutputIt out_it)
 
     const indexed_triangle_set &its = ItsWithNeighborsIndex_<Its>::get_its(m);
 
-    std::vector<bool> visited(its.indices.size(), false);
+    std::vector<char> visited(its.indices.size(), false);
 
-    const size_t        UNASSIGNED = its.vertices.size();
+    struct VertexConv {
-    std::vector<size_t> vidx_conv(its.vertices.size());
+        size_t part_id      = std::numeric_limits<size_t>::max();
+        size_t vertex_image;
+    };
+    std::vector<VertexConv> vidx_conv(its.vertices.size());
 
     const auto& neighbor_index = ItsWithNeighborsIndex_<Its>::get_index(m);
 
-    for (;;) {
+    for (size_t part_id = 0;; ++part_id) {
         std::vector<size_t> facets =
             its_find_unvisited_neighbors(its, neighbor_index, visited);
 
         if (facets.empty())
             break;
 
-        std::fill(vidx_conv.begin(), vidx_conv.end(), UNASSIGNED);
-
         // Create a new mesh for the part that was just split off.
         indexed_triangle_set mesh;
         mesh.indices.reserve(facets.size());
-        mesh.vertices.reserve(facets.size() * 3);
+        mesh.vertices.reserve(std::min(facets.size() * 3, its.vertices.size()));
 
         // Assign the facets to the new mesh.
         for (size_t face_id : facets) {
@@ -120,12 +121,12 @@ void its_split(const Its &m, OutputIt out_it)
             for (size_t v = 0; v < 3; ++v) {
                 auto vi = face(v);
 
-                if (vidx_conv[vi] == UNASSIGNED) {
+                if (vidx_conv[vi].part_id != part_id) {
-                    vidx_conv[vi] = mesh.vertices.size();
+                    vidx_conv[vi] = {part_id, mesh.vertices.size()};
                     mesh.vertices.emplace_back(its.vertices[size_t(vi)]);
                 }
 
-                new_face(v) = vidx_conv[vi];
+                new_face(v) = vidx_conv[vi].vertex_image;
             }
 
             mesh.indices.emplace_back(new_face);
@@ -150,7 +151,7 @@ template<class Its> bool its_is_splittable(const Its &m)
     const indexed_triangle_set &its = ItsWithNeighborsIndex_<Its>::get_its(m);
     const auto& neighbor_index = ItsWithNeighborsIndex_<Its>::get_index(m);
 
-    std::vector<bool> visited(its.indices.size(), false);
+    std::vector<char> visited(its.indices.size(), false);
     its_find_unvisited_neighbors(its, neighbor_index, visited);
 
     // Try finding an unvisited facet. If there are none, the mesh is not splittable.