From edbc32959d2476e3925fa9f6d163634dd9b11c28 Mon Sep 17 00:00:00 2001 From: Laurent Rineau Date: Mon, 27 Oct 2025 17:03:34 +0100 Subject: [PATCH] refactor debug API so that cdt_3_from_off can use the official API --- .clang-format | 82 ++-- .../Conforming_Delaunay_triangulation_3.h | 324 ++++++-------- ...ing_constrained_Delaunay_triangulation_3.h | 423 ++++++++++++++---- .../cdt_3_from_off.cpp | 346 ++------------ 4 files changed, 593 insertions(+), 582 deletions(-) diff --git a/.clang-format b/.clang-format index 369aadbb15f..856d54b4fbc 100644 --- a/.clang-format +++ b/.clang-format @@ -1,35 +1,63 @@ --- +# CGAL clang-format configuration +# This file defines the code formatting style for C++ files in the CGAL project. + Language: Cpp BasedOnStyle: LLVM -AccessModifierOffset: -2 -AllowShortFunctionsOnASingleLine: true -BinPackParameters: false -BreakConstructorInitializers: BeforeComma + +# Indentation +AccessModifierOffset: -2 # Indent public:/private:/protected: 2 spaces to the left + +# Function formatting +AllowShortFunctionsOnASingleLine: Inline # Allow short inline/member functions on one line, but not free functions +AlwaysBreakAfterReturnType: None # Don't force return type on separate line + +# Parameter and argument formatting +BinPackParameters: false # Put each parameter on its own line for better readability + +# Constructor formatting +BreakConstructorInitializers: BeforeComma # Put comma before each initializer: `MyClass() \n , member1(val1) \n , member2(val2)` + +# Brace wrapping configuration BreakBeforeBraces: Custom BraceWrapping: - AfterCaseLabel: false - AfterClass: true - AfterControlStatement: MultiLine - AfterEnum: false - AfterFunction: false - AfterNamespace: false - AfterObjCDeclaration: false - AfterStruct: true - AfterUnion: false - AfterExternBlock: false - BeforeCatch: false - BeforeElse: false - BeforeLambdaBody: false - BeforeWhile: false - IndentBraces: false - SplitEmptyFunction: false - SplitEmptyRecord: false - SplitEmptyNamespace: false -ColumnLimit: 120 -# Force pointers to the type for C++. + AfterCaseLabel: false # Don't put brace on new line after case labels + AfterClass: true # Put opening brace on new line after class definition + AfterControlStatement: MultiLine # Only break before braces if the control statement spans multiple lines + AfterEnum: false # Don't break after enum + AfterFunction: false # Don't break after function declaration (keep on same line) + AfterNamespace: false # Don't break after namespace + AfterObjCDeclaration: false # Objective-C related (not used in CGAL) + AfterStruct: true # Put opening brace on new line after struct definition + AfterUnion: false # Don't break after union + AfterExternBlock: false # Don't break after extern "C" blocks + BeforeCatch: false # Don't put catch on new line + BeforeElse: false # Don't put else on new line + BeforeLambdaBody: false # Don't break before lambda body + BeforeWhile: false # Don't put while on new line (do-while loops) + IndentBraces: false # Don't indent the braces themselves + SplitEmptyFunction: false # Don't split empty functions across lines + SplitEmptyRecord: false # Don't split empty classes/structs across lines + SplitEmptyNamespace: false # Don't split empty namespaces across lines + +# Line length +ColumnLimit: 120 # Maximum line length of 120 characters + +# Pointer and reference alignment +# Force pointers and references to align with the type (e.g., `int* ptr` not `int *ptr`) DerivePointerAlignment: false PointerAlignment: Left -# Control the spaces around conditionals -SpacesInConditionalStatement: false -SpaceBeforeParens: false + +# Spacing in control statements +SpacesInConditionalStatement: false # No extra spaces inside conditionals: `if(condition)` not `if( condition )` +SpaceBeforeParens: false # No space before parentheses: `if(` not `if (` + +# Include directive handling +SortIncludes: Never # Preserve the original order of #include statements (don't sort them) + +# Preprocessor directive formatting +IndentPPDirectives: None # Don't indent preprocessor directives (#ifdef, #include, etc.) + +# Blank line handling +MaxEmptyLinesToKeep: 2 # Keep up to 2 consecutive blank lines ... diff --git a/Constrained_triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h b/Constrained_triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h index f16b04ce729..f047cff574d 100644 --- a/Constrained_triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h +++ b/Constrained_triangulation_3/include/CGAL/Conforming_Delaunay_triangulation_3.h @@ -41,7 +41,96 @@ namespace CGAL { -namespace CDT_3::internal { +namespace CDT_3 { + +struct Debug_options { + enum class Flags { + Steiner_points = 0, + conforming, + input_faces, + missing_region, + regions, + copy_triangulation_into_hole, + validity, + use_older_cavity_algorithm, + debug_finite_edges_map, + use_finite_edges_map, + debug_subconstraints_to_conform, + verbose_special_cases, + debug_encroaching_vertices, + debug_conforming_validation, + debug_constraint_hierarchy, + debug_geometric_errors, + debug_polygon_insertion, + display_statistics, + nb_of_flags + }; + + bool Steiner_points() const { return flags[static_cast(Flags::Steiner_points)]; } + void Steiner_points(bool b) { flags.set(static_cast(Flags::Steiner_points), b); } + + bool input_faces() const { return flags[static_cast(Flags::input_faces)]; } + void input_faces(bool b) { flags.set(static_cast(Flags::input_faces), b); } + + bool missing_region() const { return flags[static_cast(Flags::missing_region)]; } + void missing_region(bool b) { flags.set(static_cast(Flags::missing_region), b); } + + bool regions() const { return flags[static_cast(Flags::regions)]; } + void regions(bool b) { flags.set(static_cast(Flags::regions), b); } + + bool copy_triangulation_into_hole() const { return flags[static_cast(Flags::copy_triangulation_into_hole)]; } + void copy_triangulation_into_hole(bool b) { flags.set(static_cast(Flags::copy_triangulation_into_hole), b); } + + bool validity() const { return flags[static_cast(Flags::validity)]; } + void validity(bool b) { flags.set(static_cast(Flags::validity), b); } + + bool use_older_cavity_algorithm() const { return flags[static_cast(Flags::use_older_cavity_algorithm)]; } + bool use_newer_cavity_algorithm() const { return !flags[static_cast(Flags::use_older_cavity_algorithm)]; } + void use_older_cavity_algorithm(bool b) { flags.set(static_cast(Flags::use_older_cavity_algorithm), b); } + + bool finite_edges_map() const { return flags[static_cast(Flags::debug_finite_edges_map)]; } + void finite_edges_map(bool b) { flags.set(static_cast(Flags::debug_finite_edges_map), b); } + + bool subconstraints_to_conform() const { return flags[static_cast(Flags::debug_subconstraints_to_conform)]; } + void subconstraints_to_conform(bool b) { flags.set(static_cast(Flags::debug_subconstraints_to_conform), b); } + + bool use_finite_edges_map_flag() const { return flags[static_cast(Flags::use_finite_edges_map)]; } + void use_finite_edges_map(bool b) { flags.set(static_cast(Flags::use_finite_edges_map), b); } + + bool verbose_special_cases() const { return flags[static_cast(Flags::verbose_special_cases)]; } + void verbose_special_cases(bool b) { flags.set(static_cast(Flags::verbose_special_cases), b); } + + bool encroaching_vertices() const { return flags[static_cast(Flags::debug_encroaching_vertices)]; } + void encroaching_vertices(bool b) { flags.set(static_cast(Flags::debug_encroaching_vertices), b); } + + bool conforming_validation() const { return flags[static_cast(Flags::debug_conforming_validation)]; } + void conforming_validation(bool b) { flags.set(static_cast(Flags::debug_conforming_validation), b); } + + bool constraint_hierarchy() const { return flags[static_cast(Flags::debug_constraint_hierarchy)]; } + void constraint_hierarchy(bool b) { flags.set(static_cast(Flags::debug_constraint_hierarchy), b); } + + bool geometric_errors() const { return flags[static_cast(Flags::debug_geometric_errors)]; } + void geometric_errors(bool b) { flags.set(static_cast(Flags::debug_geometric_errors), b); } + + bool polygon_insertion() const { return flags[static_cast(Flags::debug_polygon_insertion)]; } + void polygon_insertion(bool b) { flags.set(static_cast(Flags::debug_polygon_insertion), b); } + + bool display_statistics() const { return flags[static_cast(Flags::display_statistics)]; } + void display_statistics(bool b) { flags.set(static_cast(Flags::display_statistics), b); } + + double segment_vertex_epsilon() const { return segment_vertex_epsilon_; } + void set_segment_vertex_epsilon(double eps) { segment_vertex_epsilon_ = eps; } + + double vertex_vertex_epsilon() const { return vertex_vertex_epsilon_; } + void set_vertex_vertex_epsilon(double eps) { vertex_vertex_epsilon_ = eps; } + +private: + std::bitset(Flags::nb_of_flags)> flags{}; + double segment_vertex_epsilon_ = 0.0; + double vertex_vertex_epsilon_ = 0.0; +}; // end struct Debug_options + +namespace internal { auto& tasks_manager() { struct Tasks_manager { @@ -80,6 +169,12 @@ auto& tasks_manager() { __itt_task_begin(instance->cdt_3_domain, __itt_null, __itt_null, instance->task_handles[task_id]); #endif } + auto time() const { + return instance->timers[task_id].time(); + } + auto time_ms() const { + return instance->timers[task_id].time() / 1000.; + } ~Scope_guard() { instance->timers[task_id].stop(); #if CGAL_USE_ITT @@ -107,7 +202,9 @@ auto& tasks_manager() { return instance; } // end auto& tasks_manager() -} // end namespace CDT_3::internal +} // end namespace internal + +} // end namespace CDT_3 inline auto CDT_3_READ_INPUT_TASK_guard() { return CDT_3::internal::tasks_manager().READ_INPUT_TASK_guard(); @@ -155,9 +252,9 @@ public: using Line = typename T_3::Geom_traits::Line_3; using Locate_type = typename T_3::Locate_type; - inline static With_offset_tag with_offset{}; - inline static With_point_tag with_point{}; - inline static With_point_and_info_tag with_point_and_info{}; + inline static With_offset_tag with_offset{-1}; + inline static With_point_tag with_point{-1}; + inline static With_point_and_info_tag with_point_and_info{-1}; Conforming_Delaunay_triangulation_3(const Geom_traits& gt = Geom_traits()) : T_3(gt) @@ -215,7 +312,7 @@ protected: if(v1 > v2) std::swap(v1, v2); auto v1_index = v1->time_stamp(); [[maybe_unused]] auto nb_erased = self->all_finite_edges[v1_index].erase(v2); - if constexpr (cdt_3_can_use_cxx20_format()) if(self->debug_finite_edges_map() && nb_erased > 0) { + if constexpr (cdt_3_can_use_cxx20_format()) if(self->debug().finite_edges_map() && nb_erased > 0) { std::cerr << cdt_3_format("erasing edge {} {}\n", self->display_vert((std::min)(v1, v2)), self->display_vert((std::max)(v1, v2))); } @@ -274,7 +371,7 @@ protected: void add_to_subconstraints_to_conform(Vertex_handle va, Vertex_handle vb, Constrained_polyline_id id) { const auto pair = make_subconstraint(va, vb); - if(debug_subconstraints_to_conform()) { + if(debug().subconstraints_to_conform()) { std::cerr << "tr().subconstraints_to_conform.push(" << display_subcstr(pair) << ")\n"; } @@ -321,7 +418,7 @@ protected: if(tr().is_infinite(v1) || tr().is_infinite(v2)) return; [[maybe_unused]] auto [_, inserted] = all_finite_edges[v1->time_stamp()].insert(v2); - if constexpr (cdt_3_can_use_cxx20_format()) if (debug_finite_edges_map() && inserted) { + if constexpr (cdt_3_can_use_cxx20_format()) if (debug().finite_edges_map() && inserted) { if(v2 < v1) std::swap(v1, v2); std::cerr << cdt_3_format("new_edge({}, {})\n", display_vert(v1), display_vert(v2)); } @@ -372,7 +469,7 @@ protected: if(use_finite_edges_map()) { new_vertex(v); all_finite_edges.clear(); - if (debug_finite_edges_map()) std::cerr << "all_finite_edges.clear()\n"; + if (debug().finite_edges_map()) std::cerr << "all_finite_edges.clear()\n"; for(auto e: tr().all_edges()) { new_edge(e); } @@ -408,7 +505,7 @@ protected: } } - void update_max_bbox_edge_length() { + void update_max_bbox_edge_length() const { double d_x = bbox.xmax() - bbox.xmin(); double d_y = bbox.ymax() - bbox.ymin(); double d_z = bbox.zmax() - bbox.zmin(); @@ -421,137 +518,15 @@ public: segment_vertex_epsilon = epsilon; } - bool debug_Steiner_points() const { - return debug_flags[static_cast(Debug_flags::Steiner_points)]; - } + CDT_3::Debug_options& debug() { return debug_options_; } + const CDT_3::Debug_options& debug() const { return debug_options_; } - void debug_Steiner_points(bool b) { - debug_flags.set(static_cast(Debug_flags::Steiner_points), b); - } - - bool debug_input_faces() const { - return debug_flags[static_cast(Debug_flags::input_faces)]; - } - - void debug_input_faces(bool b) { - debug_flags.set(static_cast(Debug_flags::input_faces), b); - } - - bool debug_missing_region() const { - return debug_flags[static_cast(Debug_flags::missing_region)]; - } - - void debug_missing_region(bool b) { - debug_flags.set(static_cast(Debug_flags::missing_region), b); - } - - bool debug_regions() const { - return debug_flags[static_cast(Debug_flags::regions)]; - } - - void debug_regions(bool b) { - debug_flags.set(static_cast(Debug_flags::regions), b); - } - - bool debug_copy_triangulation_into_hole() const { - return debug_flags[static_cast(Debug_flags::copy_triangulation_into_hole)]; - } - - void debug_copy_triangulation_into_hole(bool b) { - debug_flags.set(static_cast(Debug_flags::copy_triangulation_into_hole), b); - } - - bool debug_validity() const { - return debug_flags[static_cast(Debug_flags::validity)]; - } - - void debug_validity(bool b) { - debug_flags.set(static_cast(Debug_flags::validity), b); - } - - bool use_older_cavity_algorithm() const { - return debug_flags[static_cast(Debug_flags::use_older_cavity_algorithm)]; - } - - bool use_newer_cavity_algorithm() const { - return !debug_flags[static_cast(Debug_flags::use_older_cavity_algorithm)]; - } - - void use_older_cavity_algorithm(bool b) { - debug_flags.set(static_cast(Debug_flags::use_older_cavity_algorithm), b); - } - - bool debug_finite_edges_map() const { - return debug_flags[static_cast(Debug_flags::debug_finite_edges_map)]; - } - - void debug_finite_edges_map(bool b) { - debug_flags.set(static_cast(Debug_flags::debug_finite_edges_map), b); - } - - bool debug_subconstraints_to_conform() const { - return debug_flags[static_cast(Debug_flags::debug_subconstraints_to_conform)]; - } - - void debug_subconstraints_to_conform(bool b) { - debug_flags.set(static_cast(Debug_flags::debug_subconstraints_to_conform), b); - } - - bool use_finite_edges_map() const { - return update_all_finite_edges_ && debug_flags[static_cast(Debug_flags::use_finite_edges_map)]; - } - - void use_finite_edges_map(bool b) { - debug_flags.set(static_cast(Debug_flags::use_finite_edges_map), b); - } - - bool debug_verbose_special_cases() const { - return debug_flags[static_cast(Debug_flags::verbose_special_cases)]; - } - - void debug_verbose_special_cases(bool b) { - debug_flags.set(static_cast(Debug_flags::verbose_special_cases), b); - } - - bool debug_encroaching_vertices() const { - return debug_flags[static_cast(Debug_flags::debug_encroaching_vertices)]; - } - - void debug_encroaching_vertices(bool b) { - debug_flags.set(static_cast(Debug_flags::debug_encroaching_vertices), b); - } - - bool debug_conforming_validation() const { - return debug_flags[static_cast(Debug_flags::debug_conforming_validation)]; - } - - void debug_conforming_validation(bool b) { - debug_flags.set(static_cast(Debug_flags::debug_conforming_validation), b); - } - - bool debug_constraint_hierarchy() const { - return debug_flags[static_cast(Debug_flags::debug_constraint_hierarchy)]; - } - - void debug_constraint_hierarchy(bool b) { - debug_flags.set(static_cast(Debug_flags::debug_constraint_hierarchy), b); - } - - bool debug_geometric_errors() const { - return debug_flags[static_cast(Debug_flags::debug_geometric_errors)]; - } - - void debug_geometric_errors(bool b) { - debug_flags.set(static_cast(Debug_flags::debug_geometric_errors), b); - } - - bool debug_polygon_insertion() const { - return debug_flags[static_cast(Debug_flags::debug_polygon_insertion)]; - } - - void debug_polygon_insertion(bool b) { - debug_flags.set(static_cast(Debug_flags::debug_polygon_insertion), b); - } + // Backward compatibility wrappers (deprecated, use debug().method() instead) + bool use_older_cavity_algorithm() const { return debug_options_.use_older_cavity_algorithm(); } + bool use_newer_cavity_algorithm() const { return debug_options_.use_newer_cavity_algorithm(); } + void use_older_cavity_algorithm(bool b) { debug_options_.use_older_cavity_algorithm(b); } + bool use_finite_edges_map() const { return update_all_finite_edges_ && debug_options_.use_finite_edges_map_flag(); } + void use_finite_edges_map(bool b) { debug_options_.use_finite_edges_map(b); } Vertex_handle insert(const Point &p, Locate_type lt, Cell_handle c, int li, int lj) @@ -579,14 +554,14 @@ public: bool is_edge(Vertex_handle va, Vertex_handle vb) const { const bool is_edge_v1 = - ((debug_finite_edges_map() && use_finite_edges_map()) || !use_finite_edges_map()) && tr().tds().is_edge(va, vb); + ((debug().finite_edges_map() && use_finite_edges_map()) || !use_finite_edges_map()) && tr().tds().is_edge(va, vb); if(use_finite_edges_map() && va > vb) std::swap(va, vb); const auto va_index = va->time_stamp(); const bool is_edge_v2 = use_finite_edges_map() && all_finite_edges[va_index].find(vb) != all_finite_edges[va_index].end(); - if(debug_finite_edges_map() && use_finite_edges_map() && is_edge_v1 != is_edge_v2) { + if(debug().finite_edges_map() && use_finite_edges_map() && is_edge_v1 != is_edge_v2) { std::cerr << "!! Inconsistent edge status\n"; std::cerr << " -> constraint " << display_vert(va) << " " << display_vert(vb) << '\n'; std::cerr << " -> edge " << (is_edge_v1 ? "is" : "is not") << " in the triangulation\n"; @@ -606,7 +581,7 @@ public: [this](const auto &sc) { const auto [va, vb] = sc; const auto is_edge = this->is_edge(va, vb); - if constexpr (cdt_3_can_use_cxx20_format()) if(debug_conforming_validation()) { + if constexpr (cdt_3_can_use_cxx20_format()) if(debug().conforming_validation()) { std::cerr << cdt_3_format("is_conforming>> Edge is 3D: {} ({} , {})\n", is_edge, CGAL::IO::oformat(va, with_point_and_info), @@ -622,7 +597,7 @@ public: Vertex_handle vb, Vertex_handle min_vertex, double min_dist, - Check_distance distance_type) + Check_distance distance_type = Check_distance::NON_SQUARED_DISTANCE) const { if(!max_bbox_edge_length) { update_max_bbox_edge_length(); @@ -643,6 +618,26 @@ public: } } + void check_vertex_vertex_distance_or_throw(Vertex_handle va, + Vertex_handle vb, + double min_dist) const + { + if(!max_bbox_edge_length) { + update_max_bbox_edge_length(); + } + if(min_dist < debug_options_.vertex_vertex_epsilon() * *max_bbox_edge_length) + { + std::stringstream ss; + ss.precision(std::cerr.precision()); + ss << "Two vertices are too close to each other.\n"; + ss << " -> vertex " << display_vert(va) << '\n'; + ss << " -> vertex " << display_vert(vb) << '\n'; + ss << " -> distance = " << min_dist << '\n'; + ss << " -> max_bbox_edge_length = " << *max_bbox_edge_length << '\n'; + CGAL_error_msg(ss.str().c_str()); + } + } + auto ancestors_of_Steiner_vertex_on_edge(Vertex_handle v) const { std::pair result; CGAL_precondition(v->ccdt_3_data().is_Steiner_vertex_on_edge()); @@ -748,7 +743,7 @@ protected: if(!constraint_hierarchy.is_subconstraint(va, vb)) { continue; } - if(debug_subconstraints_to_conform()) { + if(debug().subconstraints_to_conform()) { std::cerr << "tr().subconstraints_to_conform.pop()=" << display_subcstr(subconstraint) << "\n"; } @@ -794,7 +789,7 @@ protected: const auto& [steiner_pt, hint, ref_vertex] = construct_Steiner_point(constraint, subconstraint); [[maybe_unused]] const auto v = insert_Steiner_point_on_subconstraint(steiner_pt, hint, subconstraint, constraint, visitor); - if(debug_Steiner_points()) { + if(debug().Steiner_points()) { const auto [c_start, c_end] = constraint_extremities(constraint); std::cerr << "(" << IO::oformat(va, with_offset) << ", " << IO::oformat(vb, with_offset) << ")"; std::cerr << ": [ " << display_vert(c_start) << " - " << display_vert(c_end) << " ] "; @@ -830,7 +825,7 @@ protected: this->constraint_hierarchy.constraints_end(), c_id) != this->constraint_hierarchy.constraints_end()); CGAL_assertion(this->constraint_hierarchy.vertices_in_constraint_begin(c_id) != this->constraint_hierarchy.vertices_in_constraint_end(c_id)); - if(debug_constraint_hierarchy()) { + if(debug().constraint_hierarchy()) { std::cerr << "constraint " << static_cast(c_id.vl_ptr()) << " has " << c_id.vl_ptr()->skip_size() << " vertices\n"; } @@ -891,7 +886,7 @@ protected: encroaching_vertices.insert(v); }; auto fill_encroaching_vertices = [&](const auto simplex) { - if(debug_encroaching_vertices()) { + if(debug().encroaching_vertices()) { std::cerr << " - " << IO::oformat(simplex, With_point_tag{}) << '\n'; } auto visit_cell = [&](Cell_handle cell) { @@ -937,7 +932,7 @@ protected: std::cerr << "!! The constraint passes through a vertex!\n"; std::cerr << " -> constraint " << display_vert(va) << " " << display_vert(vb) << '\n'; std::cerr << " -> vertex " << display_vert(v) << '\n'; - if(debug_geometric_errors()) { + if(debug().geometric_errors()) { debug_dump("bug-through-vertex"); } CGAL_error(); @@ -949,7 +944,7 @@ protected: std::for_each(tr().segment_traverser_simplices_begin(va, vb), tr().segment_traverser_simplices_end(), fill_encroaching_vertices); auto vector_of_encroaching_vertices = encroaching_vertices.extract_sequence(); - if(debug_encroaching_vertices()) { + if(debug().encroaching_vertices()) { std::cerr << " -> vector_of_encroaching_vertices (before filter):\n"; std::for_each(vector_of_encroaching_vertices.begin(), vector_of_encroaching_vertices.end(), @@ -965,7 +960,7 @@ protected: this->tr().point(v), pb) == ACUTE; }); - if(debug_encroaching_vertices()) { + if(debug().encroaching_vertices()) { std::cerr << " -> vector_of_encroaching_vertices (after filter):\n"; std::for_each(vector_of_encroaching_vertices.begin(), end, [&](Vertex_handle v) { std::cerr << " " << this->display_vert(v) << " angle " << approximate_angle(pa, this->tr().point(v), pb) @@ -1001,7 +996,7 @@ protected: return {midpoint_functor(pa, pb), va->cell(), va}; } - if(debug_encroaching_vertices()) { + if(debug().encroaching_vertices()) { std::cerr << "construct_Steiner_point( " << display_vert(va) << " , " << display_vert(vb) << " )\n"; } @@ -1102,7 +1097,7 @@ protected: static_assert(CGAL::cdt_3_msvc_2019_or_older() || CGAL::is_nothrow_movable_v); Bbox_3 bbox{}; double segment_vertex_epsilon = 1e-8; - std::optional max_bbox_edge_length; + mutable std::optional max_bbox_edge_length; using Pair_of_vertex_handles = std::pair; boost::container::map pair_of_vertices_to_cid; Insert_in_conflict_visitor insert_in_conflict_visitor = {this}; @@ -1134,27 +1129,7 @@ protected: } } - enum class Debug_flags { - Steiner_points = 0, - conforming, - input_faces, - missing_region, - regions, - copy_triangulation_into_hole, - validity, - use_older_cavity_algorithm, - debug_finite_edges_map, - use_finite_edges_map, - debug_subconstraints_to_conform, - verbose_special_cases, - debug_encroaching_vertices, - debug_conforming_validation, - debug_constraint_hierarchy, - debug_geometric_errors, - debug_polygon_insertion, - nb_of_flags - }; - std::bitset(Debug_flags::nb_of_flags)> debug_flags{}; + CDT_3::Debug_options debug_options_{}; bool is_Delaunay = true; }; @@ -1162,6 +1137,5 @@ protected: #endif // not DOXYGEN_RUNNING -# #endif // CGAL_CONFORMING_DELAUNAY_TRIANGULATION_3_H diff --git a/Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h b/Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h index 26e6444e0e5..fea4d208113 100644 --- a/Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h +++ b/Constrained_triangulation_3/include/CGAL/Conforming_constrained_Delaunay_triangulation_3.h @@ -9,8 +9,8 @@ // // Author(s) : Laurent Rineau -#ifndef CGAL_CONSTRAINED_DELAUNAY_TRIANGULATION_3_H -#define CGAL_CONSTRAINED_DELAUNAY_TRIANGULATION_3_H +#ifndef CGAL_CONFORMING_CONSTRAINED_DELAUNAY_TRIANGULATION_3_H +#define CGAL_CONFORMING_CONSTRAINED_DELAUNAY_TRIANGULATION_3_H #include @@ -66,6 +66,7 @@ #include #include #include +#include #include #include #include @@ -79,6 +80,7 @@ #include #include #include +#include #include #include #include @@ -578,6 +580,9 @@ public: #endif // not DOXYGEN_RUNNING using size_type = typename Triangulation::size_type; + /// \cond SKIP_IN_MANUAL + CDT_3_impl& impl() { return cdt_impl; } + /// /endcond public: /** \name Constructors @{ @@ -647,6 +652,7 @@ public: auto mesh_vp_map = choose_parameter(get_parameter(np, internal_np::vertex_point), get_const_property_map(vertex_point, mesh)); + this->debug() = choose_parameter(get_parameter(np, internal_np::debug), CDT_3::Debug_options{}); using graph_traits = boost::graph_traits; using vertex_descriptor = typename graph_traits::vertex_descriptor; @@ -779,27 +785,57 @@ public: cdt_impl.insert_vertices_range(vertices(mesh), mesh_vp_map, tr_vertex_pmap, p::vertex_is_constrained_map(v_selected_map)); - for(auto&& edges : patch_edges) { - auto index = &edges - &patch_edges[0]; - auto face_index_opt = - cdt_impl.insert_constrained_face_defined_by_its_border_edges(edges, mesh_vp_map, tr_vertex_pmap); - for(auto vd : extra_isolated_mesh_vertices[index]) { - auto vh = get(tr_vertex_pmap, vd); - cdt_impl.mark_vertex_as_isolated_in_a_constrained_face(vh, face_index_opt); - } + cdt_impl.template validate_distances_and_report(mesh, patch_edges, tr_vertex_pmap); + + { + auto task_guard = CGAL::CDT_3_CONFORMING_TASK_guard(); + for(auto&& edges : patch_edges) { + auto index = &edges - &patch_edges[0]; + auto face_index_opt = + cdt_impl.insert_constrained_face_defined_by_its_border_edges(edges, mesh_vp_map, tr_vertex_pmap); + for(auto vd : extra_isolated_mesh_vertices[index]) { + auto vh = get(tr_vertex_pmap, vd); + cdt_impl.mark_vertex_as_isolated_in_a_constrained_face(vh, face_index_opt); + } + } } } else { cdt_impl.insert_vertices_range(vertices(mesh), mesh_vp_map, tr_vertex_pmap); - for(auto f : faces(mesh)) { - auto face_vertices = vertices_around_face(halfedge(f, mesh), mesh); - auto range_of_vertices = CGAL::make_transform_range_from_property_map(face_vertices, tr_vertex_pmap); - cdt_impl.insert_constrained_face(range_of_vertices, false); + cdt_impl.template validate_distances_and_report(mesh, patch_edges, tr_vertex_pmap); + { + auto task_guard = CGAL::CDT_3_CONFORMING_TASK_guard(); + + for(auto f : faces(mesh)) { + auto face_vertices = vertices_around_face(halfedge(f, mesh), mesh); + auto range_of_vertices = CGAL::make_transform_range_from_property_map(face_vertices, tr_vertex_pmap); + cdt_impl.insert_constrained_face(range_of_vertices, false); + } + } + } + + { + auto task_guard = CGAL::CDT_3_CONFORMING_TASK_guard(); + cdt_impl.restore_Delaunay(); + if(cdt_impl.debug().display_statistics()) { + std::cout << "[timings] restored Delaunay in " << task_guard.time_ms() << " ms\n"; + std::cout << "Number of vertices after Delaunay: " << cdt_impl.number_of_vertices() << "\n"; + } + } + + const auto& visitor = parameters::get_parameter(np, internal_np::visitor); + if constexpr (!std::is_same_v, internal_np::Param_not_found>) { + visitor(*this); + } + { + auto task_guard = CGAL::CDT_3_CDT_TASK_guard(); + cdt_impl.restore_constrained_Delaunay(); + if(cdt_impl.debug().display_statistics()) { + std::cout << "[timings] restored constrained Delaunay in " << task_guard.time_ms() << " ms\n"; + std::cout << "Number of vertices after CDT: " << cdt_impl.number_of_vertices() << "\n\n"; } } - cdt_impl.restore_Delaunay(); - cdt_impl.restore_constrained_Delaunay(); // std::cerr << cdt_3_format("cdt_impl: {} vertices, {} cells\n", cdt_impl.number_of_vertices(), // cdt_impl.number_of_cells()); } @@ -1042,15 +1078,58 @@ public: /// \name Checking /// These methods are mainly a debugging help for the users of advanced features. /// @{ - - /*! - \brief returns whether the triangulation is valid. - When `verbose` is set to `true`, messages describing the first invalidity encountered are - printed. - */ - bool is_valid(bool verbose = false) const { return cdt_impl.is_valid(verbose); } - /// @} + + // ----------------------- + // Debug and advanced control + // ----------------------- +public: + /** + * \brief Access debug options. + * + * Returns a reference to the debug options object that controls various debugging + * and algorithmic behaviors. Use this to configure debug flags: + * + * \code + * ccdt.debug().Steiner_points(true); + * ccdt.debug().input_faces(true); + * \endcode + */ + CDT_3::Debug_options& debug() { return cdt_impl.debug(); } + const CDT_3::Debug_options& debug() const { return cdt_impl.debug(); } + + // Algorithmic switches and tolerances + void use_older_cavity_algorithm(bool b) { cdt_impl.use_older_cavity_algorithm(b); } + void use_finite_edges_map(bool b) { cdt_impl.use_finite_edges_map(b); } + void set_segment_vertex_epsilon(double eps) { cdt_impl.set_segment_vertex_epsilon(eps); } + + // Validation helpers + bool is_conforming() const { return cdt_impl.is_conforming(); } + bool is_valid(bool verbose = false) const { return cdt_impl.is_valid(verbose); } + + // IO helpers for debugging + bool write_missing_segments_file(std::ostream& out) { return cdt_impl.write_missing_segments_file(out); } + void write_all_segments_file(std::ostream& out) { cdt_impl.write_all_segments_file(out); } + + // Lightweight accessors used by tests/tools + size_type number_of_vertices() const { return static_cast(cdt_impl.number_of_vertices()); } + size_type number_of_cells() const { return triangulation().number_of_cells(); } + + // Ranges and utilities for post-processing (e.g., dumping Steiner vertices) + auto finite_vertex_handles() const { return cdt_impl.finite_vertex_handles(); } + auto finite_cell_handles() const { return triangulation().finite_cell_handles(); } + auto all_cell_handles() const { return triangulation().all_cell_handles(); } + auto finite_facets() const { return triangulation().finite_facets(); } + typename Triangulation::Cell_handle infinite_cell() const { return triangulation().infinite_cell(); } + + auto ancestors_of_Steiner_vertex_on_edge(Vertex_handle v) const + { return cdt_impl.ancestors_of_Steiner_vertex_on_edge(v); } + + // Write facets (for debugging/output) + template + void write_facets(std::ostream& out, const CDT& cdt, FacetRange&& facets) const { + cdt_impl.write_facets(out, cdt, std::forward(facets)); + } }; #ifndef DOXYGEN_RUNNING @@ -1423,7 +1502,7 @@ public: this->new_vertex(p_vh); - CGAL_assume(!this->debug_validity() || this->is_valid(true)); + CGAL_assume(!this->debug().validity() || this->is_valid(true)); return p_vh; } @@ -1471,6 +1550,8 @@ public: VertexHandleMap tr_vertex_pmap, const NamedParameters& np = parameters::default_values()) { + auto task_guard = CGAL::CDT_3_INSERT_VERTICES_TASK_guard(); + using parameters::get_parameter; using parameters::is_default_parameter; @@ -1495,6 +1576,13 @@ public: } add_bbox_points_if_not_dimension_3(); + + if(this->debug().display_statistics()) { + + std::cout << "[timings] inserted vertices in " + << task_guard.time_ms() << " ms\n"; + std::cout << "Number of vertices: " << this->number_of_vertices() << "\n\n"; + } } /** @@ -1663,7 +1751,7 @@ public: restore_Delaunay = false; - if(this->debug_input_faces()) { + if(this->debug().input_faces()) { std::cerr << "insert_constrained_face (" << std::size(vertex_handles) << " vertices):\n"; for(auto v: vertex_handles) { std::cerr << " - " << this->display_vert(v) << '\n'; @@ -1703,7 +1791,7 @@ public: } } while(circ != circ_end); - if(this->debug_input_faces()) { + if(this->debug().input_faces()) { std::stringstream filename; filename << "dump-input-face-" << polygon_contraint_id << "_polygon_" << borders.size() - 1 << ".polylines.txt"; std::ofstream os(filename.str()); @@ -1748,7 +1836,7 @@ public: face_constraint_misses_subfaces.resize(face_cdt_2.size()); extra_isolated_vertices.resize(face_cdt_2.size()); } - if(this->debug_input_faces()) { + if(this->debug().input_faces()) { std::cerr << "insert_constrained_face return the polygon_contraint_id: " << polygon_contraint_id << '\n'; } return polygon_contraint_id; @@ -1817,10 +1905,181 @@ public: extra_isolated_vertices[face_index].insert(vh); } + + struct Min_distance_result { + double min_distance; + std::array vertices_of_min_edge; + }; + + Min_distance_result compute_minimum_vertex_distance() const { + const auto& tr = *this; +#if CGAL_CXX20 && __cpp_lib_concepts >= 201806L && __cpp_lib_ranges >= 201911L + auto [min_sq_distance, min_edge] = + (std::ranges::min)(tr.finite_edges() | std::views::transform([&](auto edge) { + return std::make_pair(tr.segment(edge).squared_length(), edge); + })); +#else + auto transform_fct = [&](auto edge) { return std::make_pair(tr.segment(edge).squared_length(), edge); }; + auto min_p = transform_fct(*tr.finite_edges_begin()); + for (auto ite=tr.finite_edges_begin(); ite!=tr.finite_edges_end(); ++ite) + { + auto p = transform_fct(*ite); + if (p < min_p) + p = min_p; + } + auto [min_sq_distance, min_edge] = min_p; +#endif + auto min_distance = CGAL::to_double(CGAL::approximate_sqrt(min_sq_distance)); + auto vertices_of_min_edge = this->vertices(min_edge); + + return {min_distance, vertices_of_min_edge}; + } + + static void print_minimum_distance_info(const Min_distance_result& min_dist) { + std::cout << "Min distance between vertices: " << min_dist.min_distance << '\n' + << " between vertices: : " + << CGAL::IO::oformat(min_dist.vertices_of_min_edge[0], CGAL::With_point_tag{}) << " " + << CGAL::IO::oformat(min_dist.vertices_of_min_edge[1], CGAL::With_point_tag{}) << "\n\n"; + } + + struct Constraint_distance_result { + double min_distance; + Vertex_handle min_va, min_vb, min_vertex; + }; + + template + Constraint_distance_result + compute_constraint_vertex_distances_from_patches_borders(const BordersOfPatches& patch_edges, + const VertexPointMap& tr_vertex_pmap) const + { + +#if CGAL_CXX20 && __cpp_lib_concepts >= 201806L && __cpp_lib_ranges >= 201911L + auto edge_results = patch_edges + | std::views::join + | std::views::transform([&](const auto& edge_pair) { + auto [vda, vdb] = edge_pair; + auto va = get(tr_vertex_pmap, vda); + auto vb = get(tr_vertex_pmap, vdb); + auto [min_dist, min_v] = this->min_distance_and_vertex_between_constraint_and_encroaching_vertex(va, vb); + return std::make_tuple(CGAL::to_double(min_dist), va, vb, min_v); + }); + + auto min_result = std::ranges::min_element(edge_results, {}, [](const auto& tuple) { + return std::get<0>(tuple); + }); + + auto [min_distance, min_va, min_vb, min_vertex] = *min_result; +#else + double min_distance = (std::numeric_limits::max)(); + Vertex_handle min_va, min_vb, min_vertex; + + for(const auto& edges : patch_edges) { + for(auto [vda, vdb]: edges) { + auto va = get(tr_vertex_pmap, vda); + auto vb = get(tr_vertex_pmap, vdb); + auto [min_dist, min_v] = this->min_distance_and_vertex_between_constraint_and_encroaching_vertex(va, vb); + if(min_dist < min_distance) { + min_distance = CGAL::to_double(min_dist); + min_va = va; + min_vb = vb; + min_vertex = min_v; + } + } + } +#endif + + return {min_distance, min_va, min_vb, min_vertex}; + } + + template + Constraint_distance_result + compute_constraint_vertex_distances_from_faces(const Mesh& mesh, const VertexPointMap& tr_vertex_pmap) const { + double min_distance = (std::numeric_limits::max)(); + Vertex_handle min_va, min_vb, min_vertex; + + for(auto face_descriptor : faces(mesh)) { + auto he = halfedge(face_descriptor, mesh); + const auto end = he; + do { + auto va = get(tr_vertex_pmap, source(he, mesh)); + auto vb = get(tr_vertex_pmap, target(he, mesh)); + auto [min_dist, min_v] = this->min_distance_and_vertex_between_constraint_and_encroaching_vertex(va, vb); + if(min_dist < min_distance) { + min_distance = CGAL::to_double(min_dist); + min_va = va; + min_vb = vb; + min_vertex = min_v; + } + he = next(he, mesh); + } while((he = next(he, mesh)) != end); + } + + return {min_distance, min_va, min_vb, min_vertex}; + } + + template + void validate_constraint_vertex_distances_or_throw(const Mesh& mesh, + const BordersOfPatches& patch_edges, + const VertexPointMap& tr_vertex_pmap) + { + + auto [min_distance, min_va, min_vb, min_vertex] = + merged_facets ? compute_constraint_vertex_distances_from_patches_borders(patch_edges, tr_vertex_pmap) + : compute_constraint_vertex_distances_from_faces(mesh, tr_vertex_pmap); + + if(this->debug().display_statistics()) { + std::cout << "Min distance between constraint segment and vertex: " << min_distance << '\n' + << " between segment : " + << CGAL::IO::oformat(min_va, CGAL::With_point_tag{}) << " " + << CGAL::IO::oformat(min_vb, CGAL::With_point_tag{}) << '\n' + << " and vertex : " + << CGAL::IO::oformat(min_vertex, CGAL::With_point_tag{}) << "\n\n"; + } + this->check_segment_vertex_distance_or_throw(min_va, min_vb, min_vertex, min_distance); + } + + template + void validate_vertex_vertex_distances_or_throw([[maybe_unused]] const Mesh& mesh, + [[maybe_unused]] const BordersOfPatches& patch_edges, + [[maybe_unused]] const VertexPointMap& tr_vertex_pmap) + { + + auto result = this->compute_minimum_vertex_distance(); + + if(this->debug().display_statistics()) { + std::cout << "Min distance between vertices: " << result.min_distance << '\n' + << " between vertices : " + << CGAL::IO::oformat(result.vertices_of_min_edge[0], CGAL::With_point_tag{}) << " " + << CGAL::IO::oformat(result.vertices_of_min_edge[1], CGAL::With_point_tag{}) << "\n\n"; + } + this->check_vertex_vertex_distance_or_throw(result.vertices_of_min_edge[0], + result.vertices_of_min_edge[1], + result.min_distance); + } + + template + void validate_distances_and_report(const Mesh& mesh, + const BordersOfPatches& patch_edges, + const VertexPointMap& tr_vertex_pmap) + { + auto task_guard = CGAL::CDT_3_COMPUTE_DISTANCES_TASK_guard(); + if(this->debug().segment_vertex_epsilon() > 0) { + this->validate_constraint_vertex_distances_or_throw(mesh, patch_edges, tr_vertex_pmap); + } + if(this->debug().vertex_vertex_epsilon() > 0) { + this->validate_vertex_vertex_distances_or_throw(mesh, patch_edges, tr_vertex_pmap); + } + if(this->debug().display_statistics()) { + std::cout << "[timings] compute distances in " + << task_guard.time_ms() << " ms\n"; + } + } + + private: void set_mark(Vertex_handle v, Vertex_marker m) { - if(this->debug_regions()) { + if(this->debug().regions()) { std::cerr << " set_mark(" << this->display_vert(v) << ", "; switch(m) { case Vertex_marker::CLEAR: std::cerr << "CLEAR"; break; @@ -1900,7 +2159,7 @@ private: return vec_of_handles; }); - if(this->debug_input_faces()) { + if(this->debug().input_faces()) { std::cerr << "Polygon #" << polygon_contraint_id << " normal is: " << cdt_2.geom_traits().normal() << '\n'; auto filename = "dump_cdt_2_polygons_" + std::to_string(polygon_contraint_id) + ".polylines.txt"; std::cerr << " dumping it to \"" << filename << "\".\n"; @@ -1917,7 +2176,7 @@ private: // create and fill the 2D triangulation { auto insert_constraint_in_cdt_2 = [&](const auto& va, const auto& vb) { - if(this->debug_input_faces()) { + if(this->debug().input_faces()) { std::cerr << "cdt_2.insert_constraint (" << tr().point(va->info().vertex_handle_3d) << " , " @@ -2032,7 +2291,7 @@ private: } } } // end of marking inside/outside - if(this->debug_input_faces()) { + if(this->debug().input_faces()) { int counter = 0; for(const auto fh: cdt_2.finite_face_handles()) { if(!fh->info().is_outside_the_face) ++counter; @@ -2066,7 +2325,7 @@ private: if(this->is_edge(vb_3d, vd_3d)) { // let's insert the diagonal [bd] in the CDT_2 cdt_2.insert_constraint(vb, vd); - if(this->debug_verbose_special_cases()) { + if(this->debug().verbose_special_cases()) { std::cerr << "NOTE: CDT_2 has 4 vertices. Flip the diagonal\n"; } } @@ -2146,7 +2405,7 @@ private: border_edges.emplace_back(c, i, j); } } - if(this->debug_regions()) { + if(this->debug().regions()) { std::cerr << "region size is: " << fh_region.size() << "\n"; std::cerr << "region border size is: " << border_edges.size() << "\n"; } @@ -2292,7 +2551,7 @@ private: std::set> 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()) { expensive_debug_dump_tetrahedra_intersect_region(face_index, region_index, cdt_2, fh_region); } @@ -2334,7 +2593,7 @@ private: } } - if(this->debug_regions()) { + if(this->debug().regions()) { std::cerr << "...use edges of the border facets to unify sets\n"; } for(auto facet: facets_of_border) { @@ -2350,7 +2609,7 @@ private: } if(vertices_of_cavity_union_find.number_of_sets() > 2) { - if(this->debug_regions()) { + 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) { @@ -2400,13 +2659,13 @@ private: return std::make_pair(Vertex_handle{}, Facet{}); }); CGAL_assume(vertex_above != Vertex_handle{}); - if(this->debug_regions()) { + 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 // using a DFS on the border facets if(vertices_of_cavity_union_find.number_of_sets() > 2) { - if(this->debug_regions()) { + if(this->debug().regions()) { std::cerr << "...propagate the information using a DFS on the border facets\n"; } const auto border_edges_set = std::invoke([&] { @@ -2561,7 +2820,7 @@ private: // classify the facets of the border of the cavity for(auto facet: facets_of_border) { - if(this->debug_regions()) { + if(this->debug().regions()) { debug_output_facet_vertices({facet}); } for(auto v: tr().vertices(facet)) { @@ -2579,7 +2838,7 @@ private: clear_marks(vertices_of_lower_cavity, Vertex_marker::CAVITY_BELOW); } // new algorithm - if(this->debug_regions()) { + if(this->debug().regions()) { debug_dump_cavity_outputs(face_index, region_index, intersecting_edges, facets_of_border, facets_of_upper_cavity, facets_of_lower_cavity); for(auto edge : intersecting_edges) { auto [v1, v2] = tr().vertices(edge); @@ -2617,7 +2876,7 @@ private: void restore_subface_region(CDT_3_signed_index face_index, int region_index, CDT_2& non_const_cdt_2, Fh_region& non_const_fh_region) { - if(this->debug_regions()) { + if(this->debug().regions()) { std::cerr << "restore_subface_region face index: " << face_index << ", region #" << region_index << "\n"; } const auto& cdt_2 = non_const_cdt_2; @@ -2640,7 +2899,7 @@ private: } return vertices; }); - 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 << "region_border_vertices.size() = " << region_border_vertices.size() << "\n"; for(auto v : region_border_vertices) { std::cerr << cdt_3_format(" {}\n", IO::oformat(v, with_point)); @@ -2685,7 +2944,7 @@ private: if(cdt_2.orientation(v0->point(), v1->point(), v3->point()) == CGAL::POSITIVE && cdt_2.orientation(v0->point(), v3->point(), v2->point()) == CGAL::POSITIVE) { - if(this->debug_regions()) { + if(this->debug().regions()) { std::cerr << "NOTE: the other diagonal is in the 3D triangulation: flip the edge\n"; } non_const_cdt_2.flip(non_const_fh_region[0], diagonal_index); @@ -2705,7 +2964,7 @@ private: fh->info().missing_subface = false; } return true; - } else if(this->debug_regions()) { + } else if(this->debug().regions()) { std::cerr << "NOTE: the other diagonal is in the 3D triangulation BUT the edge is not flippable!\n"; std::cerr << " The region " << region_index << " of face #F" << face_index << " has four points:\n"; std::cerr << " v0: " << v0->point() << '\n'; @@ -2715,7 +2974,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 ("NOTE: diagonal: {:.6} {:.6} {} in tr\n", IO::oformat(*diagonal.begin(), with_point), @@ -2786,7 +3045,7 @@ private: // std::to_string(region_index) + ".binary.cgal"); // CGAL::IO::save_binary_file(dump, *this); // } - if(this->debug_regions()) { + 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); @@ -2832,7 +3091,7 @@ private: return true; }; - 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 {} piercing cells and {} piercing edges, " "{} vertices in upper cavity and {} in lower, " "{} facets in upper cavity and {} in lower\n", @@ -2845,7 +3104,7 @@ private: } auto register_internal_constrained_facet = [this](Facet f) { this->register_facet_to_be_constrained(f); }; - if(this->debug_copy_triangulation_into_hole()) { + if(this->debug().copy_triangulation_into_hole()) { std::cerr << "# upper cavity\n"; } [[maybe_unused]] const auto [upper_cavity_triangulation, vertices_of_upper_cavity, @@ -2855,7 +3114,7 @@ private: const auto& upper_cavity_triangulation_ = upper_cavity_triangulation; std::for_each(interior_constrained_faces_upper.begin(), interior_constrained_faces_upper.end(), register_internal_constrained_facet); - if(this->debug_copy_triangulation_into_hole()) { + if(this->debug().copy_triangulation_into_hole()) { std::cerr << "# lower cavity\n"; } [[maybe_unused]] const auto [lower_cavity_triangulation, vertices_of_lower_cavity, @@ -2888,7 +3147,7 @@ private: << " )"; return s.str(); }; - if(this->debug_missing_region()) { + if(this->debug().missing_region()) { if(fail_upper) { std::cerr << "NOTE: Face " << display_face() << " is not a facet of the upper cavity\n"; } @@ -2900,7 +3159,7 @@ private: } return false; })) { - if(this->debug_missing_region()) { + if(this->debug().missing_region()) { // debug_region_size_4(); dump_region(face_index, region_index, cdt_2); std::for_each(fh_region.begin(), fh_region.end(), [](auto fh) { fh->info().is_in_region = 3; }); @@ -2937,7 +3196,7 @@ private: insert_in_conflict_visitor.process_cells_in_conflict(cells_of_upper_cavity.begin(), cells_of_upper_cavity.end()); insert_in_conflict_visitor.process_cells_in_conflict(cells_of_lower_cavity.begin(), cells_of_lower_cavity.end()); - if(this->debug_copy_triangulation_into_hole()) { + if(this->debug().copy_triangulation_into_hole()) { std::cerr << "# glu the upper triangulation of the cavity\n"; if(cells_of_lower_cavity.size() > original_intersecting_cells.size() || cells_of_upper_cavity.size() > original_intersecting_cells.size()) @@ -2958,7 +3217,7 @@ private: outer_map[vt] = f; CGAL_USE(this); #if CGAL_CDT_3_CAN_USE_CXX20_FORMAT - if(this->debug_copy_triangulation_into_hole()) { + if(this->debug().copy_triangulation_into_hole()) { CGAL_assertion(vt[0] != vt[1]); CGAL_assertion(vt[0] != vt[2]); CGAL_assertion(vt[1] != vt[2]); @@ -2995,7 +3254,7 @@ private: if(!is_facet) continue; // we might be in a sliver in the plane of the polygon const auto [fh_2d, reverse_orientation] = *is_facet; - if(this->debug_regions()) facets_of_polygon.push_back(f); + if(this->debug().regions()) facets_of_polygon.push_back(f); const auto vt_aux = this->make_vertex_triple(f); typename T_3::Vertex_triple vt{map_cavity_vertices_to_ambient_vertices[vt_aux[0]], map_cavity_vertices_to_ambient_vertices[vt_aux[1]], @@ -3010,7 +3269,7 @@ private: CGAL_assertion(static_cast(facet_is_facet_of_cdt_2(*this, {new_cell, 3}, cdt_2))); add_to_outer_map(vt, {new_cell, 3}, "extra "); } - if(this->debug_regions()) { + if(this->debug().regions()) { std::ofstream out(cdt_3_format("dump_{}_pseudo_cells_region_{}_{}.off", is_upper_cavity ? "upper" : "lower", face_index, region_index)); out.precision(17); @@ -3025,7 +3284,7 @@ private: const auto upper_inner_map = tr().create_triangulation_inner_map( upper_cavity_triangulation, map_upper_cavity_vertices_to_ambient_vertices, false); - 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()) { std::cerr << "upper_inner_map:\n"; for(auto [vt, _] : upper_inner_map) { std::cerr << cdt_3_format(" {:.6}, {:.6}, {:.6})\n", @@ -3034,7 +3293,7 @@ private: IO::oformat(vt[2], with_point)); } } - if(this->debug_copy_triangulation_into_hole()) { + if(this->debug().copy_triangulation_into_hole()) { std::cerr << "# glu the lower triangulation of the cavity\n"; } this->copy_triangulation_into_hole(map_upper_cavity_vertices_to_ambient_vertices, @@ -3042,7 +3301,7 @@ private: upper_inner_map, this->new_cells_output_iterator()); } - if(this->debug_copy_triangulation_into_hole()) { + if(this->debug().copy_triangulation_into_hole()) { std::cerr << "# glu the lower triangulation of the cavity\n"; } @@ -3063,7 +3322,7 @@ private: const auto lower_inner_map = tr().create_triangulation_inner_map( lower_cavity_triangulation, map_lower_cavity_vertices_to_ambient_vertices, false); #if CGAL_CDT_3_CAN_USE_CXX20_FORMAT - if(this->debug_copy_triangulation_into_hole()) { + if(this->debug().copy_triangulation_into_hole()) { std::cerr << "outer_map:\n"; for(auto [vt, _] : outer_map) { std::cerr << cdt_3_format(" {:.6}, {:.6}, {:.6})\n", @@ -3104,7 +3363,7 @@ private: set_facet_constrained(f, face_index, f2d); f2d->info().missing_subface = false; } - CGAL_assume(!this->debug_validity() || this->is_valid(true)); + CGAL_assume(!this->debug().validity() || this->is_valid(true)); }; // ------------------------- @@ -3237,7 +3496,7 @@ private: tr().is_infinite(v) ? cavity_triangulation.infinite_vertex() : cavity_triangulation.insert(this->point(v)); map_ambient_vertices_to_cavity_vertices[v] = cavity_v; map_cavity_vertices_to_ambient_vertices[cavity_v] = v; - 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("inserted {}cavity vertex {:.6} -> {:.6}\n", extra, IO::oformat(cavity_v, with_point_and_info), @@ -3337,7 +3596,7 @@ private: { const auto& cdt_2 = non_const_cdt_2; auto steiner_pt = CGAL::centroid(cdt_2.triangle(fh_2d)); - if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_verbose_special_cases()) { + if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug().verbose_special_cases()) { std::cerr << cdt_3_format("Trying to insert Steiner (centroid) point {} in non-coplanar face {}.\n", IO::oformat(steiner_pt), IO::oformat(cdt_2.triangle(fh_2d))); } @@ -3345,7 +3604,7 @@ private: if(encroached_edge_opt) { return encroached_edge_opt; } - if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_Steiner_points()) { + if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug().Steiner_points()) { std::cerr << cdt_3_format("Inserting Steiner (centroid) point {} in non-coplanar face {}: {}.\n", IO::oformat(steiner_pt), face_index, IO::oformat(cdt_2.triangle(fh_2d))); } @@ -3415,7 +3674,7 @@ private: const auto v = this->insert_in_cdt_3(steiner_pt, lt, ch, li, lj, insert_in_conflict_visitor);// TODO: use "insert in hole" // this->study_bug = false; // assert(is_valid(true)); - if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_Steiner_points()) { + if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug().Steiner_points()) { std::cerr << " -> " << IO::oformat(v, with_offset) << '\n'; } v->ccdt_3_data().set_Steiner_vertex_in_face(face_index); @@ -3438,13 +3697,13 @@ private: const auto a = this->point(va_3d); const auto b = this->point(vb_3d); const auto mid = CGAL::midpoint(a, b); - if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_Steiner_points()) { + if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug().Steiner_points()) { std::cerr << cdt_3_format("Inserting Steiner (midpoint) point {} of constrained edge ({:.6} , {:.6})\n", IO::oformat(mid), IO::oformat(va_3d, with_point_and_info), IO::oformat(vb_3d, with_point_and_info)); } auto&& contexts = this->constraint_hierarchy.contexts(va_3d, vb_3d); - if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_verbose_special_cases()) { + if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug().verbose_special_cases()) { if(std::next(contexts.begin()) != contexts.end()) { std::cerr << "ERROR: Edge is constrained by more than one constraint\n"; for(const auto& c : contexts) { @@ -3472,7 +3731,7 @@ private: [[maybe_unused]] auto v = this->insert_Steiner_point_on_subconstraint(mid, mid_c, {va_3d, vb_3d}, constrained_polyline_id, insert_in_conflict_visitor); - if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_Steiner_points()) { + if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug().Steiner_points()) { std::cerr << " -> " << IO::oformat(v, with_offset) << '\n'; } // this->study_bug = false; @@ -3483,7 +3742,7 @@ private: CDT_2& non_const_cdt_2 = face_cdt_2[face_index]; const CDT_2& cdt_2 = non_const_cdt_2; if constexpr (cdt_3_can_use_cxx20_format()) - if(this->debug_copy_triangulation_into_hole() || this->debug_verbose_special_cases()) { + 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()); } for(const auto& edge : cdt_2.finite_edges()) { @@ -3493,7 +3752,7 @@ private: const auto vb_3d = fh->vertex(cdt_2.ccw(i))->info().vertex_handle_3d; const bool is_3d = this->is_edge(va_3d, vb_3d); if constexpr(cdt_3_can_use_cxx20_format()) { - if(this->debug_copy_triangulation_into_hole() || this->debug_conforming_validation()) { + if(this->debug().copy_triangulation_into_hole() || this->debug().conforming_validation()) { std::cerr << cdt_3_format("Edge is 3D: {:6} ({} , {})\n", is_3d, IO::oformat(va_3d, with_point_and_info), IO::oformat(vb_3d, with_point_and_info)); } @@ -3529,11 +3788,11 @@ private: processed_faces.insert(fh_region.begin(), fh_region.end()); auto handle_error_with_region = [&](const char* what, CDT_2_face_handle fh_2d) { - if(this->debug_regions() || this->debug_verbose_special_cases()) { + if(this->debug().regions() || this->debug().verbose_special_cases()) { std::cerr << "NOTE: " << what << " in sub-region " << (region_index - 1) << " of face #F" << face_index << '\n'; } - if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug_verbose_special_cases()) { + if constexpr (cdt_3_can_use_cxx20_format()) if(this->debug().verbose_special_cases()) { std::cerr << " constrained edges are:\n"; for(auto [c, index]: cdt_2.constrained_edges()) { const auto va = c->vertex(cdt_2.cw(index)); @@ -3713,7 +3972,7 @@ public: fill_cdt_2(cdt_2, i); search_for_missing_subfaces(i); } - if(this->debug_input_faces()) { + if(this->debug().input_faces()) { for(CDT_3_signed_index i = 0, end = static_cast (face_constraint_misses_subfaces.size()); i < end; ++i) { dump_face(i); } @@ -3727,7 +3986,7 @@ public: if(restore_face(static_cast (i))) { face_constraint_misses_subfaces_reset(i); } else { - if(this->debug_missing_region() || this->debug_Steiner_points()) { + if(this->debug().missing_region() || this->debug().Steiner_points()) { std::cerr << "restore_face(" << i << ") incomplete, back to conforming...\n"; } Conforming_Dt::restore_Delaunay(insert_in_conflict_visitor); @@ -3739,7 +3998,7 @@ public: i = face_constraint_misses_subfaces_find_next(i); if(i == npos) { std::cerr << "ERROR: No more missing face to restore after a PLC error\n"; - if(this->debug_regions() || this->debug_verbose_special_cases()) { + if(this->debug().regions() || this->debug().verbose_special_cases()) { dump_region(e.face_index, e.region_index); } throw; @@ -3883,7 +4142,7 @@ public: for(std::size_t i = 0; i < intersecting_edges.size(); ++i) { const auto intersecting_edge = intersecting_edges[i]; const auto [v_above, v_below] = tr().vertices(intersecting_edge); - 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_edge_region_intersection(face_index, region_index, fh_region, i, v_above, v_below, intersecting_edge); } @@ -3891,7 +4150,7 @@ public: int expected) { auto value_returned = [this, v0, v1](bool b, bool not_visited) { - 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(" test_edge {} {} return {} {}\n", IO::oformat(v0, with_point_and_info), IO::oformat(v1, with_point_and_info), @@ -3952,7 +4211,7 @@ public: !test_edge(cell, v_below, index_v_below, vc, index_vc, -1) && this->use_older_cavity_algorithm()) { - if(this->debug_regions()) { + if(this->debug().regions()) { dump_triangulation(); dump_region(face_index, region_index, cdt_2); std::ofstream out(std::string("dump_two_edges_") + std::to_string(face_index) + ".polylines.txt"); @@ -3966,7 +4225,7 @@ public: } while(++facet_circ != facet_circ_end); if(this->use_newer_cavity_algorithm() && i + 1 == intersecting_edges.size()) { for(auto ch: intersecting_cells) { - if(this->debug_regions()) { + if(this->debug().regions()) { std::cerr << "tetrahedron #" << ch->time_stamp() << " intersects the region\n"; } for(int i = 0; i < 4; ++i) { @@ -3989,10 +4248,10 @@ public: })) { intersecting_cells.insert(n_ch); - if(this->debug_regions()) { + if(this->debug().regions()) { 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, new tetrahedron #" << n_ch->time_stamp() << " does not intersect the region\n"; } for(int i = 0; i < 4; ++i) { @@ -4570,4 +4829,4 @@ auto get_remeshing_triangulation(Conforming_constrained_Delaunay_triangulation_3 } // end CGAL -#endif // CGAL_CONSTRAINED_DELAUNAY_TRIANGULATION_3_H +#endif // CGAL_CONFORMING_CONSTRAINED_DELAUNAY_TRIANGULATION_3_H diff --git a/Constrained_triangulation_3/test/Constrained_triangulation_3/cdt_3_from_off.cpp b/Constrained_triangulation_3/test/Constrained_triangulation_3/cdt_3_from_off.cpp index 5359a58e5ff..86f864937a4 100644 --- a/Constrained_triangulation_3/test/Constrained_triangulation_3/cdt_3_from_off.cpp +++ b/Constrained_triangulation_3/test/Constrained_triangulation_3/cdt_3_from_off.cpp @@ -5,6 +5,7 @@ // #define CGAL_CDT_2_DEBUG_INTERSECTIONS 1 #include #include +#include #include #include #include @@ -45,12 +46,8 @@ using K = CGAL::Exact_predicates_inexact_constructions_kernel; #endif // use Epick -struct Vb : public CGAL::Conforming_constrained_Delaunay_triangulation_vertex_base_3 {}; -struct Cb : public CGAL::Conforming_constrained_Delaunay_triangulation_cell_base_3 {}; -struct Tds: public CGAL::Triangulation_data_structure_3 {}; -using Base_triantulation = CGAL::Delaunay_triangulation_3; -using CDT = CGAL::Conforming_constrained_Delaunay_triangulation_3_impl; -using Point = Base_triantulation::Point; +using CDT = CGAL::Conforming_constrained_Delaunay_triangulation_3; +using Point = K::Point_3; using Point_3 = K::Point_3; using Mesh = CGAL::Surface_mesh; @@ -258,35 +255,29 @@ CDT_options::CDT_options(int argc, char* argv[]) { } } -#if NO_TRY_CATCH -# define CDT_3_try if (true) -# define CDT_3_catch(X) if (false) -# define CDT_3_throw_exception_again -#else -// Else proceed normally. -# define CDT_3_try try -# define CDT_3_catch(X) catch(X) -# define CDT_3_throw_exception_again throw -#endif +CGAL::CDT_3::Debug_options cdt_debug_options(const CDT_options& options) { + CGAL::CDT_3::Debug_options cdt_debug; + cdt_debug.Steiner_points(options.verbose_level > 0); + cdt_debug.input_faces(options.debug_input_faces); + cdt_debug.missing_region(options.verbose_level > 1 || options.debug_missing_regions); + cdt_debug.regions(options.debug_regions); + cdt_debug.validity(options.debug_validity); + cdt_debug.finite_edges_map(options.debug_finite_edges_map); + cdt_debug.subconstraints_to_conform(options.debug_subconstraints_to_conform); + cdt_debug.verbose_special_cases(options.debug_verbose_special_cases); + cdt_debug.encroaching_vertices(options.debug_encroaching_vertices); + cdt_debug.conforming_validation(options.debug_conforming_validation); + cdt_debug.constraint_hierarchy(options.debug_constraint_hierarchy); + cdt_debug.geometric_errors(options.debug_geometric_errors); + cdt_debug.polygon_insertion(options.debug_polygon_insertion); + cdt_debug.copy_triangulation_into_hole(options.debug_copy_triangulation_into_hole); + cdt_debug.use_older_cavity_algorithm(!options.use_new_cavity_algorithm); + cdt_debug.use_finite_edges_map(options.use_finite_edges_map); + cdt_debug.display_statistics(!options.quiet); + cdt_debug.set_segment_vertex_epsilon(options.segment_vertex_epsilon); + cdt_debug.set_vertex_vertex_epsilon(options.vertex_vertex_epsilon); -void configure_cdt_debug_options(CDT& cdt, const CDT_options& options) { - cdt.debug_Steiner_points(options.verbose_level > 0); - cdt.debug_input_faces(options.debug_input_faces); - cdt.debug_missing_region(options.verbose_level > 1 || options.debug_missing_regions); - cdt.debug_regions(options.debug_regions); - cdt.debug_validity(options.debug_validity); - cdt.debug_finite_edges_map(options.debug_finite_edges_map); - cdt.debug_subconstraints_to_conform(options.debug_subconstraints_to_conform); - cdt.debug_verbose_special_cases(options.debug_verbose_special_cases); - cdt.debug_encroaching_vertices(options.debug_encroaching_vertices); - cdt.debug_conforming_validation(options.debug_conforming_validation); - cdt.debug_constraint_hierarchy(options.debug_constraint_hierarchy); - cdt.debug_geometric_errors(options.debug_geometric_errors); - cdt.debug_polygon_insertion(options.debug_polygon_insertion); - cdt.debug_copy_triangulation_into_hole(options.debug_copy_triangulation_into_hole); - cdt.use_older_cavity_algorithm(!options.use_new_cavity_algorithm); - cdt.use_finite_edges_map(options.use_finite_edges_map); - cdt.set_segment_vertex_epsilon(options.segment_vertex_epsilon); + return cdt_debug; } auto compute_bounding_box(const Mesh& mesh, const CDT_options& options) { @@ -334,7 +325,7 @@ std::function create_output_finalizer(const CDT& cdt, const CDT_options& auto& tr = cdt; - std::unordered_map cells_map; + std::unordered_map cells_map; for(auto ch : tr.all_cell_handles()) { cells_map[ch] = 1; @@ -380,13 +371,13 @@ std::function create_output_finalizer(const CDT& cdt, const CDT_options& std::ofstream dump(options.output_filename); dump.precision(17); #if CGAL_CXX20 && __cpp_lib_concepts >= 201806L && __cpp_lib_ranges >= 201911L - cdt.write_facets(dump, cdt, std::views::filter(cdt.finite_facets(), [&](auto f) { + cdt.write_facets(dump, cdt.triangulation(), std::views::filter(cdt.finite_facets(), [&](auto f) { return cdt.is_facet_constrained(f); })); #else auto is_facet_constrained = [&](auto f) { return cdt.is_facet_constrained(f); }; auto it_end = cdt.finite_facets_end(); - cdt.write_facets(dump, cdt, + cdt.write_facets(dump, cdt.triangulation(), CGAL::make_range( boost::make_filter_iterator(is_facet_constrained,cdt.finite_facets_begin(), it_end), boost::make_filter_iterator(is_facet_constrained,it_end, it_end))); @@ -396,155 +387,6 @@ std::function create_output_finalizer(const CDT& cdt, const CDT_options& }; } -struct Min_distance_result { - double min_distance; - std::array vertices_of_min_edge; -}; - -Min_distance_result compute_minimum_vertex_distance(const CDT& cdt) { -#if CGAL_CXX20 && __cpp_lib_concepts >= 201806L && __cpp_lib_ranges >= 201911L - auto [min_sq_distance, min_edge] = - (std::ranges::min)(cdt.finite_edges() | std::views::transform([&](auto edge) { - return std::make_pair(cdt.segment(edge).squared_length(), edge); - })); -#else - auto transform_fct = [&](auto edge) { return std::make_pair(cdt.segment(edge).squared_length(), edge); }; - auto min_p = transform_fct(*cdt.finite_edges_begin()); - for (auto ite=cdt.finite_edges_begin(); ite!=cdt.finite_edges_end(); ++ite) - { - auto p = transform_fct(*ite); - if (p < min_p) - p = min_p; - } - auto [min_sq_distance, min_edge] = min_p; -#endif - auto min_distance = CGAL::to_double(CGAL::approximate_sqrt(min_sq_distance)); - auto vertices_of_min_edge = cdt.vertices(min_edge); - - return {min_distance, vertices_of_min_edge}; -} - -void print_minimum_distance_info(const Min_distance_result& min_dist) { - std::cout << "Min distance between vertices: " << min_dist.min_distance << '\n' - << " between vertices: : " - << CGAL::IO::oformat(min_dist.vertices_of_min_edge[0], CGAL::With_point_tag{}) << " " - << CGAL::IO::oformat(min_dist.vertices_of_min_edge[1], CGAL::With_point_tag{}) << "\n\n"; -} - -int validate_minimum_vertex_distances(const CDT& cdt, double vertex_vertex_min_distance, const CDT_options& options) { - auto result = compute_minimum_vertex_distance(cdt); - - if(!options.quiet) { - print_minimum_distance_info(result); - } - if(result.min_distance < vertex_vertex_min_distance) { - std::cerr << "ERROR: min distance between vertices is too small\n"; - return EXIT_FAILURE; - } - return EXIT_SUCCESS; -} - -struct Constraint_distance_result { - double min_distance; - CDT::Vertex_handle min_va, min_vb, min_vertex; -}; - -template -Constraint_distance_result compute_constraint_vertex_distances_from_patches_borders( - CDT& cdt, - const BordersOfPatches& patch_edges, - const VertexPointMap& tr_vertex_pmap) { - -#if CGAL_CXX20 && __cpp_lib_concepts >= 201806L && __cpp_lib_ranges >= 201911L - auto edge_results = patch_edges - | std::views::join - | std::views::transform([&](const auto& edge_pair) { - auto [vda, vdb] = edge_pair; - auto va = get(tr_vertex_pmap, vda); - auto vb = get(tr_vertex_pmap, vdb); - auto [min_dist, min_v] = cdt.min_distance_and_vertex_between_constraint_and_encroaching_vertex(va, vb); - return std::make_tuple(CGAL::to_double(min_dist), va, vb, min_v); - }); - - auto min_result = std::ranges::min_element(edge_results, {}, [](const auto& tuple) { - return std::get<0>(tuple); - }); - - auto [min_distance, min_va, min_vb, min_vertex] = *min_result; -#else - double min_distance = (std::numeric_limits::max)(); - CDT::Vertex_handle min_va, min_vb, min_vertex; - - for(const auto& edges : patch_edges) { - for(auto [vda, vdb]: edges) { - auto va = get(tr_vertex_pmap, vda); - auto vb = get(tr_vertex_pmap, vdb); - auto [min_dist, min_v] = cdt.min_distance_and_vertex_between_constraint_and_encroaching_vertex(va, vb); - if(min_dist < min_distance) { - min_distance = CGAL::to_double(min_dist); - min_va = va; - min_vb = vb; - min_vertex = min_v; - } - } - } -#endif - - return {min_distance, min_va, min_vb, min_vertex}; -} - -template -Constraint_distance_result compute_constraint_vertex_distances_from_faces( - CDT& cdt, - const Mesh& mesh, - const VertexPointMap& tr_vertex_pmap) { - - double min_distance = (std::numeric_limits::max)(); - CDT::Vertex_handle min_va, min_vb, min_vertex; - - for(auto face_descriptor : faces(mesh)) { - auto he = halfedge(face_descriptor, mesh); - const auto end = he; - do { - auto va = get(tr_vertex_pmap, source(he, mesh)); - auto vb = get(tr_vertex_pmap, target(he, mesh)); - auto [min_dist, min_v] = cdt.min_distance_and_vertex_between_constraint_and_encroaching_vertex(va, vb); - if(min_dist < min_distance) { - min_distance = CGAL::to_double(min_dist); - min_va = va; - min_vb = vb; - min_vertex = min_v; - } - he = next(he, mesh); - } while((he = next(he, mesh)) != end); - } - - return {min_distance, min_va, min_vb, min_vertex}; -} - -template -void validate_constraint_vertex_distances_or_throw( - CDT& cdt, - const Mesh& mesh, - const CDT_options& options, - const BordersOfPatches& patch_edges, - const VertexPointMap& tr_vertex_pmap) { - - auto [min_distance, min_va, min_vb, min_vertex] = - options.merge_facets ? compute_constraint_vertex_distances_from_patches_borders(cdt, patch_edges, tr_vertex_pmap) - : compute_constraint_vertex_distances_from_faces(cdt, mesh, tr_vertex_pmap); - - if(!options.quiet) { - std::cout << "Min distance between constraint segment and vertex: " << min_distance << '\n' - << " between segment : " - << CGAL::IO::oformat(min_va, CDT::Conforming_Dt::with_point) << " " - << CGAL::IO::oformat(min_vb, CDT::Conforming_Dt::with_point) << '\n' - << " and vertex : " - << CGAL::IO::oformat(min_vertex, CDT::Conforming_Dt::with_point) << "\n\n"; - } - cdt.check_segment_vertex_distance_or_throw(min_va, min_vb, min_vertex, min_distance, - CDT::Check_distance::NON_SQUARED_DISTANCE); -} template struct Mesh_property_maps { @@ -766,9 +608,9 @@ void dump_patches_borders(const BordersOfPatches& patch_edges, } int go(Mesh mesh, CDT_options options) { - CDT cdt; - configure_cdt_debug_options(cdt, options); + CGAL::CDT_3::Debug_options cdt_debug = cdt_debug_options(options); + // Compute bbox (used for distance checks scaling) auto bbox = compute_bounding_box(mesh, options); auto bbox_max_span = (std::max)(bbox.x_span(), (std::max)(bbox.y_span(), bbox.z_span())); @@ -781,99 +623,18 @@ int go(Mesh mesh, CDT_options options) { dump_patches_borders(patch_edges, pmaps, options.dump_patches_borders_prefix); } - auto mesh_descriptor_to_vertex_handle_pmap = get(CGAL::dynamic_vertex_property_t(), mesh); - - { - auto _ = CGAL::CDT_3_INSERT_VERTICES_TASK_guard(); - auto start_time = std::chrono::high_resolution_clock::now(); - if(options.merge_facets) { - cdt.insert_vertices_range(vertices(mesh), pmaps.mesh_vertex_point_map, mesh_descriptor_to_vertex_handle_pmap, - CGAL::parameters::vertex_is_constrained_map(pmaps.v_selected_map)); - } else { - cdt.insert_vertices_range(vertices(mesh), pmaps.mesh_vertex_point_map, mesh_descriptor_to_vertex_handle_pmap); + auto dump_mesh_with_steiner_points = [&](const auto& cdt) { + if(options.dump_after_conforming_filename.empty()) { + return; } - CDT::Cell_handle hint{}; - for(auto v: vertices(mesh)) { - if(options.merge_facets && false == get(pmaps.v_selected_map, v)) continue; - auto vh = cdt.insert(get(pmaps.mesh_vertex_point_map, v), hint, false); - hint = vh->cell(); - put(mesh_descriptor_to_vertex_handle_pmap, v, vh); - } - if(!options.quiet) { - auto timing = std::chrono::high_resolution_clock::now() - start_time; - std::cout << "[timings] inserted vertices in " - << std::chrono::duration_cast(timing).count() << " ms\n"; - std::cout << "Number of vertices: " << cdt.number_of_vertices() << "\n\n"; - } - cdt.add_bbox_points_if_not_dimension_3(); - } - - { - auto _ = CGAL::CDT_3_COMPUTE_DISTANCES_TASK_guard(); - auto start_time = std::chrono::high_resolution_clock::now(); - - auto exit_code = validate_minimum_vertex_distances(cdt, options.vertex_vertex_epsilon * bbox_max_span, options); - if(exit_code != EXIT_SUCCESS) { - return exit_code; - } - - validate_constraint_vertex_distances_or_throw(cdt, mesh, options, patch_edges, mesh_descriptor_to_vertex_handle_pmap); - - if(!options.quiet) { - auto timing = std::chrono::high_resolution_clock::now() - start_time; - std::cout << "[timings] compute distances on " - << std::chrono::duration_cast(timing).count() << " ms\n"; - } - } - - { - auto _ = CGAL::CDT_3_CONFORMING_TASK_guard(); - int poly_id = 0; - auto output_on_exit_scope_guard = CGAL::make_scope_exit(create_output_finalizer(cdt, options)); - auto start_time = std::chrono::high_resolution_clock::now(); - if(options.merge_facets) { - for(auto& edges: patch_edges) { - cdt.insert_constrained_face_defined_by_its_border_edges(std::move(edges), pmaps.mesh_vertex_point_map, - mesh_descriptor_to_vertex_handle_pmap); - } - } else { - for(auto face_descriptor : faces(mesh)) { - std::vector polygon; - const auto he = halfedge(face_descriptor, mesh); - for(auto vertex_it : CGAL::vertices_around_face(he, mesh)) { - polygon.push_back(get(pmaps.mesh_vertex_point_map, vertex_it)); - } - if(cdt.debug_polygon_insertion()) { - std::cerr << "NEW POLYGON #" << poly_id << '\n'; - } - [[maybe_unused]] auto id = cdt.insert_constrained_polygon(polygon, false); - assert(id == poly_id); - ++poly_id; - // std::ofstream dump("dump.binary.cgal"); - // CGAL::Mesh_3::save_binary_file(dump, cdt); - } - } // not merge_facets - if(!options.quiet) { - auto timing = std::chrono::high_resolution_clock::now() - start_time; - std::cout << "[timings] registered facets in " - << std::chrono::duration_cast(timing).count() << " ms\n"; - } - start_time = std::chrono::high_resolution_clock::now(); - cdt.restore_Delaunay(); - if(!options.quiet) { - std::cout << "[timings] restored Delaunay (conforming of facets borders) in " << std::chrono::duration_cast( - std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; - } - } - - if(!options.dump_after_conforming_filename.empty()) { auto _ = CGAL::CDT_3_OUTPUT_TASK_guard(); using Vertex_index = Mesh::Vertex_index; [[maybe_unused]] std::size_t time_stamp_counter = 0u; - for(auto v: cdt.finite_vertex_handles()) { + for(auto v : cdt.finite_vertex_handles()) { [[maybe_unused]] const auto time_stamp = v->time_stamp(); assert(++time_stamp_counter == time_stamp); - if(!v->ccdt_3_data().is_Steiner_vertex_on_edge()) continue; + if(!v->ccdt_3_data().is_Steiner_vertex_on_edge()) + continue; const auto [va, vb] = cdt.ancestors_of_Steiner_vertex_on_edge(v); const auto index_va = Vertex_index{static_cast(va->time_stamp() - 1)}; const auto index_vb = Vertex_index{static_cast(vb->time_stamp() - 1)}; @@ -889,36 +650,28 @@ int go(Mesh mesh, CDT_options options) { out_mesh.precision(17); out_mesh << mesh; out_mesh.close(); - } + }; + // Build CDT directly from the polygon mesh using the official API + auto build_start = std::chrono::high_resolution_clock::now(); + CDT cdt{std::invoke([&]() { + auto np = CGAL::parameters::debug(cdt_debug).visitor(std::ref(dump_mesh_with_steiner_points)); + return options.merge_facets ? CDT(mesh, np.plc_face_id(pmaps.patch_id_map)) : CDT(mesh, np); + })}; if(!options.quiet) { - std::cout << "Number of vertices after conforming: " << cdt.number_of_vertices() << "\n\n"; + auto timing = std::chrono::high_resolution_clock::now() - build_start; + std::cout << "[timings] built CDT from mesh in " + << std::chrono::duration_cast(timing).count() << " ms\n"; + std::cout << "Number of vertices: " << cdt.number_of_vertices() << "\n\n"; } + // Validation: check conforming status and validity { auto _ = CGAL::CDT_3_VALIDATION_TASK_guard(); - CGAL_assertion(!options.call_is_valid || cdt.Base_triantulation::is_valid(true)); CGAL_assertion(!options.call_is_valid || cdt.is_valid(true)); CGAL_assertion(!options.call_is_valid || cdt.is_conforming()); } - { - auto _ = CGAL::CDT_3_CDT_TASK_guard(); - auto start_time = std::chrono::high_resolution_clock::now(); - cdt.restore_constrained_Delaunay(); - if(!options.quiet) { - std::cout << "[timings] restored constrained Delaunay in " << std::chrono::duration_cast( - std::chrono::high_resolution_clock::now() - start_time).count() << " ms\n"; - std::cout << "Number of vertices after CDT: " << cdt.number_of_vertices() << "\n\n"; - } - } - - { - auto _ = CGAL::CDT_3_VALIDATION_TASK_guard(); - CGAL_assertion(!options.call_is_valid || cdt.is_conforming()); - CGAL_assertion(!options.call_is_valid || cdt.is_valid(true)); - } - return EXIT_SUCCESS; } @@ -1008,9 +761,6 @@ int bisect_errors(Mesh mesh, CDT_options options) { } int main(int argc, char* argv[]) { - CDT::Conforming_Dt::with_offset.offset = -1; - CDT::Conforming_Dt::with_point.offset = -1; - CDT::Conforming_Dt::with_point_and_info.offset = -1; std::cerr.precision(17); std::cout.precision(17);