reorganize some of the debugging code

This commit is contained in:
Laurent Rineau 2025-10-08 15:59:00 +02:00
parent 89393e1b7c
commit c26c013b5a
1 changed files with 68 additions and 80 deletions

View File

@ -1914,7 +1914,7 @@ private:
}; };
template <typename Fh_region> template <typename Fh_region>
int does_edge_intersect_region(Cell_handle cell, int index_vc, int index_vd, int does_edge_interior_intersect_region(Cell_handle cell, int index_vc, int index_vd,
const CDT_2& cdt_2, const Fh_region& fh_region) const CDT_2& cdt_2, const Fh_region& fh_region)
{ {
const auto vc = cell->vertex(index_vd); const auto vc = cell->vertex(index_vd);
@ -1938,7 +1938,7 @@ private:
if(opc == CGAL::COPLANAR || opd == CGAL::COPLANAR || opc == opd) { if(opc == CGAL::COPLANAR || opd == CGAL::COPLANAR || opc == opd) {
continue; continue;
} else { } else {
// otherwise the segment intersects the plane of the triangle // otherwise the segment interior intersects the plane of the triangle
if(CGAL::orientation(pc, pd, t0, t1) != opc && if(CGAL::orientation(pc, pd, t0, t1) != opc &&
CGAL::orientation(pc, pd, t1, t2) != opc && CGAL::orientation(pc, pd, t1, t2) != opc &&
CGAL::orientation(pc, pd, t2, t0) != opc) CGAL::orientation(pc, pd, t2, t0) != opc)
@ -1987,7 +1987,7 @@ private:
//write_segment(dump_edges_around, cell_circ->vertex(index_vc), cell_circ->vertex(index_vd)); //write_segment(dump_edges_around, cell_circ->vertex(index_vc), cell_circ->vertex(index_vd));
if(is_marked(cell_circ->vertex(index_vc), Vertex_marker::REGION_BORDER)) continue; if(is_marked(cell_circ->vertex(index_vc), Vertex_marker::REGION_BORDER)) continue;
if(is_marked(cell_circ->vertex(index_vd), Vertex_marker::REGION_BORDER)) continue; if(is_marked(cell_circ->vertex(index_vd), Vertex_marker::REGION_BORDER)) continue;
int cd_intersects_region = does_edge_intersect_region(cell_circ, index_vc, index_vd, cdt_2, fh_region); int cd_intersects_region = does_edge_interior_intersect_region(cell_circ, index_vc, index_vd, cdt_2, fh_region);
if(cd_intersects_region == 1) { if(cd_intersects_region == 1) {
return Search_first_intersection_result_type{ Edge{cell_circ, index_vc, index_vd}, border_edge }; return Search_first_intersection_result_type{ Edge{cell_circ, index_vc, index_vd}, border_edge };
} }
@ -2044,7 +2044,7 @@ private:
std::set<std::pair<Vertex_handle, Vertex_handle>> non_intersecting_edges_set; std::set<std::pair<Vertex_handle, Vertex_handle>> non_intersecting_edges_set;
if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_regions()) { if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_regions()) {
debug_dump_tetrahedra_intersect_region(face_index, region_index, cdt_2, fh_region); expensive_debug_dump_tetrahedra_intersect_region(face_index, region_index, cdt_2, fh_region);
} }
detect_edges_and_cells_intersecting_region(face_index, region_index, cdt_2, fh_region, region_border_vertices, detect_edges_and_cells_intersecting_region(face_index, region_index, cdt_2, fh_region, region_border_vertices,
@ -2085,7 +2085,9 @@ private:
} }
} }
// use edges of the border facets to unify sets if(this->debug_regions()) {
std::cerr << "...use edges of the border facets to unify sets\n";
}
for(auto facet: facets_of_border) { for(auto facet: facets_of_border) {
auto vertices = tr().vertices(facet); auto vertices = tr().vertices(facet);
for(int i = 0; i < 3; ++i) { for(int i = 0; i < 3; ++i) {
@ -2098,8 +2100,10 @@ private:
} }
} }
// use non-intersecting edges to unify sets, until we have at most 2 sets
if(vertices_of_cavity_union_find.number_of_sets() > 2) { if(vertices_of_cavity_union_find.number_of_sets() > 2) {
if(this->debug_regions()) {
std::cerr << "...use non-intersecting edges to unify sets, until we have at most 2 sets\n";
}
for(auto c : cr_intersecting_cells) { for(auto c : cr_intersecting_cells) {
for(int i = 0; i < 4; ++i) { for(int i = 0; i < 4; ++i) {
@ -2147,10 +2151,15 @@ private:
return std::make_pair(Vertex_handle{}, Facet{}); return std::make_pair(Vertex_handle{}, Facet{});
}); });
CGAL_assume(vertex_above != Vertex_handle{}); CGAL_assume(vertex_above != Vertex_handle{});
if(this->debug_regions()) {
std::cerr << "The vertex above the region is " << IO::oformat(vertex_above, with_point_and_info) << "\n";
}
// if there are still more than 2 sets, we need to propagate the information // if there are still more than 2 sets, we need to propagate the information
// using a DFS on the border facets // using a DFS on the border facets
if(vertices_of_cavity_union_find.number_of_sets() > 2) { if(vertices_of_cavity_union_find.number_of_sets() > 2) {
if(this->debug_regions()) {
std::cerr << "...propagate the information using a DFS on the border facets\n";
}
const auto border_edges_set = std::invoke([&] { const auto border_edges_set = std::invoke([&] {
using Ordered_pair = std::pair<Vertex_handle, Vertex_handle>; using Ordered_pair = std::pair<Vertex_handle, Vertex_handle>;
using Hash = boost::hash<Ordered_pair>; using Hash = boost::hash<Ordered_pair>;
@ -2527,8 +2536,12 @@ private:
// std::ofstream dump(std::string("dump_no_segment_found_") + std::to_string(face_index) + "_" + // std::ofstream dump(std::string("dump_no_segment_found_") + std::to_string(face_index) + "_" +
// std::to_string(region_index) + ".binary.cgal"); // std::to_string(region_index) + ".binary.cgal");
// CGAL::IO::save_binary_file(dump, *this); // CGAL::IO::save_binary_file(dump, *this);
dump_region(face_index, region_index, cdt_2);
// } // }
if(this->debug_regions()) {
std::cerr << "ERROR: No first segment found intersecting region " << region_index
<< " of face #" << face_index << "\n";
dump_region(face_index, region_index, cdt_2);
}
throw Next_region{"No segment found", fh_region[0]}; throw Next_region{"No segment found", fh_region[0]};
} }
CGAL_assertion(found_edge_opt != std::nullopt); CGAL_assertion(found_edge_opt != std::nullopt);
@ -2571,7 +2584,7 @@ private:
}; };
if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_regions()) { if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_regions()) {
std::cerr << cdt_3_format("Cavity has {} cells and {} edges, " std::cerr << cdt_3_format("Cavity has {} piercing cells and {} piercing edges, "
"{} vertices in upper cavity and {} in lower, " "{} vertices in upper cavity and {} in lower, "
"{} facets in upper cavity and {} in lower\n", "{} facets in upper cavity and {} in lower\n",
original_intersecting_cells.size(), original_intersecting_cells.size(),
@ -2580,20 +2593,6 @@ private:
original_vertices_of_lower_cavity.size(), original_vertices_of_lower_cavity.size(),
original_facets_of_upper_cavity.size(), original_facets_of_upper_cavity.size(),
original_facets_of_lower_cavity.size()); original_facets_of_lower_cavity.size());
if(original_intersecting_cells.size() > 3 || intersecting_edges.size() > 1) {
std::cerr << "!! Interesting case !!\n";
// dump_region(face_index, region_index, cdt_2);
// {
// std::ofstream out(std::string("dump_intersecting_edges_") + std::to_string(face_index) + "_" +
// std::to_string(region_index) + ".polylines.txt");
// out.precision(17);
// for(auto edge: intersecting_edges) {
// write_segment(out, edge);
// }
// }
// dump_facets_of_cavity(face_index, region_index, "lower", original_facets_of_lower_cavity);
// dump_facets_of_cavity(face_index, region_index, "upper", original_facets_of_upper_cavity);
}
} }
auto register_internal_constrained_facet = [this](Facet f) { this->register_facet_to_be_constrained(f); }; auto register_internal_constrained_facet = [this](Facet f) { this->register_facet_to_be_constrained(f); };
@ -3234,7 +3233,8 @@ private:
bool restore_face(CDT_3_signed_index face_index) { bool restore_face(CDT_3_signed_index face_index) {
CDT_2& non_const_cdt_2 = face_cdt_2[face_index]; CDT_2& non_const_cdt_2 = face_cdt_2[face_index];
const CDT_2& cdt_2 = non_const_cdt_2; const CDT_2& cdt_2 = non_const_cdt_2;
if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_copy_triangulation_into_hole()) { if constexpr (cdt_3_can_use_cxx20_format())
if(this->debug_copy_triangulation_into_hole() || this->debug_verbose_special_cases()) {
std::cerr << cdt_3_format("restore_face({}): CDT_2 has {} vertices\n", face_index, cdt_2.number_of_vertices()); std::cerr << cdt_3_format("restore_face({}): CDT_2 has {} vertices\n", face_index, cdt_2.number_of_vertices());
} }
for(const auto& edge : cdt_2.finite_edges()) { for(const auto& edge : cdt_2.finite_edges()) {
@ -3280,7 +3280,7 @@ private:
processed_faces.insert(fh_region.begin(), fh_region.end()); processed_faces.insert(fh_region.begin(), fh_region.end());
auto handle_error_with_region = [&](const char* what, CDT_2_face_handle fh_2d) { auto handle_error_with_region = [&](const char* what, CDT_2_face_handle fh_2d) {
if(this->debug_regions()) { if(this->debug_regions() || this->debug_verbose_special_cases()) {
std::cerr << "NOTE: " << what << " in sub-region " << (region_index - 1) std::cerr << "NOTE: " << what << " in sub-region " << (region_index - 1)
<< " of face #F" << face_index << '\n'; << " of face #F" << face_index << '\n';
} }
@ -3490,7 +3490,9 @@ public:
i = face_constraint_misses_subfaces_find_next(i); i = face_constraint_misses_subfaces_find_next(i);
if(i == npos) { if(i == npos) {
std::cerr << "ERROR: No more missing face to restore after a PLC error\n"; std::cerr << "ERROR: No more missing face to restore after a PLC error\n";
if(this->debug_regions() || this->debug_verbose_special_cases()) {
dump_region(e.face_index, e.region_index); dump_region(e.face_index, e.region_index);
}
throw; throw;
} }
std::cerr << "Next face is face #F " << i << '\n'; std::cerr << "Next face is face #F " << i << '\n';
@ -3625,9 +3627,6 @@ public:
return visited_edges.emplace(CGAL::make_sorted_pair(v0, v1), does_intersect); return visited_edges.emplace(CGAL::make_sorted_pair(v0, v1), does_intersect);
}; };
// Alias for C++20 extension warning workaround
auto& intersecting_edges_ = intersecting_edges;
intersecting_edges.push_back(first_intersecting_edge); intersecting_edges.push_back(first_intersecting_edge);
const auto [v0, v1] = tr().vertices(first_intersecting_edge); const auto [v0, v1] = tr().vertices(first_intersecting_edge);
(void)new_edge(v0, v1, true); (void)new_edge(v0, v1, true);
@ -3640,24 +3639,25 @@ public:
} }
auto test_edge = [&](Cell_handle cell, Vertex_handle v0, int index_v0, Vertex_handle v1, int index_v1, auto test_edge = [&](Cell_handle cell, Vertex_handle v0, int index_v0, Vertex_handle v1, int index_v1,
int expected) { int expected)
auto value_returned = [this](bool b) { {
auto value_returned = [this, v0, v1](bool b, bool not_visited) {
CGAL_USE(this); CGAL_USE(this);
if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_regions()) { if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_regions()) {
std::cerr << cdt_3_format(" return {}\n", b); std::cerr << cdt_3_format(" test_edge {} {} return {} {}\n",
IO::oformat(v0, with_point_and_info),
IO::oformat(v1, with_point_and_info),
b,
not_visited ? "(new)" : "(cached)");
} }
return b; return b;
}; };
if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_regions()) {
std::cerr << cdt_3_format(" test_edge {} {} ", IO::oformat(v0, with_point_and_info),
IO::oformat(v1, with_point_and_info));
}
auto [cached_value_it, not_visited] = new_edge(v0, v1, false); auto [cached_value_it, not_visited] = new_edge(v0, v1, false);
if(!not_visited) return value_returned(cached_value_it->second); if(!not_visited) return value_returned(cached_value_it->second, not_visited);
int v0v1_intersects_region = (is_marked(v0, Vertex_marker::REGION_INSIDE) || int v0v1_intersects_region =
is_marked(v1, Vertex_marker::REGION_INSIDE)) (is_marked(v0, Vertex_marker::REGION_INSIDE) || is_marked(v1, Vertex_marker::REGION_INSIDE))
? expected ? expected
: does_edge_intersect_region(cell, index_v0, index_v1, cdt_2, fh_region); : does_edge_interior_intersect_region(cell, index_v0, index_v1, cdt_2, fh_region);
if(v0v1_intersects_region != 0) { if(v0v1_intersects_region != 0) {
if(this->use_older_cavity_algorithm()) { if(this->use_older_cavity_algorithm()) {
if(v0v1_intersects_region != expected) { if(v0v1_intersects_region != expected) {
@ -3669,13 +3669,13 @@ public:
if(v0v1_intersects_region < 0) { if(v0v1_intersects_region < 0) {
std::swap(index_v0, index_v1); std::swap(index_v0, index_v1);
} }
intersecting_edges_.emplace_back(cell, index_v0, index_v1); intersecting_edges.emplace_back(cell, index_v0, index_v1);
cached_value_it->second = true; cached_value_it->second = true;
return value_returned(true); return value_returned(true, not_visited);
} else { } else {
non_intersecting_edges_set.insert(make_sorted_pair(v0, v1)); non_intersecting_edges_set.insert(make_sorted_pair(v0, v1));
cached_value_it->second = false; cached_value_it->second = false;
return value_returned(false); return value_returned(false, not_visited);
} }
}; };
@ -3741,10 +3741,10 @@ public:
{ {
intersecting_cells.insert(n_ch); intersecting_cells.insert(n_ch);
if(this->debug_regions()) { if(this->debug_regions()) {
std::cerr << "tetrahedron #" << n_ch->time_stamp() << " intersects the region\n"; std::cerr << "new tetrahedron #" << n_ch->time_stamp() << " intersects the region\n";
} }
} else if(this->debug_regions()) { } else if(this->debug_regions()) {
std::cerr << "NO, tetrahedron #" << n_ch->time_stamp() << " does not intersect the region\n"; std::cerr << "NO, new tetrahedron #" << n_ch->time_stamp() << " does not intersect the region\n";
} }
for(int i = 0; i < 4; ++i) { for(int i = 0; i < 4; ++i) {
for(int j = i + 1; j < 4; ++j) { for(int j = i + 1; j < 4; ++j) {
@ -3820,12 +3820,13 @@ public:
{ {
using EK = CGAL::Exact_predicates_exact_constructions_kernel; using EK = CGAL::Exact_predicates_exact_constructions_kernel;
const auto to_exact = CGAL::Cartesian_converter<Geom_traits, EK>(); const auto to_exact = CGAL::Cartesian_converter<Geom_traits, EK>();
const auto& cdt_2 = this->face_cdt_2[face_index];
std::cerr << cdt_3_format("restore_subface_region face index: {}, region #{}, intersecting edge #{}: ({} {})\n", std::cerr << cdt_3_format("restore_subface_region face index: {}, region #{}, intersecting edge #{}: ({} {})\n",
face_index, region_index, edge_index, face_index, region_index, edge_index,
IO::oformat(v_above, with_point_and_info), IO::oformat(v_above, with_point_and_info),
IO::oformat(v_below, with_point_and_info)); IO::oformat(v_below, with_point_and_info));
dump_region(face_index, region_index, this->face_cdt_2[face_index]); dump_region(face_index, region_index, cdt_2);
const auto p_above = this->point(v_above); const auto p_above = this->point(v_above);
const auto p_below = this->point(v_below); const auto p_below = this->point(v_below);
@ -3851,30 +3852,17 @@ public:
intersect_out.close(); intersect_out.close();
auto cells_around_intersecting_edge = Container_from_circulator{this->incident_cells(intersecting_edge)}; auto cells_around_intersecting_edge = Container_from_circulator{this->incident_cells(intersecting_edge)};
for(const auto& cell: cells_around_intersecting_edge) { for(const auto& ch: make_prevent_deref_range(cells_around_intersecting_edge)) {
CGAL_assertion(!cell.has_vertex(tr().infinite_vertex())); CGAL_assertion(!ch->has_vertex(tr().infinite_vertex()));
auto tetrahedron = auto tetrahedron = tr().tetrahedron(ch.current_circulator());
typename Geom_traits::Tetrahedron_3{tr().point(cell.vertex(0)), tr().point(cell.vertex(1)),
tr().point(cell.vertex(2)), tr().point(cell.vertex(3))};
for(auto fh: fh_region) {
auto v0 = fh->vertex(0)->info().vertex_handle_3d;
auto v1 = fh->vertex(1)->info().vertex_handle_3d;
auto v2 = fh->vertex(2)->info().vertex_handle_3d;
auto triangle = typename Geom_traits::Triangle_3{tr().point(v0), tr().point(v1), tr().point(v2)};
auto exact_triangle = to_exact(triangle);
}
std::cerr << cdt_3_format("Test tetrahedron (#{}):\n {}\n {}\n {}\n {}\n", std::cerr << cdt_3_format("Test tetrahedron (#{}):\n {}\n {}\n {}\n {}\n",
cell.time_stamp(), ch->time_stamp(),
IO::oformat(cell.vertex(0), with_point_and_info), IO::oformat(ch->vertex(0), with_point_and_info),
IO::oformat(cell.vertex(1), with_point_and_info), IO::oformat(ch->vertex(1), with_point_and_info),
IO::oformat(cell.vertex(2), with_point_and_info), IO::oformat(ch->vertex(2), with_point_and_info),
IO::oformat(cell.vertex(3), with_point_and_info)); IO::oformat(ch->vertex(3), with_point_and_info));
if(!std::any_of(fh_region.begin(), fh_region.end(), [&](const auto fh) { if(!std::any_of(fh_region.begin(), fh_region.end(), [&](const auto fh) {
auto v0 = fh->vertex(0)->info().vertex_handle_3d; auto triangle = cdt_2.triangle(fh);
auto v1 = fh->vertex(1)->info().vertex_handle_3d;
auto v2 = fh->vertex(2)->info().vertex_handle_3d;
auto triangle = typename Geom_traits::Triangle_3{tr().point(v0), tr().point(v1), tr().point(v2)};
bool b = does_tetrahedron_intersect_triangle_interior(tetrahedron, triangle, tr().geom_traits()); bool b = does_tetrahedron_intersect_triangle_interior(tetrahedron, triangle, tr().geom_traits());
if(b) { if(b) {
std::cerr << " intersects the region\n"; std::cerr << " intersects the region\n";
@ -3884,15 +3872,15 @@ public:
{ {
std::cerr << cdt_3_format( std::cerr << cdt_3_format(
"ERROR: The following tetrahedron (#{}) does not intersect the region:\n {}\n {}\n {}\n {}\n", "ERROR: The following tetrahedron (#{}) does not intersect the region:\n {}\n {}\n {}\n {}\n",
cell.time_stamp(), ch->time_stamp(),
IO::oformat(cell.vertex(0), with_point_and_info), IO::oformat(cell.vertex(1), with_point_and_info), IO::oformat(ch->vertex(0), with_point_and_info), IO::oformat(ch->vertex(1), with_point_and_info),
IO::oformat(cell.vertex(2), with_point_and_info), IO::oformat(cell.vertex(3), with_point_and_info)); IO::oformat(ch->vertex(2), with_point_and_info), IO::oformat(ch->vertex(3), with_point_and_info));
} }
} }
} }
template <typename Fh_region> template <typename Fh_region>
void debug_dump_tetrahedra_intersect_region(CDT_3_signed_index face_index, void expensive_debug_dump_tetrahedra_intersect_region(CDT_3_signed_index face_index,
int region_index, int region_index,
const CDT_2& cdt_2, const CDT_2& cdt_2,
const Fh_region& fh_region) const Fh_region& fh_region)
@ -3922,7 +3910,7 @@ public:
bool intersects = false; bool intersects = false;
for(int i = 0; i < 4; ++i) { for(int i = 0; i < 4; ++i) {
for(int j = i + 1; j < 4; ++j) { for(int j = i + 1; j < 4; ++j) {
int intersects_region = does_edge_intersect_region(ch, i, j, cdt_2, fh_region); int intersects_region = does_edge_interior_intersect_region(ch, i, j, cdt_2, fh_region);
if(intersects_region != 0) { if(intersects_region != 0) {
intersects = true; intersects = true;
} }