diff --git a/Tetrahedral_remeshing/include/CGAL/Tetrahedral_remeshing/internal/flip_edges.h b/Tetrahedral_remeshing/include/CGAL/Tetrahedral_remeshing/internal/flip_edges.h index 3e9f1902764..ed33ce8a714 100644 --- a/Tetrahedral_remeshing/include/CGAL/Tetrahedral_remeshing/internal/flip_edges.h +++ b/Tetrahedral_remeshing/include/CGAL/Tetrahedral_remeshing/internal/flip_edges.h @@ -25,6 +25,7 @@ #include #include +#include #include #include @@ -88,7 +89,7 @@ void update_c3t3_facets(C3t3& c3t3, template Sliver_removal_result flip_3_to_2(typename C3t3::Edge& edge, C3t3& c3t3, - const std::vector& vertices_around_edge, + const boost::container::small_vector& vertices_around_edge, const Flip_Criterion& criterion, IncCellsVectorMap& inc_cells, CellSelector& cell_selector) @@ -578,23 +579,36 @@ void find_best_flip_to_improve_dh(C3t3& c3t3, Facet_circulator curr_fcirc = tr.incident_facets(edge); Facet_circulator curr_fdone = curr_fcirc; - //Only keep the possible flips - std::vector opposite_vertices; - int nb_cells_around_edge = 0; + //Collect the vertex opposite to the edge in each facet around it, in + //circulation order. The chord test below used to re-circulate the facets + //around the edge once per such vertex to enumerate the other apices; those + //are the same vertices collected here, so the tests can be indexed on the + //ring instead. is_edge_uv is read-only, so stopping at the first chord + //found yields the same verdict as running the ring to its end. + boost::container::small_vector ring_apices; do { - Vertex_handle vh; - //Get the ids of the opposite vertices + //Get the id of the opposite vertex for (int i = 0; i < 3; ++i) { Vertex_handle curr_vertex = curr_fcirc->first->vertex( indices(curr_fcirc->second, i)); if (curr_vertex != vh0 && curr_vertex != vh1) { - vh = curr_vertex; + ring_apices.push_back(curr_vertex); break; } } + } + while (++curr_fcirc != curr_fdone); + + //Only keep the possible flips + std::vector opposite_vertices; + int nb_cells_around_edge = 0; + const int n_apices = static_cast(ring_apices.size()); + for (int p = 0; p < n_apices; ++p) + { + const Vertex_handle vh = ring_apices[p]; if(tr.is_infinite(vh)) continue; @@ -603,37 +617,23 @@ void find_best_flip_to_improve_dh(C3t3& c3t3, if (o_inc_vh.empty()) tr.incident_cells(vh, std::back_inserter(o_inc_vh)); - Facet_circulator facet_circulator = curr_fcirc; - Facet_circulator facet_done = curr_fcirc; - - facet_done--; - facet_circulator++; - facet_circulator++; + //a chord is an edge joining vh to an apex that is not one of its two + //neighbors on the ring (positions p-1 and p+1) bool is_edge = false; - do + for (int j = p + 2; j <= p + n_apices - 2; ++j) { - //Get the ids of the opposite vertices - for (int i = 0; i < 3; ++i) + if (is_edge_uv(vh, ring_apices[j % n_apices], o_inc_vh)) { - Vertex_handle curr_vertex = facet_circulator->first->vertex( - indices(facet_circulator->second, i)); - if (curr_vertex != vh0 && curr_vertex != vh1) - { - if (is_edge_uv(vh, curr_vertex, o_inc_vh)) - { - is_edge = true; - break; - } - } + is_edge = true; + break; } - } while (++facet_circulator != facet_done); + } if (!is_edge) opposite_vertices.push_back(vh); nb_cells_around_edge++; } - while (++curr_fcirc != curr_fdone); if (nb_cells_around_edge < 4) return; @@ -1115,8 +1115,12 @@ Sliver_removal_result find_best_flip(typename C3t3::Edge& edge, Facet_circulator circ = tr.incident_facets(edge); Facet_circulator done = circ; - //Identify the vertices around this edge - std::unordered_set vertices_around_edge; + //Identify the vertices around this edge. The ring of apices around an edge + //holds a handful of distinct vertices, so they are kept inline and scanned + //rather than hashed. flip_3_to_2 does not depend on their order : it picks + //vh2/vh3 by testing each vertex against ch0/ch1 individually, and is_facet + //is symmetric in its three vertices. + boost::container::small_vector vertices_around_edge; bool boundary_edge = false; bool hull_edge = false; @@ -1130,7 +1134,9 @@ Sliver_removal_result find_best_flip(typename C3t3::Edge& edge, Vertex_handle vi = circ->first->vertex(indices(circ->second, i)); if (vi != v0 && vi != v1) { - vertices_around_edge.insert(vi); + if (std::find(vertices_around_edge.begin(), vertices_around_edge.end(), vi) + == vertices_around_edge.end()) + vertices_around_edge.push_back(vi); if ( circ->first->subdomain_index() != circ->first->neighbor(circ->second)->subdomain_index()) @@ -1161,9 +1167,7 @@ Sliver_removal_result find_best_flip(typename C3t3::Edge& edge, { if (!boundary_edge && !hull_edge) { - std::vector vertices; - vertices.insert(vertices.end(), vertices_around_edge.begin(), vertices_around_edge.end()); - res = flip_3_to_2(edge, c3t3, vertices, criterion, inc_cells, cell_selector); + res = flip_3_to_2(edge, c3t3, vertices_around_edge, criterion, inc_cells, cell_selector); } } else diff --git a/Tetrahedral_remeshing/include/CGAL/Tetrahedral_remeshing/internal/tetrahedral_remeshing_helpers.h b/Tetrahedral_remeshing/include/CGAL/Tetrahedral_remeshing/internal/tetrahedral_remeshing_helpers.h index 15f8a868ca3..56074da2274 100644 --- a/Tetrahedral_remeshing/include/CGAL/Tetrahedral_remeshing/internal/tetrahedral_remeshing_helpers.h +++ b/Tetrahedral_remeshing/include/CGAL/Tetrahedral_remeshing/internal/tetrahedral_remeshing_helpers.h @@ -976,6 +976,63 @@ std::size_t nb_incident_subdomains(const typename C3t3::Vertex_handle v, return indices.size(); } +// `nb_incident_subdomains(v, c3t3) > 1`, without counting the whole star. The +// traversal is the one `TDS_3::incident_cells_3()` performs - from v's cell, +// across the facets that contain v, marking cells as it goes - so it sees the +// same cells in the same order, and stops as soon as a second index appears. +// The marks are the TDS's own conflict flags and are cleared before +// returning, so this must not be called from inside another marking +// traversal; topology_test, which calls it, is not. +template +bool has_several_incident_subdomains(const typename C3t3::Vertex_handle v, + const C3t3& c3t3) +{ + typedef typename C3t3::Triangulation::Cell_handle Cell_handle; + CGAL_USE(c3t3); + + const Cell_handle start = v->cell(); + const auto si0 = start->subdomain_index(); + + boost::container::small_vector marked; + boost::container::small_vector to_visit; + + start->tds_data().mark_in_conflict(); + marked.push_back(start); + to_visit.push_back(start); + + bool several = false; + while (!to_visit.empty()) + { + const Cell_handle c = to_visit.back(); + to_visit.pop_back(); + + if (c->subdomain_index() != si0) + { + several = true; + break; + } + + for (int i = 0; i < 4; ++i) + { + if (c->vertex(i) == v) + continue; + + const Cell_handle n = c->neighbor(i); + if (!n->tds_data().is_clear()) + continue; + + n->tds_data().mark_in_conflict(); + marked.push_back(n); + to_visit.push_back(n); + } + } + + for (const Cell_handle c : marked) + c->tds_data().clear(); + + return several; +} + template std::size_t nb_incident_subdomains(const typename C3t3::Edge& e, const C3t3& c3t3) @@ -1312,10 +1369,15 @@ bool topology_test(const typename C3t3::Edge& edge, for (int i = 1; i < 4; i++) { Vertex_handle vi = f.first->vertex((f.second + i) % 4); - if (vi != v0 && vi != v1 && nb_incident_subdomains(vi, c3t3) > 1) + if (vi != v0 && vi != v1) { - if (is_edge_in_complex(v0, vi, c3t3) - && is_edge_in_complex(v1, vi, c3t3)) + //(v0,vi) and (v1,vi) are edges of f, so testing them for the + //complex needs no is_edge() star walk, and feature edges are rare + //enough that the subdomain star walk is skipped almost always. + //The three tests are pure, so the conjunction is unchanged. + if (c3t3.is_in_complex(v0, vi) + && c3t3.is_in_complex(v1, vi) + && has_several_incident_subdomains(vi, c3t3)) return false; } }