Speed up island ordering on layers with thousands of islands

src/libslic3r/GCode/OrderingStrategies.cpp

@@ -121,9 +121,9 @@ bool tsp_remove_crossings(std::vector<size_t>& path, const Points& centers)
     // TSP helpers (2-opt, closing-edge rotation) that operate on the full cycle.
     size_t n_edges = pn;
 
-    // Scan for the first pair of intersecting edges worth reversing; returns {i, j} or {npos, npos} if none.
-    auto find_crossing = [&]() -> std::pair<size_t, size_t> {
-        for (size_t i = 0; i < n_edges; ++i) {
+    // Scan from edge `from` for the first pair of intersecting edges worth reversing; returns {i, j} or {npos, npos} if none.
+    auto find_crossing = [&](size_t from) -> std::pair<size_t, size_t> {
+        for (size_t i = from; i < n_edges; ++i) {
             const Point& ai = centers[path[i]];
             const Point& bi = centers[path[(i + 1) % pn]];
 
@@ -147,16 +147,27 @@ bool tsp_remove_crossings(std::vector<size_t>& path, const Points& centers)
         return {std::numeric_limits<size_t>::max(), std::numeric_limits<size_t>::max()};
     };
 
-    // Process crossings one at a time: find first, reverse it, restart scan.
-    // Every reversal shortens the cycle, so no tour repeats and the scan ends once no pair is worth
-    // reversing; the cap is a safety net.
-    int max_iters = static_cast<int>(pn * pn);
+    // Process crossings one at a time. Reversing path[ci+1..cj] replaces only edges ci and cj; the edges between
+    // them are the same segments traversed backwards. The scan therefore resumes at edge ci rather than at edge 0,
+    // so a pass over the tour costs one quadratic scan however many reversals it makes. Pairs before ci that
+    // involve a replaced edge are new too, so a scan that reaches the end starts over from edge 0, and the loop
+    // ends when a scan from edge 0 finds nothing.
+    // Every reversal shortens the cycle, so no tour repeats and the loop ends; the cap is a safety net.
+    size_t max_reversals = pn * pn;
     bool improved = false;
-    while (max_iters-- > 0) {
-        auto [ci, cj] = find_crossing();
-        if (ci == std::numeric_limits<size_t>::max()) break;
+    for (size_t from = 0;;) {
+        auto [ci, cj] = find_crossing(from);
+        if (ci == std::numeric_limits<size_t>::max()) {
+            if (from == 0)
+                break;
+            from = 0;
+            continue;
+        }
         improved = true;
         std::reverse(path.begin() + ci + 1, path.begin() + cj + 1);
+        if (--max_reversals == 0)
+            break;
+        from = ci;
     }
     return improved;
 }

tests/libslic3r/test_ordering_strategies.cpp

@@ -1,6 +1,7 @@
 #include <catch2/catch_test_macros.hpp>
 #include <catch2/generators/catch_generators.hpp>
 #include <catch2/catch_all.hpp>
+#include <catch2/generators/catch_generators.hpp>
 #include <catch2/matchers/catch_matchers.hpp>
 #include <catch2/matchers/catch_matchers_floating_point.hpp>
 #include <vector>
@@ -17,8 +18,10 @@
 #include "libslic3r/GCode/OrderingStrategies.hpp"
 #include "libslic3r/ShortestPath.hpp"
 #include "libslic3r/Geometry.hpp"
+#include "libslic3r/ShortestPath.hpp"
 
 #include <algorithm>
+#include <random>
 #include <unordered_set>
 #include <numeric>
 #include <random>
@@ -53,6 +56,27 @@ static bool has_crossings(const std::vector<size_t>& path, const Points& centers
     return false;
 }
 
+// Intersecting edge pairs whose reversal would shorten the cycle by more than SCALED_EPSILON, which
+// tsp_remove_crossings must leave none of.
+static size_t count_shortening_crossings(const std::vector<size_t>& path, const Points& centers)
+{
+    size_t pn = path.size();
+    size_t count = 0;
+    for (size_t i = 0; i < pn; ++i) {
+        for (size_t j = i + 2; j < pn; ++j) {
+            if (i == 0 && j == pn - 1) continue;
+            const Point& ai = centers[path[i]];
+            const Point& bi = centers[path[(i + 1) % pn]];
+            const Point& aj = centers[path[j]];
+            const Point& bj = centers[path[(j + 1) % pn]];
+            if (Geometry::segments_intersect(ai, bi, aj, bj) &&
+                ai.distance_to(bi) + aj.distance_to(bj) - ai.distance_to(aj) - bi.distance_to(bj) > SCALED_EPSILON)
+                ++count;
+        }
+    }
+    return count;
+}
+
 static bool is_permutation(const std::vector<size_t>& path, size_t n)
 {
     if (path.size() != n) return false;
@@ -185,6 +209,29 @@ TEST_CASE("tsp_remove_crossings reverses overlapping collinear edges when that s
     CHECK_THAT(tsp_cycle_path_length(path, centers), Catch::Matchers::WithinAbs(6000000., 1.));
 }
 
+TEST_CASE("tsp_remove_crossings leaves no shortening crossing in a tour of thousands of grid-snapped points", "[TSPPostProcessing]") {
+    // 4000 island centroids on a 170 mm plate, snapped to 1 mm as GCode::process_layer does, so many of them share
+    // a position or a line. The nearest-neighbour tour the island ordering starts from has crossings to remove.
+    std::mt19937 rng(GENERATE(1u, 2u, 3u));
+    Points centers;
+    for (int i = 0; i < 4000; ++i) {
+        const double x = double(rng() % 171);
+        const double y = double(rng() % 171);
+        centers.push_back(Point::new_scale(x, y));
+    }
+    std::vector<size_t> path = chain_points(centers);
+    const size_t first = path.front();
+    const double length_before = tsp_cycle_path_length(path, centers);
+    REQUIRE(count_shortening_crossings(path, centers) > 0);
+
+    CHECK(tsp_remove_crossings(path, centers));
+    REQUIRE(is_permutation(path, centers.size()));
+    CHECK(count_shortening_crossings(path, centers) == 0);
+    CHECK(tsp_cycle_path_length(path, centers) < length_before);
+    // Callers that start the tour near a given point rely on the first node staying in place.
+    CHECK(path.front() == first);
+}
+
 
 
 TEST_CASE("tsp_rotate_minimize_closing shortens closing edge", "[TSPPostProcessing]") {