diff --git a/docs/HLSD/travel-optimization.md b/docs/HLSD/travel-optimization.md new file mode 100644 index 0000000000..119133cee2 --- /dev/null +++ b/docs/HLSD/travel-optimization.md @@ -0,0 +1,198 @@ +# Travel optimization — High Level Design + +## Purpose and scope + +Inside an island, the infill paths of each region are chained into one print +order by a greedy pass that joins the closest free path ends first. That pass +leaves its longest joins for last, and they often cross the whole island. Two +things shorten the travel between infill paths: + +- the fill generators join their lines into as few paths as they can, since + every separate path needs a travel to it; +- the G-code export refines the greedy chain with a local search that shortens + the straight-line travel between the paths. + +The refinement covers every chain that `chain_and_reorder_extrusion_entities()` +builds from a start point: the infill and ironing of `GCode::extrude_infill()` +and the support of `GCode::extrude_support()`, as well as the order of objects +and islands under `print_order` Default. It changes only the order and the +direction of the paths, not the paths or the walls. Neither part has a setting. + +## The chain + +`ObjectByExtruder::Island::Region::append()` flattens every sortable fill +collection into its paths, so the paths of all the fills of a region in one +island are chained together. A collection that must keep its order (`no_sort`: +monotonic, concentric, a forced fill order, the flow calibration) stays a single +entry that cannot be reversed. Paths that must keep their direction cannot be +reversed either: Grid's, since it alternates their direction between layers for +strength, and those of a fill with a forced order. A loop starts and ends at the +same point, so reversing it changes nothing. + +`extrude_infill()` chains each region from the current nozzle position with +`chain_and_reorder_extrusion_entities()`, and chains the paths inside a nested +sortable collection the same way; `extrude_support()` chains the support paths +of a layer so too, unless they keep their order. The overload without a start +point, used to chain collections while slicing, keeps the greedy chain. The +order the fill generators give their paths does not reach the G-code for a +sortable fill, which is chained again. + +## Connecting the lines + +Most patterns cut their lines at the infill boundary, and `Fill::connect_infill()` +joins pairs of line ends that are neighbors on a boundary contour by the arc of +wall between them, shortest arcs first, as long as neither end is joined yet and +the join closes no loop. An arc is joined when it is shorter than +`infill_anchor_max`, or than two anchors and a line width if that is longer: +two ends further apart each get an anchor toward the other, and where those +would meet, the join prints about as much wall as one path instead of two. + +Regular patterns leave many arcs of the same length, and the order among them +decides the pairing. In an arbitrary order the pairs along one wall can fall out +of step, such as lines 1–2 and 3–4, then 6–7 and 8–9, which leaves the end of +line 5 loose, and every loose end starts or ends a separate path. So the arcs +are also sorted in steps of half the line spacing and, within a step, along each +contour forwards and backwards, which pairs the ends of a wall in one phase. +Each order is replayed on the flags and a union-find of the lines, without +building any path, and the one that joins the most lines is used; a tie keeps +the shortest-first order. + +The replay includes the last pass, which joins any two open ends that are +neighbors on a contour. With multiline infill the main pass leaves the copies of +one line to that pass, so that each copy reaches a wall of its own. Two +neighboring groups cut at both walls can then be joined at either wall by arcs +of the same length: at one, their open ends stay next to their own copies, which +the last pass joins, so each line prints out and back; at the other, the copies +end up in separate paths printed apart. + +A line end left unpaired gets an anchor, a piece of the wall at most +`infill_anchor` long. It runs toward the nearest open end of another path, +looking up to 64 ends along the contour, since the next path usually starts +there; without one it takes the side with more free wall. With multiline infill +an anchor toward a copy of the same line stops short of the copy rather than +printing over it. A copy is a neighboring end less than the group width +(`multiline` spacings) across the line, at any wall angle, while the next group +lies a whole line spacing away. Meeting the wall at an angle θ, the copy covers +half a width over sin θ of the wall on each side of its end, so the anchor stops +that far and half its own width before it, or is left out; the end still counts +as anchored, so the last pass pairs the other ends as before. + +Last, `Fill::fill_surface_extrusion()` joins the sparse infill lines whose ends +are at most two line spacings apart with a straight piece. With the anchors off +(`infill_anchor_max` 0) nothing above runs and this is the only joining; with +them on, it takes what the steps above leave next to each other, such as two +anchors that stop just short of each other or a short piece cut off by a corner. +A gap that short leaves no room to print anything in between, and as a travel it +is the one hop a reordering can trade for a longer travel that retracts. +Multiline copies leave the most such ends: one spacing apart across the line, +their ends where it meets a wall at an angle θ are one spacing over sin θ apart +along it, so two spacings reach walls met at 30° and steeper. + +Ends that meet are joined first, then the gaps between the leading ends of their +lines, along each line's overall direction, then the rest, shortest first within +each, without closing a loop or crossing a line or another join. Taking the +leading ends first joins every pair of copies at the same side, which leaves +their open ends next to each other; by length alone, neighboring pairs were +joined at different walls and the travel between them crossed the island. Lines +that cannot reverse, such as Grid's, only join the end of one to the start of the +next, and fills that keep their order are left alone. + +## The search + +A chain of up to eight paths is found exactly, by dynamic programming over the +subsets of paths: a state is the paths visited so far, the last one and its +direction, about four thousand states for eight paths. The shared greedy +chaining, `chain_segments_greedy_constrained_reversals_()`, does the same for two +to four segments, which also orders up to four islands of a layer and the loops +of a perimeter collection exactly. Either way the first path is no longer +always the one nearest the start; with a free end, a path near the start can be +worth leaving for last. + +Objects (`print_order` Default) and the tour of islands across objects are +visited in the same order on every layer, so `chain_points_with_postprocessing()` +orders them as a closed tour rather than an open chain. Up to ten points the +shortest tour is found exactly by the same dynamic programming with the first +point fixed; with a start point the tour then begins at the point nearest it. +Larger sets take the greedy chain and alternate 2-opt, or-opt and crossing +removal; or-opt moves runs of up to three points next to one of the eight points +nearest an end of the run, found on the same grid as the chain's candidates. + +For longer chains, the search refines the greedy chain. On an island cut into many short pieces, +such as infill without anchors, that chain can sweep away from the nozzle and +come back for the rest with one long jump, which no local move below removes. +So a nearest-next chain, going from the nozzle to the nearest free path end each +time, is built too; when it is already shorter than the refined greedy chain, it +is refined and used instead. Elsewhere it costs one grid query per path. + +`improve_chain_by_local_search()` sees the chain as a sequence of segments with +an entry and an exit point each. Its cost is the straight distance from the +nozzle to the first entry plus the distance from every exit to the next entry; +the end of the chain is free. It applies three kinds of moves: + +- 2-opt reverses a run of segments, which flips each of them and replaces the + two links at the ends of the run. +- Or-opt moves a run of up to three segments between two others, flipped or + not. +- A run of any length that starts or ends next to the nozzle moves to the front. + Local search can settle on a chain that starts with a block of paths far from + the nozzle while the paths next to it sit in the middle, and moving that block + takes more than three segments at once. + +A move that would flip a segment that cannot be reversed is skipped, so such +segments are only reordered. + +A move is tried only if every link it creates joins an end point to one of the +32 end points nearest to it, which covers every move on a chain of up to 17 +segments, as most islands are; with 8 candidates the search missed most of the +gain on grid-like patterns such as Triangles. The candidates come from a uniform +grid of about two end points per cell, searched in rings of cells around the +query. `KDTreeIndirect` visits its subtrees in a fixed order rather than the +nearer one first, so a query for 32 points walks most of the tree, and on long +chains finding the candidates would cost more than the search itself. + +The candidates of an end point are tried nearest first, and only while they are +closer to it than the link the move would replace there. A 2-opt move that pays +off creates a link shorter than the one it replaces at one of its two ends, so +this misses no 2-opt move. An or-opt move is tried only while its new link is +shorter than what taking the run out saves, which skips a few or-opt moves that +would still pay off but cuts the search on long chains to a fraction. + +The search works through a queue of segments and applies the first move that +shortens the chain around the segment it examines. After a move, only the +segments next to a changed link are queued again. A new link can also make room +for a run elsewhere that no segment next to it would try, so when the queue +empties after any move, every segment is queued once more, and the search ends +when such a sweep finds no move. + +Every move shortens the chain by more than `SCALED_EPSILON`, so the result is +never longer than the greedy or the nearest-next chain and the search ends; a +cap of 50 moves per segment only bounds its time. The search has no randomness +and visits the segments in a fixed order, so the same input always gives the +same G-code. + +## Cost + +Finding the candidates of an end point visits a few cells, and the nearest-next +chain makes one such query per path, which only grows where the paths around the +nozzle are used up and the chain jumps. Examining a segment tries at most a few +hundred candidate links, usually only the nearest few, plus a scan of the chain +for the move to the front when it is the first segment, and a move rewrites only +the range of the chain it changes, so the search stays close to linear in the +number of segments. Most islands hold a few to a few dozen paths; the cost shows +on infill cut into thousands of pieces per island, such as dense infill with the +anchors off. The straight joins only look at the line segments near a candidate +join for crossings. + +## Limitations + +The search minimizes straight-line distance. It does not know which travels +leave the infill area and retract under `reduce_infill_retraction`, nor the +detours of `reduce_crossing_wall`, and a travel's time is not proportional to +its length: a short travel is bound by acceleration, so its time grows about +with the square root of its length. A chain that is shorter in total but trades +very short travels for longer ones, or for ones past `retraction_minimum_travel`, +can occasionally take a little longer to print. + +The end of the chain is free: it is not steered toward the next island or the +start of the next layer. The regions of an island are chained one after the +other, in region order, each from where the previous one ended. diff --git a/src/libslic3r/Fill/FillBase.cpp b/src/libslic3r/Fill/FillBase.cpp index a852f2c42b..63ebff49c2 100644 --- a/src/libslic3r/Fill/FillBase.cpp +++ b/src/libslic3r/Fill/FillBase.cpp @@ -11,6 +11,10 @@ #include #include #include +#include +#include +#include +#include "../AABBTreeLines.hpp" #include "../ClipperUtils.hpp" #include "../EdgeGrid.hpp" #include "../Geometry.hpp" @@ -206,6 +210,9 @@ void Fill::fill_surface_extrusion(const Surface* surface, const FillParams& para thick_polylines.clear(); } else { + // Orca: sparse infill lines whose ends are next to each other print as one path, without a travel between them. + if (params.extrusion_role == erInternalInfill && ! eec->no_sort) + join_close_ends(polylines, scale_(2. * this->spacing), params.can_reverse); extrusion_entities_append_paths( eec->entities, std::move(polylines), params.extrusion_role, @@ -1644,6 +1651,8 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vector SCALED_EPSILON ? std::max(anchor_length_max, 2. * (anchor_length + line_half_width)) : anchor_length_max; #if 0 // Connection from end of one infill line to the start of another infill line. @@ -1720,6 +1729,10 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vector open_end(graph.map_infill_end_point_to_boundary.size(), 1); + auto end_index = [&graph](const ContourIntersectionPoint *cp) { return size_t(cp - graph.map_infill_end_point_to_boundary.data()); }; + struct Arc { ContourIntersectionPoint *intersection; double arc_length; @@ -1731,6 +1744,81 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vector &order) { + std::vector consumed(graph.map_infill_end_point_to_boundary.size()); + for (size_t i = 0; i < consumed.size(); ++ i) + consumed[i] = graph.map_infill_end_point_to_boundary[i].consumed; + std::vector parent(infill_ordered.size()); + std::iota(parent.begin(), parent.end(), 0); + auto find = [&parent](size_t i) { + while (parent[i] != i) + i = parent[i] = parent[parent[i]]; + return i; + }; + size_t joins = 0; + for (const Arc &arc : order) { + const size_t end1 = end_index(arc.intersection), end2 = end_index(arc.intersection->next_on_contour); + if (consumed[end1] || consumed[end2] || (params.multiline > 1 && arc.arc_length < scale_(spacing) * params.multiline)) + continue; + const size_t polyline1 = find(end1 / 2), polyline2 = find(end2 / 2); + if (polyline1 == polyline2) + continue; + if (arc.arc_length < join_length) { + parent[polyline2] = polyline1; + ++ joins; + } else if (anchor_length <= SCALED_EPSILON) + continue; + consumed[end1] = consumed[end2] = 1; + } + // The last pass joins the open ends next to each other, also those of a multiline group. + for (size_t end = 0; end < consumed.size(); ++ end) { + const ContourIntersectionPoint &cp = graph.map_infill_end_point_to_boundary[end]; + if (consumed[end] || cp.contour_idx == boundary_idx_unconnected) + continue; + const double length = graph.boundary_params[cp.contour_idx].back(); + auto open = [&](const ContourIntersectionPoint *other, bool trimmed) { + return other != &cp && ! consumed[end_index(other)] && ! trimmed; + }; + const double lprev = open(cp.prev_on_contour, cp.prev_trimmed || cp.prev_on_contour->next_trimmed) ? + path_length_along_contour_ccw(cp.prev_on_contour, &cp, length) : std::numeric_limits::max(); + const double lnext = open(cp.next_on_contour, cp.next_trimmed || cp.next_on_contour->prev_trimmed) ? + path_length_along_contour_ccw(&cp, cp.next_on_contour, length) : std::numeric_limits::max(); + consumed[end] = anchor_length > SCALED_EPSILON; + for (const double l : { std::min(lprev, lnext), std::max(lprev, lnext) }) { + if (l > join_length) + break; + const size_t other = end_index(l == lprev ? cp.prev_on_contour : cp.next_on_contour); + const size_t polyline1 = find(end / 2), polyline2 = find(other / 2); + if (polyline1 == polyline2) + continue; + parent[polyline2] = polyline1; + consumed[end] = consumed[other] = 1; + ++ joins; + break; + } + } + return joins; + }; + // Equal arcs taken along the contour pair a wall in one phase; keep the order that joins the most lines. + size_t best_joins = count_joins(arches); + const double bucket = 0.5 * scale_(spacing); + for (const bool backwards : { false, true }) { + std::vector along_contour = arches; + std::sort(along_contour.begin(), along_contour.end(), [bucket, backwards](const Arc &l, const Arc &r) { + const long long bucket_l = std::llround(l.arc_length / bucket), bucket_r = std::llround(r.arc_length / bucket); + if (bucket_l != bucket_r) + return bucket_l < bucket_r; + if (l.intersection->contour_idx != r.intersection->contour_idx) + return l.intersection->contour_idx < r.intersection->contour_idx; + return backwards ? l.intersection->param > r.intersection->param : l.intersection->param < r.intersection->param; + }); + if (const size_t joins = count_joins(along_contour); joins > best_joins) { + best_joins = joins; + arches = std::move(along_contour); + } + } } //FIXME improve the Traveling Salesman problem with 2-opt and 3-opt local optimization. @@ -1752,7 +1840,7 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vectorpoint_idx] == polyline1.points.front() || contour[cp1->point_idx] == polyline1.points.back()); if (contour[cp1->point_idx] == polyline1.points.front()) @@ -1761,6 +1849,7 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vectorpoint_idx] == polyline2.points.back()) polyline2.reverse(); take(polyline1, polyline2, contour, cp1, cp2, false); + open_end[end_index(cp1)] = open_end[end_index(cp2)] = 0; // Mark the second polygon as merged with the first one. if (polyline_idx2 < polyline_idx1) { polyline2 = std::move(polyline1); @@ -1778,6 +1867,23 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vectorprev_on_contour : cp->next_on_contour; + if (other == nullptr || other == &from) + break; + distance += backwards ? closed_contour_distance_cw(cp->param, other->param, length) : closed_contour_distance_ccw(cp->param, other->param, length); + if (open_end[end_index(other)] && get_and_update_merged_with(end_index(other) / 2) != polyline_idx) + return distance; + cp = other; + } + return std::numeric_limits::max(); + }; + // Connect the remaining open infill lines to the perimeter lines if possible. for (ContourIntersectionPoint &contour_point : graph.map_infill_end_point_to_boundary) if (! contour_point.consumed && contour_point.contour_idx != boundary_idx_unconnected) { @@ -1796,7 +1902,7 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vector::max() || l > anchor_length_max) + if (l == std::numeric_limits::max() || l > join_length) break; // Take the complete contour. bool reversed = l == lprev; @@ -1815,6 +1921,7 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vectorpoint_idx] == polyline2.points.back()) polyline2.reverse(); take(polyline, polyline2, contour, &contour_point, cp2, reversed); + open_end[end_index(&contour_point)] = open_end[end_index(cp2)] = 0; if (polyline_idx < polyline_idx2) { // Mark the second polyline as merged with the first one. merged_with[polyline_idx2] = polyline_idx; @@ -1829,18 +1936,35 @@ void Fill::connect_infill(Polylines &&infill_ordered, const std::vector SCALED_EPSILON) { - // Which to take? One could optimize for: - // 1) Shortest path - // 2) Hook length - // ... - // Let's take the longer now, as this improves the chance of another hook to be placed on the other side of this contour point. - double l = std::max(contour_point.contour_not_taken_length_prev, contour_point.contour_not_taken_length_next); - if (l > SCALED_EPSILON) { - if (contour_point.contour_not_taken_length_prev > contour_point.contour_not_taken_length_next) - take_limited(polyline, contour, contour_params, &contour_point, contour_point.prev_on_contour, true, anchor_length, line_half_width); - else - take_limited(polyline, contour, contour_params, &contour_point, contour_point.next_on_contour, false, anchor_length, line_half_width); + // Hook toward the nearest open end of another polyline, where the next path likely starts. + const double free_prev = contour_point.contour_not_taken_length_prev; + const double free_next = contour_point.contour_not_taken_length_next; + const double to_prev = distance_to_open_end(contour_point, polyline_idx, true); + const double to_next = distance_to_open_end(contour_point, polyline_idx, false); + bool backwards = free_prev > free_next; + if (to_prev < to_next && free_prev > SCALED_EPSILON) + backwards = true; + else if (to_next < to_prev && free_next > SCALED_EPSILON) + backwards = false; + // Orca: with multiline, a hook toward a copy stops short of it; the copy takes half a width over the sine of its wall angle. + double length = anchor_length; + if (const ContourIntersectionPoint *other = backwards ? contour_point.prev_on_contour : contour_point.next_on_contour; + params.multiline > 1 && other != &contour_point) { + const Point &end = contour[contour_point.point_idx]; + Vec2d dir = Vec2d::Zero(); + for (size_t i = 1; i < polyline.size() && dir.squaredNorm() < sqr(2. * scale_(spacing) * params.multiline); ++ i) + dir = (end - polyline.points[polyline.points.front() == end ? i : polyline.size() - 1 - i]).cast(); + const Vec2d to_copy = (contour[other->point_idx] - end).cast(); + const double across = dir.squaredNorm() > 0. ? std::abs(cross2(dir.normalized(), to_copy)) : 0.; + if (across > SCALED_EPSILON && across < scale_(spacing) * params.multiline) { + const double gap = backwards ? path_length_along_contour_ccw(other, &contour_point, contour_params.back()) : + path_length_along_contour_ccw(&contour_point, other, contour_params.back()); + length = std::min(length, gap - line_half_width * (1. + to_copy.norm() / across)); + } } + if ((backwards ? free_prev : free_next) > SCALED_EPSILON) + take_limited(polyline, contour, contour_params, &contour_point, backwards ? contour_point.prev_on_contour : contour_point.next_on_contour, + backwards, std::max(0., length), line_half_width); } } @@ -2798,4 +2922,123 @@ void multiline_fill(Polylines& polylines, const FillParams& params, float spacin polylines = std::move(all_polylines); } +void join_close_ends(Polylines &polylines, double max_gap, bool can_reverse) +{ + const size_t n = polylines.size(); + if (n < 2 || max_gap <= 0.) + return; + static constexpr size_t npos = size_t(-1); + // End 2 * i is the first point of polyline i, 2 * i + 1 its last. + auto end_point = [&polylines](size_t e) -> const Point & { return (e & 1) ? polylines[e >> 1].points.back() : polylines[e >> 1].points.front(); }; + // Ends in cells of max_gap, so the ends within max_gap of one are in its cell or the eight around it. + auto cell_of = [max_gap](const Point &p) { return Vec2i64(int64_t(std::floor(p.x() / max_gap)), int64_t(std::floor(p.y() / max_gap))); }; + auto key = [](const Vec2i64 &c) { return (uint64_t(c.x()) << 32) ^ uint64_t(uint32_t(c.y())); }; + std::unordered_map> grid; + for (size_t e = 0; e < 2 * n; ++ e) + if (polylines[e >> 1].size() >= 2) + grid[key(cell_of(end_point(e)))].push_back(e); + + // The end of a polyline that comes last along its overall direction, taken along the dominant axis. + std::vector backwards(n, 0); + for (size_t i = 0; i < n; ++ i) + if (polylines[i].size() >= 2) { + const Vec2d dir = (polylines[i].points.back() - polylines[i].points.front()).cast(); + backwards[i] = std::abs(dir.x()) >= std::abs(dir.y()) ? dir.x() < 0. : dir.y() < 0.; + } + auto trailing = [&backwards](size_t e) { return bool(e & 1) != bool(backwards[e >> 1]); }; + // Ends that meet first, then shortest gaps between leading ends, then the others: the copies of neighboring multiline + // groups are all joined at the same side, which leaves their open ends next to each other rather than on alternating walls. + std::vector> gaps; + for (const auto &[cell_key, ends] : grid) { + const Vec2i64 c = cell_of(end_point(ends.front())); + for (int64_t dx = -1; dx <= 1; ++ dx) + for (int64_t dy = -1; dy <= 1; ++ dy) + if (auto it = grid.find(key(c + Vec2i64(dx, dy))); it != grid.end()) + for (size_t a : ends) + for (size_t b : it->second) + if (a < b && (a >> 1) != (b >> 1) && (can_reverse || (a & 1) != (b & 1))) + if (const double d = (end_point(a) - end_point(b)).cast().norm(); d <= max_gap) + gaps.emplace_back(d > SCALED_EPSILON, d > SCALED_EPSILON && (trailing(a) || trailing(b)), d, a, b); + } + if (gaps.empty()) + return; + std::sort(gaps.begin(), gaps.end()); + + std::vector partner(2 * n, npos), parent(n); + std::iota(parent.begin(), parent.end(), 0); + auto find = [&parent](size_t i) { + while (parent[i] != i) + i = parent[i] = parent[parent[i]]; + return i; + }; + // A join may not cross a line or another join; touching them at its own ends is fine, as where a line ends on another. + // Only the segments through a cell next to a candidate join can cross one. + std::unordered_set near_gaps; + for (const auto &[apart, trail, d, a, b] : gaps) + for (size_t e : { a, b }) + for (int64_t dx = -1; dx <= 1; ++ dx) + for (int64_t dy = -1; dy <= 1; ++ dy) + near_gaps.insert(key(cell_of(end_point(e)) + Vec2i64(dx, dy))); + Lines near_lines; + for (const Polyline &pl : polylines) + for (size_t i = 1; i < pl.size(); ++ i) { + const Vec2i64 c0 = cell_of(pl.points[i - 1]), c1 = cell_of(pl.points[i]); + const Vec2i64 lo = c0.cwiseMin(c1), hi = c0.cwiseMax(c1); + bool near = (hi - lo).maxCoeff() > 3; + for (int64_t x = lo.x(); x <= hi.x() && ! near; ++ x) + for (int64_t y = lo.y(); y <= hi.y() && ! near; ++ y) + near = near_gaps.count(key(Vec2i64(x, y))) > 0; + if (near) + near_lines.emplace_back(pl.points[i - 1], pl.points[i]); + } + const AABBTreeLines::LinesDistancer lines(std::move(near_lines)); + auto crosses = [&](size_t a, size_t b) { + const Vec2d p = end_point(a).cast(), q = end_point(b).cast(); + const double length = (q - p).norm(); + if (length <= 2. * SCALED_EPSILON) + return false; + const Vec2d dir = (q - p) / length; + const Line inner(Point((p + SCALED_EPSILON * dir).cast()), Point((q - SCALED_EPSILON * dir).cast())); + if (! lines.get_lines().empty() && ! lines.intersections_with_line(inner).empty()) + return true; + const Vec2i64 c = cell_of(end_point(a)); + for (int64_t dx = -3; dx <= 3; ++ dx) + for (int64_t dy = -3; dy <= 3; ++ dy) + if (auto it = grid.find(key(c + Vec2i64(dx, dy))); it != grid.end()) + for (size_t e : it->second) + if (partner[e] != npos && Geometry::segments_intersect(inner.a, inner.b, end_point(e), end_point(partner[e]))) + return true; + return false; + }; + for (const auto &[apart, trail, d, a, b] : gaps) + if (partner[a] == npos && partner[b] == npos && find(a >> 1) != find(b >> 1) && ! crosses(a, b)) { + partner[a] = b; + partner[b] = a; + parent[find(a >> 1)] = find(b >> 1); + } + + // Walk each group from a polyline with a free end, entering the next polyline at the end joined to the last one. + Polylines out; + std::vector done(n, 0); + auto walk = [&](size_t first_end) { + Polyline joined; + for (size_t e = first_end; e != npos;) { + Polyline &pl = polylines[e >> 1]; + done[e >> 1] = 1; + if (e & 1) + pl.reverse(); + joined.points.insert(joined.points.end(), pl.points.begin() + (! joined.empty() && joined.points.back() == pl.points.front()), pl.points.end()); + e = partner[e ^ 1]; + } + out.emplace_back(std::move(joined)); + }; + for (size_t i = 0; i < n; ++ i) + if (! done[i] && partner[2 * i] == npos) + walk(2 * i); + for (size_t i = 0; i < n; ++ i) + if (! done[i] && partner[2 * i + 1] == npos) + walk(2 * i + 1); + polylines = std::move(out); +} + } // namespace Slic3r diff --git a/src/libslic3r/Fill/FillBase.hpp b/src/libslic3r/Fill/FillBase.hpp index d1f07af77c..f60fde6d1d 100644 --- a/src/libslic3r/Fill/FillBase.hpp +++ b/src/libslic3r/Fill/FillBase.hpp @@ -258,8 +258,11 @@ public: static coord_t _adjust_solid_spacing(const coord_t width, const coord_t distance); }; - //Fill Multiline + //Fill Multiline void multiline_fill(Polylines& polylines, const FillParams& params, float spacing); + // Joins polylines whose ends are at most max_gap apart, shortest gaps first and without closing a loop, so each such + // group prints as one path. Without can_reverse only the end of one polyline joins the start of another. + void join_close_ends(Polylines &polylines, double max_gap, bool can_reverse); } // namespace Slic3r #endif // slic3r_FillBase_hpp_ diff --git a/src/libslic3r/GCode/OrderingStrategies.cpp b/src/libslic3r/GCode/OrderingStrategies.cpp index 43bf9ecf87..98b8dfa37c 100644 --- a/src/libslic3r/GCode/OrderingStrategies.cpp +++ b/src/libslic3r/GCode/OrderingStrategies.cpp @@ -8,6 +8,7 @@ #include "libslic3r/PrintConfig.hpp" #include +#include #include #include #include @@ -153,6 +154,41 @@ bool tsp_remove_crossings(std::vector& path, const Points& centers) return improved; } +std::vector tsp_exact_cycle(const Points& centers) +{ + const size_t n = centers.size(); + std::vector path(n); + std::iota(path.begin(), path.end(), 0); + if (n <= 3) return path; + // State mask * m + j: the points 1..m of mask visited after point 0, the last one j + 1. + const size_t m = n - 1; + assert(m <= 16); + std::vector cost((size_t(1) << m) * m, std::numeric_limits::max()); + std::vector from(cost.size(), size_t(-1)); + auto dist = [¢ers](size_t a, size_t b) { return (centers[a].cast() - centers[b].cast()).norm(); }; + for (size_t j = 0; j < m; ++j) + cost[(size_t(1) << j) * m + j] = dist(0, j + 1); + for (size_t mask = 1; mask < (size_t(1) << m); ++mask) + for (size_t j = 0; j < m; ++j) + if (const double c = cost[mask * m + j]; c < std::numeric_limits::max()) + for (size_t k = 0; k < m; ++k) + if (!(mask >> k & 1)) { + const size_t next = (mask | size_t(1) << k) * m + k; + if (const double nc = c + dist(j + 1, k + 1); nc < cost[next]) { + cost[next] = nc; + from[next] = mask * m + j; + } + } + const size_t full = (size_t(1) << m) - 1; + size_t best = full * m; + for (size_t j = 1; j < m; ++j) + if (cost[full * m + j] + dist(j + 1, 0) < cost[best] + dist(best % m + 1, 0)) + best = full * m + j; + for (size_t s = best, k = n; k-- > 1; s = from[s]) + path[k] = s % m + 1; + return path; +} + void tsp_rotate_minimize_closing(std::vector& path, const Points& centers) { size_t pn = path.size(); diff --git a/src/libslic3r/GCode/OrderingStrategies.hpp b/src/libslic3r/GCode/OrderingStrategies.hpp index a42a338a52..95a6e77bcf 100644 --- a/src/libslic3r/GCode/OrderingStrategies.hpp +++ b/src/libslic3r/GCode/OrderingStrategies.hpp @@ -29,6 +29,9 @@ bool tsp_2opt_improve(std::vector& path, const Points& centers, int max_ // Returns true if any crossing was removed. bool tsp_remove_crossings(std::vector& path, const Points& centers); +// The shortest cycle through the points, starting at point 0, by dynamic programming over subsets; for a few points only. +std::vector tsp_exact_cycle(const Points& centers); + // Rotate the cycle so the closing edge (last -> first) is minimized. void tsp_rotate_minimize_closing(std::vector& path, const Points& centers); diff --git a/src/libslic3r/ShortestPath.cpp b/src/libslic3r/ShortestPath.cpp index bbea1a6900..31e02f35bc 100644 --- a/src/libslic3r/ShortestPath.cpp +++ b/src/libslic3r/ShortestPath.cpp @@ -26,8 +26,10 @@ #include "Print.hpp" #include "GCode/OrderingStrategies.hpp" +#include #include #include +#include #include "MultiMaterialSegmentation.hpp" namespace Slic3r { @@ -96,6 +98,46 @@ std::vector> chain_segments_closest_point(std::vector> chain_segments_exact(const std::vector &end_points, const std::vector &could_reverse, const Vec2d *start) +{ + const size_t n = could_reverse.size(); + assert(n >= 1 && n <= 12 && end_points.size() == 2 * n); + // State (mask * n + j) * 2 + o: the segments of mask visited, the last one j, flipped when o. + const size_t num_states = (size_t(1) << n) * n * 2; + std::vector cost(num_states, std::numeric_limits::max()); + std::vector from(num_states, size_t(-1)); + auto entry = [&end_points](size_t j, size_t o) -> const Vec2d & { return end_points[2 * j + o]; }; + auto exit = [&end_points](size_t j, size_t o) -> const Vec2d & { return end_points[2 * j + 1 - o]; }; + for (size_t j = 0; j < n; ++ j) + for (size_t o = 0; o <= size_t(could_reverse[j] != 0); ++ o) + cost[((size_t(1) << j) * n + j) * 2 + o] = start == nullptr ? 0. : (entry(j, o) - *start).norm(); + for (size_t mask = 1; mask < (size_t(1) << n); ++ mask) + for (size_t s = mask * n * 2; s < (mask + 1) * n * 2; ++ s) + if (cost[s] < std::numeric_limits::max()) + for (size_t k = 0; k < n; ++ k) + if (! (mask >> k & 1)) + for (size_t o = 0; o <= size_t(could_reverse[k] != 0); ++ o) { + const size_t next = (((mask | size_t(1) << k) * n) + k) * 2 + o; + const double c = cost[s] + (entry(k, o) - exit(s / 2 % n, s & 1)).norm(); + if (c < cost[next]) { + cost[next] = c; + from[next] = s; + } + } + const size_t full = (size_t(1) << n) - 1; + size_t best = full * n * 2; + for (size_t s = best; s < num_states; ++ s) + if (cost[s] < cost[best]) + best = s; + std::vector> out(n); + for (size_t s = best, k = n; k -- > 0; s = from[s]) + out[k] = { s / 2 % n, (s & 1) != 0 }; + return out; +} + // Chain perimeters (always closed) and thin fills (closed or open) using a greedy algorithm. // Solving a Traveling Salesman Problem (TSP) with the modification, that the sites are not always points, but points and segments. // Solving using a greedy algorithm, where a shortest edge is added to the solution if it does not produce a bifurcation or a cycle. @@ -104,6 +146,7 @@ std::vector> chain_segments_closest_point(std::vector std::vector> chain_segments_greedy_constrained_reversals_(SegmentEndPointFunc end_point_func, CouldReverseFunc could_reverse_func, size_t num_segments, const PointType *start_near) { @@ -118,6 +161,19 @@ std::vector> chain_segments_greedy_constrained_reversals out.emplace_back(0, could_reverse_func(0) && start_near != nullptr && (end_point_func(0, false) - *start_near).template cast().squaredNorm() < (end_point_func(0, true) - *start_near).template cast().squaredNorm()); } + else if (num_segments <= 4) + { + // Orca: so few segments are cheap to chain exactly. The first one is then not always the one nearest start_near. + std::vector end_points(2 * num_segments); + std::vector could_reverse(num_segments); + for (size_t i = 0; i < num_segments; ++ i) { + end_points[2 * i] = end_point_func(i, true).template cast(); + end_points[2 * i + 1] = end_point_func(i, false).template cast(); + could_reverse[i] = could_reverse_func(i); + } + const Vec2d start = start_near != nullptr ? Vec2d(start_near->template cast()) : Vec2d::Zero(); + out = chain_segments_exact(end_points, could_reverse, start_near != nullptr ? &start : nullptr); + } else { // End points of segments for the KD tree closest point search. @@ -1039,12 +1095,370 @@ std::vector> chain_segments_greedy2(SegmentEndPointFunc return chain_segments_greedy_constrained_reversals2_(end_point_func, could_reverse_func, num_segments, start_near); } +// The K points nearest to a query, nearest first, ties by index, from a uniform grid of about two points per cell. +template +class NearestPoints +{ +public: + static constexpr size_t npos = size_t(-1); + + explicit NearestPoints(const std::vector &points) : m_points(points) + { + BoundingBoxf bbox; + for (const Vec2d &p : points) + bbox.merge(p); + const Vec2d size = bbox.size(); + const double count = double(std::max(points.size(), 1)); + m_cell = std::max(std::sqrt(2. * size.x() * size.y() / count), std::max(size.x(), size.y()) / count); + if (m_cell <= 0.) + m_cell = 1.; + m_origin = bbox.min; + m_cols = int(size.x() / m_cell) + 1; + m_rows = int(size.y() / m_cell) + 1; + m_first.assign(size_t(m_cols) * m_rows + 1, 0); + for (const Vec2d &p : points) + ++ m_first[cell(p) + 1]; + for (size_t c = 1; c < m_first.size(); ++ c) + m_first[c] += m_first[c - 1]; + m_index.resize(points.size()); + std::vector next(m_first.begin(), m_first.end() - 1); + for (size_t i = 0; i < points.size(); ++ i) + m_index[next[cell(points[i])] ++] = i; + } + + // Indices of the nearest points that filter accepts, npos past the last one. + template + std::array query(const Vec2d &q, Filter filter) const + { + std::array, K> best; + size_t found = 0; + const int cx = int(std::clamp(std::floor((q.x() - m_origin.x()) / m_cell), 0., double(m_cols - 1))); + const int cy = int(std::clamp(std::floor((q.y() - m_origin.y()) / m_cell), 0., double(m_rows - 1))); + auto visit = [&](int x, int y) { + if (x < 0 || y < 0 || x >= m_cols || y >= m_rows) + return; + const size_t c = size_t(y) * m_cols + x; + for (size_t k = m_first[c]; k < m_first[c + 1]; ++ k) { + const size_t i = m_index[k]; + if (! filter(i)) + continue; + const std::pair candidate((m_points[i] - q).squaredNorm(), i); + if (found == K && ! (candidate < best[K - 1])) + continue; + auto it = std::upper_bound(best.begin(), best.begin() + found, candidate); + std::move_backward(it, best.begin() + std::min(found, K - 1), best.begin() + std::min(found + 1, K)); + *it = candidate; + found = std::min(found + 1, K); + } + }; + // Points beyond ring r are at least r cells away. + for (int r = 0; r <= std::max(m_cols, m_rows); ++ r) { + if (r == 0) + visit(cx, cy); + for (int d = -r; d <= r && r > 0; ++ d) { + visit(cx + d, cy - r); + visit(cx + d, cy + r); + if (d != -r && d != r) { + visit(cx - r, cy + d); + visit(cx + r, cy + d); + } + } + if (found == K && best[K - 1].first <= sqr(r * m_cell)) + break; + } + std::array out; + out.fill(npos); + for (size_t k = 0; k < found; ++ k) + out[k] = best[k].second; + return out; + } + +private: + size_t cell(const Vec2d &p) const + { + const int x = std::clamp(int((p.x() - m_origin.x()) / m_cell), 0, m_cols - 1); + const int y = std::clamp(int((p.y() - m_origin.y()) / m_cell), 0, m_rows - 1); + return size_t(y) * m_cols + x; + } + + const std::vector &m_points; + Vec2d m_origin; + double m_cell; + int m_cols; + int m_rows; + std::vector m_first; + std::vector m_index; +}; + +// Shortens the travel of a chain with 2-opt moves (reversing a run of segments) and or-opt moves (moving a run of up to +// three segments elsewhere, or a longer one to the front, flipped or not). end_points holds the first and the last point +// of each segment. The chain starts at start_near when given, and its end is free. A move that would flip a segment that +// cannot reverse is skipped. Moves only create links between nearest neighbor end points, and after a move only the +// segments next to a changed link are examined again, until a sweep of all segments finds no move. This keeps the search +// close to linear. +static void improve_chain_by_local_search(std::vector> &chain, const std::vector &end_points, + const std::vector &could_reverse, const Vec2d *start_near) +{ + const size_t n = chain.size(); + if (n < 2) + return; + assert(end_points.size() == 2 * n && could_reverse.size() == n); + + std::vector order(n), pos(n), fixed_before(n + 1, 0); + std::vector flipped(n); + for (size_t k = 0; k < n; ++ k) { + order[k] = chain[k].first; + flipped[order[k]] = chain[k].second; + } + // Updates positions lo..hi after a move, which keeps the number of segments that cannot reverse past hi. + auto refresh = [&](size_t lo, size_t hi) { + for (size_t k = lo; k <= hi; ++ k) { + pos[order[k]] = k; + fixed_before[k + 1] = fixed_before[k] + ! could_reverse[order[k]]; + } + }; + refresh(0, n - 1); + auto entry_idx = [&](size_t k) { return 2 * order[k] + flipped[order[k]]; }; + auto exit_idx = [&](size_t k) { return 2 * order[k] + 1 - flipped[order[k]]; }; + auto entry = [&](size_t k) -> const Vec2d& { return end_points[entry_idx(k)]; }; + auto exit = [&](size_t k) -> const Vec2d& { return end_points[exit_idx(k)]; }; + auto is_exit = [&](size_t e) { return char(e & 1) != flipped[e >> 1]; }; + auto reversible = [&](size_t i, size_t j) { return fixed_before[j + 1] == fixed_before[i]; }; + auto dist = [](const Vec2d *a, const Vec2d *b) { return a == nullptr || b == nullptr ? 0. : (*a - *b).norm(); }; + // Travel end points around position k, null where the chain is open. + auto prev_exit = [&](size_t k) -> const Vec2d* { return k > 0 ? &exit(k - 1) : start_near; }; + auto next_entry = [&](size_t k) -> const Vec2d* { return k + 1 < n ? &entry(k + 1) : nullptr; }; + + // The 32 nearest end points are the candidates, which covers every move on short chains; 8 missed much of the gain. + static constexpr size_t K = 32; + static constexpr size_t npos = NearestPoints::npos; + const NearestPoints nearest(end_points); + std::vector> neighbors(end_points.size()); + for (size_t e = 0; e < end_points.size(); ++ e) + neighbors[e] = nearest.query(end_points[e], [e](size_t idx) { return (idx >> 1) != (e >> 1); }); + std::array start_neighbors; + start_neighbors.fill(npos); + if (start_near != nullptr) + start_neighbors = nearest.query(*start_near, [](size_t) { return true; }); + + // Segments to examine, by id, each queued once. + std::deque queue; + std::vector queued(n, 0); + auto enqueue = [&](size_t k) { + if (k < n && ! queued[order[k]]) { + queued[order[k]] = 1; + queue.push_back(order[k]); + } + }; + + const double min_gain = SCALED_EPSILON; + // Reverses positions i..j if that shortens the chain. + auto try_2opt = [&](size_t i, size_t j) { + if (! reversible(i, j)) + return false; + const Vec2d *a = prev_exit(i), *b = next_entry(j); + if (dist(a, &exit(j)) + dist(&entry(i), b) > dist(a, &entry(i)) + dist(&exit(j), b) - min_gain) + return false; + std::reverse(order.begin() + i, order.begin() + j + 1); + for (size_t k = i; k <= j; ++ k) + flipped[order[k]] ^= 1; + refresh(i, j); + for (size_t k : { i - 1, i, j, j + 1 }) + enqueue(k); + return true; + }; + // Moves positions i..last after position gap (size_t(-1) for the front) if that shortens the chain. + auto try_oropt = [&](size_t i, size_t last, size_t gap) { + if (gap + 1 == i || (gap >= i && gap <= last)) + return false; + const Vec2d *a = gap == size_t(-1) ? start_near : &exit(gap); + const Vec2d *b = gap + 1 < n ? &entry(gap + 1) : nullptr; + const double removal_gain = dist(prev_exit(i), &entry(i)) + dist(&exit(last), next_entry(last)) - dist(prev_exit(i), next_entry(last)); + const double base = dist(a, b); + const double cost_kept = dist(a, &entry(i)) + dist(&exit(last), b) - base; + const double cost_flipped = reversible(i, last) ? dist(a, &exit(last)) + dist(&entry(i), b) - base : std::numeric_limits::max(); + if (std::min(cost_kept, cost_flipped) > removal_gain - min_gain) + return false; + const size_t len = last - i + 1; + // The run lands on first..first+len-1; lo..hi spans every position that changed. + const bool forward = gap > last && gap != size_t(-1); + const size_t first = forward ? gap + 1 - len : gap + 1; + const size_t lo = forward ? i : first; + const size_t hi = forward ? gap : last; + if (forward) + std::rotate(order.begin() + i, order.begin() + last + 1, order.begin() + gap + 1); + else + std::rotate(order.begin() + first, order.begin() + i, order.begin() + last + 1); + if (cost_flipped < cost_kept) { + std::reverse(order.begin() + first, order.begin() + first + len); + for (size_t k = first; k < first + len; ++ k) + flipped[order[k]] ^= 1; + } + refresh(lo, hi); + // The links around the run's new place, and the link closing its old place. + for (size_t k : { first - 1, first, first + len - 1, first + len, forward ? i - 1 : last, forward ? i : last + 1 }) + enqueue(k); + return true; + }; + // The gap before position k of the chain with positions i..last taken out. + auto gap_before = [](size_t k, size_t i, size_t last) { return k == last + 1 ? i - 1 : k - 1; }; + + // Applies the first move found that breaks a link next to position p. + // Candidates are tried nearest first, while closer than the link the move would replace at that end point. + auto closer = [&](const Vec2d &at, size_t q, double limit) { return q != npos && (end_points[q] - at).norm() < limit; }; + auto improve_at = [&](size_t p) { + // 2-opt breaking the link into p, as the first or the last link of the reversed run. + if (const Vec2d *a = prev_exit(p); a != nullptr) { + const double old_link = dist(a, &entry(p)); + for (size_t q : p > 0 ? neighbors[exit_idx(p - 1)] : start_neighbors) { + if (! closer(*a, q, old_link)) + break; + if (is_exit(q) && pos[q >> 1] >= p && try_2opt(p, pos[q >> 1])) + return true; + } + if (p > 0) + for (size_t q : neighbors[entry_idx(p)]) { + if (! closer(entry(p), q, old_link)) + break; + if (! is_exit(q) && pos[q >> 1] < p && try_2opt(pos[q >> 1], p - 1)) + return true; + } + } + // 2-opt breaking the link out of p. + if (p + 1 < n) { + const double old_link = dist(&exit(p), &entry(p + 1)); + for (size_t q : neighbors[exit_idx(p)]) { + if (! closer(exit(p), q, old_link)) + break; + if (is_exit(q) && pos[q >> 1] > p && try_2opt(p + 1, pos[q >> 1])) + return true; + } + for (size_t q : neighbors[entry_idx(p + 1)]) { + if (! closer(entry(p + 1), q, old_link)) + break; + if (! is_exit(q) && pos[q >> 1] <= p && try_2opt(pos[q >> 1], p)) + return true; + } + } + // A run of any length that starts or ends next to the start point moves to the front. + if (p == 0) + for (size_t q : start_neighbors) { + if (q == npos) + break; + const size_t m = pos[q >> 1]; + if (m == 0) + continue; + // Runs m..k enter next to the start, runs k..m leave there and go flipped. + const bool exit_near = is_exit(q); + for (size_t k = exit_near ? 1 : m; k < (exit_near ? m + 1 : n); ++ k) + if (try_oropt(exit_near ? k : m, exit_near ? m : k, size_t(-1))) + return true; + } + // Or-opt of the runs starting at p: next to an exit a run goes after that segment, next to an entry before it. + for (size_t last = p; last < std::min(p + 3, n); ++ last) { + if (try_oropt(p, last, size_t(-1))) + return true; + const double removal_gain = dist(prev_exit(p), &entry(p)) + dist(&exit(last), next_entry(last)) - dist(prev_exit(p), next_entry(last)); + for (size_t e : { entry_idx(p), exit_idx(last) }) + for (size_t q : neighbors[e]) { + if (! closer(end_points[e], q, removal_gain)) + break; + const size_t m = pos[q >> 1]; + if ((m < p || m > last) && try_oropt(p, last, is_exit(q) ? m : gap_before(m, p, last))) + return true; + } + } + return false; + }; + + // A new link can make room for a run that sits elsewhere, hence the sweeps until one finds no move. Every move shortens + // the chain, so the search ends; the cap only bounds its time. + for (size_t k = 0; k < n; ++ k) + enqueue(k); + bool moved_since_sweep = false; + for (size_t moves = 0; moves < 50 * n;) { + if (queue.empty()) { + if (! moved_since_sweep) + break; + moved_since_sweep = false; + for (size_t k = 0; k < n; ++ k) + enqueue(k); + } + const size_t segment = queue.front(); + queue.pop_front(); + queued[segment] = 0; + if (improve_at(pos[segment])) { + ++ moves; + moved_since_sweep = true; + enqueue(pos[segment]); + } + } + + for (size_t k = 0; k < n; ++ k) + chain[k] = { order[k], flipped[order[k]] != 0 }; +} + +// Chain that always goes on to the nearest free end point, from start. +static std::vector> chain_nearest_next(const std::vector &end_points, const std::vector &could_reverse, const Vec2d &start) +{ + const NearestPoints<1> nearest(end_points); + std::vector taken(could_reverse.size(), 0); + std::vector> out; + out.reserve(could_reverse.size()); + for (Vec2d at = start; out.size() < could_reverse.size();) { + const size_t e = nearest.query(at, [&](size_t idx) { return ! taken[idx >> 1] && ((idx & 1) == 0 || could_reverse[idx >> 1]); })[0]; + taken[e >> 1] = 1; + out.emplace_back(e >> 1, (e & 1) != 0); + at = end_points[e ^ 1]; + } + return out; +} + +static double chain_travel(const std::vector> &chain, const std::vector &end_points, const Vec2d &start) +{ + double travel = 0.; + Vec2d at = start; + for (const std::pair &segment : chain) { + travel += (end_points[2 * segment.first + segment.second] - at).norm(); + at = end_points[2 * segment.first + 1 - segment.second]; + } + return travel; +} + template -static std::vector> chain_extrusion_entities_impl(const std::vector &entities, const Point *start_near) +static std::vector> chain_extrusion_entities_impl(const std::vector &entities, const Point *start_near, bool improve = false) { auto segment_end_point = [&entities](size_t idx, bool first_point) -> Point { return first_point ? entities[idx]->first_point() : entities[idx]->last_point(); }; auto could_reverse = [&entities](size_t idx) { const ExtrusionEntity *ee = entities[idx]; return ee->is_loop() || ee->can_reverse(); }; std::vector> out = chain_segments_greedy_constrained_reversals(segment_end_point, could_reverse, entities.size(), start_near); + if (improve && out.size() > 1) { + std::vector end_points; + std::vector reversible; + end_points.reserve(2 * entities.size()); + reversible.reserve(entities.size()); + for (size_t idx = 0; idx < entities.size(); ++ idx) { + end_points.emplace_back(segment_end_point(idx, true).template cast()); + end_points.emplace_back(segment_end_point(idx, false).template cast()); + reversible.emplace_back(could_reverse(idx)); + } + Vec2d start = Vec2d::Zero(); + if (start_near != nullptr) + start = start_near->cast(); + if (entities.size() <= 8) { + // Few paths are chained exactly; the local search below can stop short of the shortest chain. + out = chain_segments_exact(end_points, reversible, start_near != nullptr ? &start : nullptr); + } else { + improve_chain_by_local_search(out, end_points, reversible, start_near != nullptr ? &start : nullptr); + if (start_near != nullptr) { + // The greedy chain can settle into a sweep that local moves cannot undo, which a nearest-next chain may avoid. + std::vector> nearest_next = chain_nearest_next(end_points, reversible, start); + if (chain_travel(nearest_next, end_points, start) < chain_travel(out, end_points, start) - SCALED_EPSILON) { + improve_chain_by_local_search(nearest_next, end_points, reversible, &start); + out = std::move(nearest_next); + } + } + } + } for (std::pair &segment : out) { const ExtrusionEntity *ee = entities[segment.first]; if (ee->is_loop()) @@ -1086,7 +1500,8 @@ void reorder_extrusion_entities(std::vector &entities, const s void chain_and_reorder_extrusion_entities(std::vector &entities, const Point &start_near) { - chain_and_reorder_extrusion_entities(entities, &start_near); + remove_entities_without_endpoints(entities); + reorder_extrusion_entities(entities, chain_extrusion_entities_impl(entities, &start_near, true)); } void chain_and_reorder_extrusion_entities(std::vector &entities, const Point *start_near) @@ -1101,7 +1516,7 @@ void chain_and_reorder_extrusion_entities(std::vector &e remove_entities_without_endpoints(entities); std::vector out; out.reserve(entities.size()); - for (const auto &[idx, reverse] : chain_extrusion_entities_impl(entities, &start_near)) { + for (const auto &[idx, reverse] : chain_extrusion_entities_impl(entities, &start_near, true)) { if (reverse) { ExtrusionEntity *clone = reversed_clones.emplace_back(entities[idx]->clone()).get(); clone->reverse(); @@ -1159,16 +1574,74 @@ std::vector chain_points(const Points &points, const Point *start_near) return out; } +// Moves runs of up to three points of the cycle next to one of the 8 points nearest to an end of the run, reversed or not, +// while that shortens the cycle. path[0] stays first. Returns true if any move was made. +static bool improve_cycle_by_or_opt(std::vector &path, const Points &points) +{ + const size_t pn = path.size(); + if (pn < 5) + return false; + std::vector coords(points.size()); + for (size_t i = 0; i < points.size(); ++ i) + coords[i] = points[i].cast(); + const NearestPoints<8> grid(coords); + std::vector pos(points.size()); + auto refresh = [&]() { for (size_t k = 0; k < pn; ++ k) pos[path[k]] = k; }; + refresh(); + auto dist = [&](size_t a, size_t b) { return (coords[path[a % pn]] - coords[path[b % pn]]).norm(); }; + bool improved = false; + for (bool moved = true; moved;) { + moved = false; + for (size_t len = 1; len <= 3; ++ len) + for (size_t i = 1; i + len <= pn; ++ i) { + // The run path[i..last] goes between path[k] and path[k + 1]. + const size_t last = i + len - 1; + const double removal = dist(i - 1, i) + dist(last, last + 1) - dist(i - 1, last + 1); + bool applied = removal <= SCALED_EPSILON; + for (size_t end : { path[i], path[last] }) + for (size_t q : grid.query(coords[end], [end](size_t idx) { return idx != end; })) + for (size_t k : { q == NearestPoints<8>::npos ? pn : (pos[q] + pn - 1) % pn, q == NearestPoints<8>::npos ? pn : pos[q] }) { + if (applied || k == pn || (k + 1 >= i && k <= last)) + continue; + const double base = dist(k, k + 1); + const double kept = dist(k, i) + dist(last, k + 1) - base; + const double reversed = dist(k, last) + dist(i, k + 1) - base; + if (std::min(kept, reversed) > removal - SCALED_EPSILON) + continue; + std::vector run(path.begin() + i, path.begin() + last + 1); + if (reversed < kept) + std::reverse(run.begin(), run.end()); + path.erase(path.begin() + i, path.begin() + last + 1); + path.insert(path.begin() + (k < i ? k + 1 : k + 1 - len), run.begin(), run.end()); + refresh(); + applied = moved = improved = true; + } + } + } + return improved; +} + std::vector chain_points_with_postprocessing(const Points &points, const Point *start_near) { - std::vector path = chain_points(points, start_near); - // Alternate 2-opt and crossing removal until convergence. - // 2-opt can create new crossings, and crossing removal can create new - // opportunities for 2-opt improvement. Break early if neither improves. - for (int iter = 0; iter < 3; ++iter) { - bool improved = tsp_2opt_improve(path, points); - improved |= tsp_remove_crossings(path, points); - if (!improved) break; + std::vector path; + if (points.size() <= 10) { + // Orca: few points are toured exactly, from the point nearest start_near toward its nearer neighbor. + path = tsp_exact_cycle(points); + if (start_near != nullptr && path.size() > 2) { + auto dist = [&points](size_t a, const Point &b) { return (points[a] - b).cast().squaredNorm(); }; + std::rotate(path.begin(), std::min_element(path.begin(), path.end(), [&](size_t a, size_t b) { return dist(a, *start_near) < dist(b, *start_near); }), path.end()); + if (dist(path.back(), points[path.front()]) < dist(path[1], points[path.front()])) + std::reverse(path.begin() + 1, path.end()); + } + } else { + path = chain_points(points, start_near); + // Alternate 2-opt, or-opt and crossing removal: each can open moves for the others. Break early if none improves. + for (int iter = 0; iter < 3; ++iter) { + bool improved = tsp_2opt_improve(path, points); + improved |= improve_cycle_by_or_opt(path, points); + improved |= tsp_remove_crossings(path, points); + if (!improved) break; + } } if (start_near == nullptr) tsp_rotate_minimize_closing(path, points); diff --git a/src/libslic3r/ShortestPath.hpp b/src/libslic3r/ShortestPath.hpp index 53ec00f9e4..34127fdb9a 100644 --- a/src/libslic3r/ShortestPath.hpp +++ b/src/libslic3r/ShortestPath.hpp @@ -23,9 +23,12 @@ std::vector chain_expolygons(const ExPolygons &input_exploy); std::vector> chain_extrusion_entities(std::vector &entities, const Point *start_near = nullptr); void reorder_extrusion_entities(std::vector &entities, const std::vector> &chain); +// From start_near, the greedy chain is refined by 2-opt and or-opt moves, which shortens the travel at some extra cost. void chain_and_reorder_extrusion_entities(std::vector &entities, const Point &start_near); +// The greedy chain only. void chain_and_reorder_extrusion_entities(std::vector &entities, const Point *start_near = nullptr); -// Each entity the chain reverses is replaced by a reversed clone that reversed_clones owns, so the originals stay unchanged. +// Refined as above. Each entity the chain reverses is replaced by a reversed clone that reversed_clones owns, so the +// originals stay unchanged. void chain_and_reorder_extrusion_entities(std::vector &entities, const Point &start_near, std::vector> &reversed_clones); diff --git a/tests/fff_print/test_fill.cpp b/tests/fff_print/test_fill.cpp index 7ca8487965..b2d0ef9abe 100644 --- a/tests/fff_print/test_fill.cpp +++ b/tests/fff_print/test_fill.cpp @@ -1719,9 +1719,10 @@ TEST_CASE("Smoothing multiline lightning infill keeps its outlines connected", " const SparseInfillShape smooth = shape_for("100%"); REQUIRE(sharp.path_count > 0); - // Smoothing should leave the loop count about where it was; uncapping the smoothing reach, the - // regression this guards against, adds about 10%. - const size_t allowed_extra = sharp.path_count / 50; // 2% + // Smoothing should leave the path count about where it was; uncapping the smoothing reach, the + // regression this guards against, adds a quarter to a third. Outline ends next to each other are + // joined into one path, so a few loops more or less are a larger share of the count. + const size_t allowed_extra = sharp.path_count / 10; // 10% REQUIRE(smooth.path_count <= sharp.path_count + allowed_extra); // The outlines are still rounded. REQUIRE(smooth.point_count > sharp.point_count); @@ -2353,3 +2354,66 @@ TEST_CASE("Adaptive TPMS anchors match the printed infill", "[Fill][InternalBrid max_distance = std::max(max_distance, printed_tree.distance_from_lines(point)); CHECK(unscale(max_distance) <= 0.012); } + +TEST_CASE("Infill lines whose ends are next to each other join into one path", "[Fill]") +{ + // An out-and-back pair 1 mm apart, as multiline infill leaves, and a line 19 mm away. + Polylines lines { Polyline(Point::new_scale(0., 0.), Point::new_scale(10., 0.)), Polyline(Point::new_scale(0., 1.), Point::new_scale(10., 1.)), + Polyline(Point::new_scale(0., 20.), Point::new_scale(10., 20.)) }; + join_close_ends(lines, scale_(1.5), true); + + REQUIRE(lines.size() == 2); + const auto joined = std::find_if(lines.begin(), lines.end(), [](const Polyline &pl) { return pl.size() == 4; }); + REQUIRE(joined != lines.end()); + // Both lines plus the 1 mm between them, joined at one end only, as joining the other end too would close a loop. + CHECK_THAT(unscale(joined->length()), Catch::Matchers::WithinAbs(21., 1e-6)); + CHECK(joined->first_point() != joined->last_point()); +} + +TEST_CASE("Joined infill lines print each line once", "[Fill]") +{ + // The middle line comes first and both of its ends get joined, so the group has to be walked from an outer line. + Polylines lines { Polyline(Point::new_scale(0., 0.), Point::new_scale(10., 0.)), Polyline(Point::new_scale(0., 1.), Point::new_scale(10., 1.)), + Polyline(Point::new_scale(10., -1.), Point::new_scale(0., -1.)) }; + join_close_ends(lines, scale_(1.5), true); + + REQUIRE(lines.size() == 1); + CHECK(lines.front().size() == 6); + // Three 10 mm lines and two 1 mm joins. + CHECK_THAT(unscale(lines.front().length()), Catch::Matchers::WithinAbs(32., 1e-6)); +} + +TEST_CASE("Infill line ends are not joined across another line", "[Fill]") +{ + // A pair 1 mm apart whose left ends have a third line running between them, which a join there would cross. + Polylines lines { Polyline(Point::new_scale(0., 0.), Point::new_scale(10., 0.)), Polyline(Point::new_scale(0., 1.), Point::new_scale(10., 1.)), + Polyline(Point::new_scale(-5., 0.5), Point::new_scale(5., 0.5)) }; + join_close_ends(lines, scale_(1.5), true); + + REQUIRE(lines.size() == 2); + const auto joined = std::find_if(lines.begin(), lines.end(), [](const Polyline &pl) { return pl.size() == 4; }); + REQUIRE(joined != lines.end()); + // Joined on the right, where nothing is in between, so the left ends stay open. + CHECK(joined->first_point().x() == 0); + CHECK(joined->last_point().x() == 0); +} + +TEST_CASE("Infill lines that cannot reverse join only end to start", "[Fill]") +{ + const Point a(Point::new_scale(0., 0.)), b(Point::new_scale(10., 0.)), c(Point::new_scale(10., 1.)), d(Point::new_scale(0., 1.)); + SECTION("Ends that would need a reversal stay apart") { + Polylines lines { Polyline(a, b), Polyline(d, c) }; + join_close_ends(lines, scale_(1.5), false); + CHECK(lines.size() == 2); + } + SECTION("The end of one line joins the start of the next, both kept in their direction") { + Polylines lines { Polyline(c, d), Polyline(a, b) }; + join_close_ends(lines, scale_(1.5), false); + REQUIRE(lines.size() == 1); + const Points &points = lines.front().points; + REQUIRE(points.size() == 4); + auto position = [&points](const Point &p) { return std::find(points.begin(), points.end(), p) - points.begin(); }; + CHECK(position(a) < position(b)); + CHECK(position(c) < position(d)); + } +} diff --git a/tests/libslic3r/CMakeLists.txt b/tests/libslic3r/CMakeLists.txt index 5c34b151b2..518a336bc5 100644 --- a/tests/libslic3r/CMakeLists.txt +++ b/tests/libslic3r/CMakeLists.txt @@ -67,6 +67,7 @@ add_executable(${_TEST_NAME}_tests test_wipe_path.cpp test_optimizers.cpp test_ordering_strategies.cpp + test_shortest_path.cpp # test_png_io.cpp test_indexed_triangle_set.cpp test_connected_bodies.cpp diff --git a/tests/libslic3r/test_ordering_strategies.cpp b/tests/libslic3r/test_ordering_strategies.cpp index 97e5a83f60..f4df7e8b5e 100644 --- a/tests/libslic3r/test_ordering_strategies.cpp +++ b/tests/libslic3r/test_ordering_strategies.cpp @@ -1,4 +1,5 @@ #include +#include #include #include #include @@ -12,10 +13,13 @@ #include "libslic3r/Point.hpp" #include "libslic3r/GCode/OrderingStrategies.hpp" +#include "libslic3r/ShortestPath.hpp" #include "libslic3r/Geometry.hpp" #include #include +#include +#include #include using namespace Slic3r; @@ -304,3 +308,60 @@ TEST_CASE("snake_core handles three points", "[Snake]") { REQUIRE(is_permutation(p2, 3)); } + +// The shortest cycle through the points, trying every order. +static double shortest_cycle_length(const Points& centers) +{ + std::vector path(centers.size()); + std::iota(path.begin(), path.end(), 0); + double best = std::numeric_limits::max(); + do { + best = std::min(best, tsp_cycle_path_length(path, centers)); + } while (std::next_permutation(path.begin() + 1, path.end())); + return best; +} + +static Points random_centers(unsigned seed, size_t count) +{ + std::mt19937 rng(seed); + std::uniform_real_distribution<> coordinate(0., 200.); + Points centers; + for (size_t i = 0; i < count; ++i) { + // One draw per statement, as the evaluation order of function arguments is unspecified. + const double x = coordinate(rng); + const double y = coordinate(rng); + centers.emplace_back(Point::new_scale(x, y)); + } + return centers; +} + +TEST_CASE("tsp_exact_cycle finds the shortest cycle", "[TSPPostProcessing]") { + const size_t count = GENERATE(4, 5, 6, 8); + for (unsigned seed = 1; seed <= 10; ++seed) { + const Points centers = random_centers(seed, count); + const std::vector path = tsp_exact_cycle(centers); + REQUIRE(is_permutation(path, count)); + CHECK(euclidean_path_length(path, centers) <= shortest_cycle_length(centers) + SCALED_EPSILON); + } +} + +TEST_CASE("Default ordering tours up to ten points along the shortest cycle", "[TSPPostProcessing]") { + const bool from_start = GENERATE(false, true); + const Point start = Point::new_scale(-20., -20.); + for (unsigned seed = 1; seed <= 5; ++seed) { + const Points centers = random_centers(seed, 9); + const std::vector path = chain_points_with_postprocessing(centers, from_start ? &start : nullptr); + REQUIRE(is_permutation(path, centers.size())); + CHECK(euclidean_path_length(path, centers) <= shortest_cycle_length(centers) + SCALED_EPSILON); + } +} + +TEST_CASE("Default ordering of many points is never longer than the greedy cycle", "[TSPPostProcessing]") { + for (unsigned seed = 1; seed <= 10; ++seed) { + const Points centers = random_centers(seed, 150); + const std::vector path = chain_points_with_postprocessing(centers, nullptr); + REQUIRE(is_permutation(path, centers.size())); + CHECK(euclidean_path_length(path, centers) <= euclidean_path_length(chain_points(centers, nullptr), centers) + SCALED_EPSILON); + CHECK_FALSE(has_crossings(path, centers)); + } +} diff --git a/tests/libslic3r/test_shortest_path.cpp b/tests/libslic3r/test_shortest_path.cpp new file mode 100644 index 0000000000..e4b8c60ddc --- /dev/null +++ b/tests/libslic3r/test_shortest_path.cpp @@ -0,0 +1,202 @@ +#include +#include +#include +#include + +#include "libslic3r/ExtrusionEntity.hpp" +#include "libslic3r/ExtrusionEntityCollection.hpp" +#include "libslic3r/Point.hpp" +#include "libslic3r/Polyline.hpp" +#include "libslic3r/ShortestPath.hpp" +#include "libslic3r/libslic3r.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using namespace Slic3r; + +namespace { + +// Infill-like paths: lines of 1 to 30 mm scattered over a 100 mm square. The collection owns them, and entities hands +// them out as the raw pointers chain_and_reorder_extrusion_entities() reorders. +struct Paths +{ + ExtrusionEntityCollection collection; + std::vector entities; + + Paths(unsigned seed, size_t count, bool can_reverse) + { + std::mt19937 rng(seed); + std::uniform_real_distribution<> coordinate(0., 100.), length(1., 30.), angle(0., 2. * M_PI); + Polylines lines; + for (size_t i = 0; i < count; ++ i) { + // One draw per statement, as the evaluation order of function arguments is unspecified. + const double x = coordinate(rng); + const double y = coordinate(rng); + const double l = length(rng); + const double t = angle(rng); + const Vec2d a(x, y); + const Vec2d b = a + l * Vec2d(std::cos(t), std::sin(t)); + lines.emplace_back(Point::new_scale(a), Point::new_scale(b)); + } + extrusion_entities_append_paths(collection.entities, std::move(lines), erInternalInfill, 0.1, 0.45f, 0.2f, can_reverse); + entities = collection.entities; + } +}; + +double travel(const std::vector &entities, const Point &start) +{ + double length = 0.; + Point last = start; + for (const ExtrusionEntity *ee : entities) { + length += (ee->first_point() - last).cast().norm(); + last = ee->last_point(); + } + return length; +} + +// Travel of the chain that always goes on to the nearest free path end. +double nearest_next_travel(const std::vector &entities, bool can_reverse, const Point &from) +{ + std::vector taken(entities.size(), 0); + double length = 0.; + Point last = from; + for (size_t step = 0; step < entities.size(); ++ step) { + double best = std::numeric_limits::max(); + size_t best_idx = 0; + bool best_reversed = false; + for (size_t i = 0; i < entities.size(); ++ i) + for (const bool reversed : { false, true }) + if (! taken[i] && (can_reverse || ! reversed)) { + const double d = ((reversed ? entities[i]->last_point() : entities[i]->first_point()) - last).cast().norm(); + if (d < best) { + best = d; + best_idx = i; + best_reversed = reversed; + } + } + taken[best_idx] = 1; + length += best; + last = best_reversed ? entities[best_idx]->first_point() : entities[best_idx]->last_point(); + } + return length; +} + +// The shortest travel over every order and allowed reversal of the paths, from start. +double shortest_travel(const std::vector &entities, bool can_reverse, const Point &from) +{ + std::vector order(entities.size()); + std::iota(order.begin(), order.end(), 0); + double best = std::numeric_limits::max(); + do { + for (size_t mask = 0; mask < (size_t(1) << order.size()) && (can_reverse || mask == 0); ++ mask) { + double length = 0.; + Point last = from; + for (size_t k = 0; k < order.size(); ++ k) { + const ExtrusionEntity *ee = entities[order[k]]; + const bool reversed = (mask >> k & 1) != 0; + length += ((reversed ? ee->last_point() : ee->first_point()) - last).cast().norm(); + last = reversed ? ee->first_point() : ee->last_point(); + } + best = std::min(best, length); + } + } while (std::next_permutation(order.begin(), order.end())); + return best; +} + +const Point start = Point::new_scale(-10., -10.); + +} // namespace + +TEST_CASE("An improved infill chain is never longer than the greedy or nearest-next one and keeps every path", "[ShortestPath]") +{ + const size_t count = GENERATE(2, 3, 10, 60, 300); + const bool can_reverse = GENERATE(true, false); + double greedy_total = 0., improved_total = 0.; + for (unsigned seed = 1; seed <= 20; ++ seed) { + Paths greedy(seed, count, can_reverse), improved(seed, count, can_reverse); + const double nearest_travel = nearest_next_travel(improved.entities, can_reverse, start); + chain_and_reorder_extrusion_entities(greedy.entities, &start); + chain_and_reorder_extrusion_entities(improved.entities, start); + + std::vector sorted = improved.entities; + std::sort(sorted.begin(), sorted.end()); + REQUIRE(std::unique(sorted.begin(), sorted.end()) == sorted.end()); + REQUIRE(sorted.size() == count); + + const double greedy_travel = travel(greedy.entities, start), improved_travel = travel(improved.entities, start); + CHECK(improved_travel <= greedy_travel + SCALED_EPSILON); + CHECK(improved_travel <= nearest_travel + SCALED_EPSILON); + greedy_total += greedy_travel; + improved_total += improved_travel; + } + // Twenty random sets of ten or more paths always leave the greedy chain something to shorten. + if (count >= 10) + CHECK(improved_total < 0.95 * greedy_total); +} + +TEST_CASE("Improving an infill chain keeps the direction of paths that cannot reverse", "[ShortestPath]") +{ + Paths paths(7, 100, false); + std::map first_points; + for (const ExtrusionEntity *ee : paths.entities) + first_points.emplace(ee, ee->first_point()); + + chain_and_reorder_extrusion_entities(paths.entities, start); + + for (const ExtrusionEntity *ee : paths.entities) + CHECK(ee->first_point() == first_points.at(ee)); +} + +TEST_CASE("Chains of a few paths are the shortest possible", "[ShortestPath]") +{ + const size_t count = GENERATE(2, 3, 4, 5, 6); + const bool can_reverse = GENERATE(true, false); + for (unsigned seed = 1; seed <= 10; ++ seed) { + Paths paths(seed, count, can_reverse); + const double shortest = shortest_travel(paths.entities, can_reverse, start); + chain_and_reorder_extrusion_entities(paths.entities, start); + CHECK(travel(paths.entities, start) <= shortest + SCALED_EPSILON); + } +} + +TEST_CASE("The greedy chain of two to four paths is the shortest possible", "[ShortestPath]") +{ + const size_t count = GENERATE(2, 3, 4); + const bool can_reverse = GENERATE(true, false); + for (unsigned seed = 1; seed <= 10; ++ seed) { + Paths paths(seed, count, can_reverse); + const double shortest = shortest_travel(paths.entities, can_reverse, start); + chain_and_reorder_extrusion_entities(paths.entities, &start); + CHECK(travel(paths.entities, start) <= shortest + SCALED_EPSILON); + } +} + +TEST_CASE("A short chain need not start with the path nearest to its start", "[ShortestPath]") +{ + // The nearest end leads 10 mm away from the two other paths, so the shortest chain leaves that path for last. + const Point from = Point::new_scale(0., 0.); + for (const bool greedy : { false, true }) { + ExtrusionEntityCollection collection; + extrusion_entities_append_paths(collection.entities, + Polylines{ Polyline(Point::new_scale(-0.9, 0.), Point::new_scale(-10., 0.)), Polyline(Point::new_scale(1., 0.), Point::new_scale(2., 0.)), + Polyline(Point::new_scale(3., 1.), Point::new_scale(4., 1.)) }, + erInternalInfill, 0.1, 0.45f, 0.2f); + std::vector entities = collection.entities; + const double shortest = shortest_travel(entities, true, from); + if (greedy) + chain_and_reorder_extrusion_entities(entities, &from); + else + chain_and_reorder_extrusion_entities(entities, from); + CHECK(entities.back() == collection.entities.front()); + CHECK_THAT(travel(entities, from), Catch::Matchers::WithinAbs(shortest, SCALED_EPSILON)); + } +}