Lagrange
chain_edges_into_simple_loops.h
1/*
2 * Copyright 2020 Adobe. All rights reserved.
3 * This file is licensed to you under the Apache License, Version 2.0 (the "License");
4 * you may not use this file except in compliance with the License. You may obtain a copy
5 * of the License at http://www.apache.org/licenses/LICENSE-2.0
6 *
7 * Unless required by applicable law or agreed to in writing, software distributed under
8 * the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR REPRESENTATIONS
9 * OF ANY KIND, either express or implied. See the License for the specific language
10 * governing permissions and limitations under the License.
11 */
12#pragma once
13
14#include <lagrange/Logger.h>
15#include <lagrange/common.h>
16#include <lagrange/legacy/inline.h>
17
18#include <Eigen/Dense>
19
20#include <fstream>
21#include <stack>
22#include <string>
23#include <vector>
24
25namespace lagrange {
26LAGRANGE_LEGACY_INLINE
27namespace legacy {
28
45template <typename DerivedE, typename DerivedO>
46bool chain_edges_into_simple_loops(
47 const Eigen::MatrixBase<DerivedE>& edges,
48 std::vector<std::vector<typename DerivedE::Scalar>>& output,
49 Eigen::PlainObjectBase<DerivedO>& remaining_edges)
50{
51 using Index = typename DerivedE::Scalar;
52 using IndexArray = Eigen::Matrix<Index, Eigen::Dynamic, 1>;
53
54 output.clear();
55
56 if (edges.rows() == 0) {
57 // Empty graph
58 return true;
59 }
60
61 Index num_vertices = edges.maxCoeff() + 1;
62 Index num_edges = (Index)edges.rows();
63
64 // Count degree_in and degree_out, and make sure that it matches
65 IndexArray degree_in(num_vertices);
66 IndexArray degree_out(num_vertices);
67 degree_in.setZero();
68 degree_out.setZero();
69 for (Index e = 0; e < num_edges; ++e) {
70 Index v0 = edges(e, 0);
71 Index v1 = edges(e, 1);
72 degree_out[v0]++;
73 degree_in[v1]++;
74 }
75 if (!(degree_in.array() == degree_out.array()).all()) {
76 logger().debug("Input digraph has dangling vertices.");
77 }
78
79 // path -> first edge in the path
80 std::vector<Index> path_to_first_edge;
81
82 // vertex -> single outgoing edge along path (or INVALID if degree_out is > 1)
83 std::vector<Index> vertex_to_outgoing_edge(num_vertices, invalid<Index>());
84
85 // edge -> next edge along path (or INVALID if last edge along path)
86 std::vector<Index> next_edge_along_path(num_edges, invalid<Index>());
87
88 // Chain edges into paths
89 for (Index e = 0; e < num_edges; ++e) {
90 Index v0 = edges(e, 0);
91 if (degree_out[v0] == 1 && degree_in[v0] == 1) {
92 // v0 is a mid-path vertex. There is only one possibility for the next edge along a path
93 // going through v0.
94 vertex_to_outgoing_edge[v0] = e;
95 } else {
96 // v0 is a junction vertex. We start one path for each outgoing edge e starting from
97 // vertex v0.
98 path_to_first_edge.push_back(e);
99 }
100 }
101 // Chain edges together based on the vertex -> outgoing edge mapping we previously computed
102 for (Index e = 0; e < num_edges; ++e) {
103 Index v1 = edges(e, 1);
104 next_edge_along_path[e] = vertex_to_outgoing_edge[v1];
105 }
106
107 // Follow each path until we reach the last edge
108 std::vector<Index> edge_label(num_edges, invalid<Index>());
109 std::vector<Index> path_to_last_edge(path_to_first_edge.size());
110 std::stack<Index> ears;
111 std::vector<std::vector<Index>> paths_in(num_vertices);
112 std::vector<std::vector<Index>> paths_out(num_vertices);
113 std::vector<bool> path_is_pending(path_to_first_edge.size(), false);
114
115 auto first_vertex_in_path = [&](Index a) { return edges(path_to_first_edge[a], 0); };
116
117 auto last_vertex_in_path = [&](Index a) { return edges(path_to_last_edge[a], 1); };
118
119 // For each path that we have started at a junction vertex, follow edges along the path until
120 // each edge has been labeled as belonging to the path. Note that our paths do not contain
121 // junction vertices by construction. So the only case where the path can contain a cycle is
122 // when the path starts and ends at the same vertex.
123 for (Index a = 0; a < (Index)path_to_first_edge.size(); ++a) {
124 // Follow edges along the path and label them as belonging to the path
125 for (Index e = path_to_first_edge[a];
126 e != invalid<Index>() && edge_label[e] == invalid<Index>();
127 e = next_edge_along_path[e]) {
128 edge_label[e] = a;
129 path_to_last_edge[a] = e;
130 }
131
132 const Index v_first = first_vertex_in_path(a);
133 const Index v_last = last_vertex_in_path(a);
134
135 // Compute outgoing and ingoing paths per vertices of degree > 1
136 paths_out[v_first].push_back(a);
137 paths_in[v_last].push_back(a);
138
139 if (v_first == v_last) {
140 // Path is an ear (simple loop), will be popped next
141 la_debug_assert(path_is_pending[a] == false);
142 path_is_pending[a] = true;
143 ears.push(a);
144 }
145 }
146
147 // For paths which are isolated cycles, there is no "starting vertex" (each vertex has total
148 // degree 2). We do an additional pass on each edge and start a new path for each unlabeled
149 // edge.
150 for (Index e = 0; e < num_edges; ++e) {
151 if (edge_label[e] == invalid<Index>()) {
152 // Add a new path starting here
153 const Index a = (Index)path_to_first_edge.size();
154 path_to_first_edge.emplace_back(e);
155 path_to_last_edge.emplace_back(e);
156 path_is_pending.push_back(false);
157
158 // Follow edges along the path and record the last edge on the path
159 edge_label[e] = a;
160 for (Index ei = next_edge_along_path[e];
161 ei != invalid<Index>() && edge_label[ei] == invalid<Index>();
162 ei = next_edge_along_path[ei]) {
163 edge_label[ei] = a;
164 path_to_last_edge[a] = ei;
165 }
166 la_debug_assert(next_edge_along_path[path_to_last_edge[a]] == e);
167 la_debug_assert(first_vertex_in_path(a) == last_vertex_in_path(a));
168 next_edge_along_path[path_to_last_edge[a]] = invalid<Index>();
169
170 // Path is an isolated cycle
171 la_debug_assert(path_is_pending[a] == false);
172 path_is_pending[a] = true;
173 ears.push(a);
174 }
175 }
176
177 // Swap/pop_back trick to remove an element from a std::vector
178 auto remove_from_vector = [](std::vector<Index>& vec, size_t idx) {
179 la_debug_assert(idx < vec.size());
180 std::swap(vec[idx], vec.back());
181 vec.pop_back();
182 };
183
184 Index num_edges_removed = 0;
185 std::vector<bool> edge_is_removed(num_edges, false);
186
187 // Pop ears repeatedly
188 std::vector<bool> path_is_removed(path_to_first_edge.size(), false);
189 while (!ears.empty()) {
190 const Index a = ears.top();
191 ears.pop();
192 la_debug_assert(!path_is_removed[a]);
193 path_is_removed[a] = true;
194
195 // path starts and ends on the same vertex, compute corresponding (simple) loop
196 la_debug_assert(first_vertex_in_path(a) == last_vertex_in_path(a));
197 la_debug_assert(path_to_first_edge[a] != invalid<Index>());
198 std::vector<Index> loop;
199 for (Index e = path_to_first_edge[a]; e != invalid<Index>(); e = next_edge_along_path[e]) {
200 loop.push_back(e);
201 la_debug_assert(edge_is_removed[e] == false);
202 edge_is_removed[e] = true;
203 ++num_edges_removed;
204 }
205 output.emplace_back(std::move(loop));
206
207 // Remove current path from the in/out paths of the endpoint vertex v
208 const Index v = first_vertex_in_path(a);
209 la_debug_assert(v == last_vertex_in_path(a));
210 for (size_t i = 0; i < paths_out[v].size();) {
211 if (paths_out[v][i] == a) {
212 remove_from_vector(paths_out[v], i);
213 } else {
214 la_debug_assert(!path_is_removed[paths_out[v][i]]);
215 ++i;
216 }
217 }
218 for (size_t i = 0; i < paths_in[v].size();) {
219 if (paths_in[v][i] == a) {
220 remove_from_vector(paths_in[v], i);
221 } else {
222 la_debug_assert(!path_is_removed[paths_in[v][i]]);
223 ++i;
224 }
225 }
226
227 // If there are only 1 remaining in/out paths, join them
228 if (paths_in[v].size() == 1 && paths_out[v].size() == 1) {
229 const Index a_in = paths_in[v].front();
230 const Index a_out = paths_out[v].front();
231 la_debug_assert(last_vertex_in_path(a_in) == v);
232 la_debug_assert(first_vertex_in_path(a_out) == v);
233 if (a_in != a_out) {
234 // If the paths are different, then join them
235 la_debug_assert(path_to_first_edge[a_in] != invalid<Index>());
236 la_debug_assert(path_to_last_edge[a_in] != invalid<Index>());
237 la_debug_assert(path_to_first_edge[a_out] != invalid<Index>());
238 la_debug_assert(path_to_last_edge[a_out] != invalid<Index>());
239 la_debug_assert(next_edge_along_path[path_to_last_edge[a_in]] == invalid<Index>());
241 edges(path_to_last_edge[a_in], 1) == edges(path_to_first_edge[a_out], 0));
242
243 // Replace a_out by a_in in last_vertex_in_path(a_out)
244 for (auto& ai : paths_in[last_vertex_in_path(a_out)]) {
245 if (ai == a_out) {
246 ai = a_in;
247 }
248 }
249
250 // Update chain to joint a_in --> a_out
251 next_edge_along_path[path_to_last_edge[a_in]] = path_to_first_edge[a_out];
252 path_to_last_edge[a_in] = path_to_last_edge[a_out];
253
254 // Cleaning up
255 path_to_first_edge[a_out] = invalid<Index>();
256 path_to_last_edge[a_out] = invalid<Index>();
257 path_is_removed[a_out] = true;
258 paths_in[v].clear();
259 paths_out[v].clear();
260 }
261 const Index v_first = first_vertex_in_path(a_in);
262 const Index v_last = last_vertex_in_path(a_in);
263 la_debug_assert(!path_is_removed[a_in]);
264 if (v_first == v_last && !path_is_pending[a_in]) {
265 path_is_pending[a_in] = true;
266 ears.push(a_in);
267 }
268 }
269 }
270
271 if (num_edges_removed != num_edges) {
272 logger().warn(
273 "Removing ears didn't result in an empty graph, number of edges remaining: {}",
274 num_edges - num_edges_removed);
275 Index true_num_removed = 0;
276 for (auto b : edge_is_removed) {
277 if (b) {
278 ++true_num_removed;
279 }
280 }
281 la_debug_assert(true_num_removed == num_edges_removed);
282 remaining_edges.resize(num_edges - num_edges_removed, 2);
283 for (Index e = 0, cnt = 0; e < num_edges; ++e) {
284 if (!edge_is_removed[e]) {
285 la_debug_assert(cnt < remaining_edges.rows());
286 remaining_edges.row(cnt++) = edges.row(e);
287 }
288 }
289 return false;
290 }
291
292 return true;
293}
294
295} // namespace legacy
296} // namespace lagrange
LA_CORE_API spdlog::logger & logger()
Retrieves the current logger.
Definition: Logger.cpp:40
#define la_debug_assert(...)
Debug assertion check.
Definition: assert.h:189
constexpr void swap(function_ref< R(Args...)> &lhs, function_ref< R(Args...)> &rhs) noexcept
Swaps the referred callables of lhs and rhs.
Definition: function_ref.h:114
Main namespace for Lagrange.
Definition: AABBIGL.h:30