From e7117e08d136708a12abdaf54590fdd3292b5506 Mon Sep 17 00:00:00 2001 From: Sascha Hofstetter Date: Sun, 29 Jun 2025 21:10:06 +0200 Subject: [PATCH] resulting CGAL Polygon is now returned instead of given as argument by reference in all functions. Mapping is now also argument fot dealii_tria_to_cgal_polygon but the actual changes to that function are not done yet --- include/deal.II/cgal/polygon.h | 25 +++++++++------------ source/cgal/polygon.cc | 37 ++++++++++++++++--------------- source/cgal/polygon.inst.in | 22 +++++++++---------- tests/cgal/cgal_polygon_01.cc | 12 +++------- tests/cgal/cgal_polygon_02.cc | 14 +++++------- tests/cgal/cgal_polygon_03.cc | 40 +++++++++++++++------------------- 6 files changed, 68 insertions(+), 82 deletions(-) diff --git a/include/deal.II/cgal/polygon.h b/include/deal.II/cgal/polygon.h index 0b2c810e9d..030e0ddb7b 100644 --- a/include/deal.II/cgal/polygon.h +++ b/include/deal.II/cgal/polygon.h @@ -47,14 +47,12 @@ namespace CGALWrappers * * @param[in] cell The input deal.II cell iterator * @param[in] mapping The mapping used to map the vertices of the cell - * @param[out] polygon The output CGAL::Polygon_2 */ template - void + CGAL::Polygon_2 dealii_cell_to_cgal_polygon( const typename Triangulation<2, 2>::cell_iterator &cell, - const Mapping<2, 2> &mapping, - CGAL::Polygon_2 &polygon); + const Mapping<2, 2> &mapping); @@ -64,12 +62,12 @@ namespace CGALWrappers * Triangulations that have holes are not supported. * * @param[in] tria The input deal.II triangulation - * @param[out] polygon The output CGAL::Polygon_2 + * @param[in] mapping The mapping used to map the vertices of cells */ template - void - dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria, - CGAL::Polygon_2 &polygon); + CGAL::Polygon_2 + dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria, + const Mapping<2, 2> &mapping); @@ -91,15 +89,12 @@ namespace CGALWrappers * @param[in] polygon_1 The first input CGAL::Polygon_2 * @param[in] polygon_2 The second input CGAL::Polygon_2 * @param[in] boolean_operation The input BooleanOperation - * @param[out] polygon_out The output CGAL::Polygon_2_with_holes */ template - void - compute_boolean_operation( - const CGAL::Polygon_2 &polygon_1, - const CGAL::Polygon_2 &polygon_2, - const BooleanOperation &boolean_operation, - std::vector> &polygon_out); + std::vector> + compute_boolean_operation(const CGAL::Polygon_2 &polygon_1, + const CGAL::Polygon_2 &polygon_2, + const BooleanOperation &boolean_operation); } // namespace CGALWrappers DEAL_II_NAMESPACE_CLOSE diff --git a/source/cgal/polygon.cc b/source/cgal/polygon.cc index be023257da..c450f1fd68 100644 --- a/source/cgal/polygon.cc +++ b/source/cgal/polygon.cc @@ -31,14 +31,14 @@ DEAL_II_NAMESPACE_OPEN namespace CGALWrappers { template - void + CGAL::Polygon_2 dealii_cell_to_cgal_polygon( const typename Triangulation<2, 2>::cell_iterator &cell, - const Mapping<2, 2> &mapping, - CGAL::Polygon_2 &polygon) + const Mapping<2, 2> &mapping) { - const auto &vertices = mapping.get_vertices(cell); - polygon.clear(); + CGAL::Polygon_2 polygon; + const auto &vertices = mapping.get_vertices(cell); + if (cell->reference_cell() == ReferenceCells::Triangle) { polygon.push_back( @@ -70,14 +70,15 @@ namespace CGALWrappers { DEAL_II_ASSERT_UNREACHABLE(); } + return polygon; } template - void - dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria, - CGAL::Polygon_2 &polygon) + CGAL::Polygon_2 + dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria, + const Mapping<2, 2> &mapping) { // This map holds the two vertex indices of each face. // Counterclockwise first vertex index on first position, @@ -110,8 +111,8 @@ namespace CGALWrappers } } - const auto &vertices = tria.get_vertices(); - polygon.clear(); + CGAL::Polygon_2 polygon; + const auto &vertices = tria.get_vertices(); // Vertex to start counterclockwise insertion const unsigned int start_index = face_vertex_indices.begin()->first; @@ -141,22 +142,22 @@ namespace CGALWrappers ExcMessage( "This should not occur, reasons might be a non closed boundary or a bug in this function")); } + + return polygon; } template - void - compute_boolean_operation( - const CGAL::Polygon_2 &polygon_1, - const CGAL::Polygon_2 &polygon_2, - const BooleanOperation &boolean_operation, - std::vector> &polygon_out) + std::vector> + compute_boolean_operation(const CGAL::Polygon_2 &polygon_1, + const CGAL::Polygon_2 &polygon_2, + const BooleanOperation &boolean_operation) { Assert(!(boolean_operation == BooleanOperation::compute_corefinement), ExcMessage("Corefinement has no usecase for 2D polygons")); - polygon_out.clear(); + std::vector> polygon_out; if (boolean_operation == BooleanOperation::compute_intersection) { @@ -173,6 +174,8 @@ namespace CGALWrappers polygon_out.resize(1); CGAL::join(polygon_1, polygon_2, polygon_out[0]); } + + return polygon_out; } // Explicit instantiations. diff --git a/source/cgal/polygon.inst.in b/source/cgal/polygon.inst.in index 9a62764b36..2622f01c67 100644 --- a/source/cgal/polygon.inst.in +++ b/source/cgal/polygon.inst.in @@ -16,17 +16,17 @@ for (cgal_kernel : CGAL_KERNELS) { - template void dealii_cell_to_cgal_polygon( - const typename Triangulation<2, 2>::cell_iterator &cell, - const Mapping<2, 2> &mapping, - CGAL::Polygon_2 &polygon); + template CGAL::Polygon_2 dealii_cell_to_cgal_polygon< + cgal_kernel>(const typename Triangulation<2, 2>::cell_iterator &cell, + const Mapping<2, 2> &mapping); - template void dealii_tria_to_cgal_polygon( - const Triangulation<2, 2> &tria, CGAL::Polygon_2 &polygon); + template CGAL::Polygon_2 + dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria, + const Mapping<2, 2> &mapping); - template void compute_boolean_operation( - const CGAL::Polygon_2 &polygon_1, - const CGAL::Polygon_2 &polygon_2, - const BooleanOperation &boolean_operation, - std::vector> &polygon_out); + template std::vector> + compute_boolean_operation( + const CGAL::Polygon_2 &polygon_1, + const CGAL::Polygon_2 &polygon_2, + const BooleanOperation &boolean_operation); } diff --git a/tests/cgal/cgal_polygon_01.cc b/tests/cgal/cgal_polygon_01.cc index 0a23fbe1e6..4c95dbfe3a 100644 --- a/tests/cgal/cgal_polygon_01.cc +++ b/tests/cgal/cgal_polygon_01.cc @@ -26,9 +26,7 @@ using namespace CGALWrappers; -using K = CGAL::Exact_predicates_exact_constructions_kernel; -using CGALPoint2 = CGAL::Point_2; -using CGALPolygon = CGAL::Polygon_2; +using K = CGAL::Exact_predicates_exact_constructions_kernel; void test_quadrilateral() @@ -41,10 +39,8 @@ test_quadrilateral() GridGenerator::reference_cell(tria, r_cell); - CGALPolygon poly; - const auto cell = tria.begin_active(); - dealii_cell_to_cgal_polygon(cell, *mapping, poly); + auto poly = dealii_cell_to_cgal_polygon(cell, *mapping); deallog << "Quadrilateral: " << std::endl; deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() @@ -76,10 +72,8 @@ test_triangle() GridGenerator::reference_cell(tria, r_cell); - CGALPolygon poly; - const auto cell = tria.begin_active(); - dealii_cell_to_cgal_polygon(cell, *mapping, poly); + auto poly = dealii_cell_to_cgal_polygon(cell, *mapping); deallog << "Triangle: " << std::endl; deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() diff --git a/tests/cgal/cgal_polygon_02.cc b/tests/cgal/cgal_polygon_02.cc index a3792ab9e9..18c09a20a2 100644 --- a/tests/cgal/cgal_polygon_02.cc +++ b/tests/cgal/cgal_polygon_02.cc @@ -26,15 +26,14 @@ using namespace CGALWrappers; -using K = CGAL::Exact_predicates_exact_constructions_kernel; -using CGALPoint2 = CGAL::Point_2; -using CGALPolygon = CGAL::Polygon_2; +using K = CGAL::Exact_predicates_exact_constructions_kernel; +using CGALPoint2 = CGAL::Point_2; void test_quadrilaterals() { Triangulation<2, 2> tria_in; - CGALPolygon poly; + MappingQ<2, 2> mapping(1); std::vector> names_and_args; @@ -53,7 +52,7 @@ test_quadrilaterals() deallog << "name: " << name << std::endl; GridGenerator::generate_from_name_and_arguments(tria_in, name, args); tria_in.refine_global(2); - dealii_tria_to_cgal_polygon(tria_in, poly); + auto poly = dealii_tria_to_cgal_polygon(tria_in, mapping); deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() << std::endl; @@ -76,7 +75,6 @@ test_quadrilaterals() } tria_in.clear(); - poly.clear(); deallog << std::endl; } } @@ -88,13 +86,13 @@ test_triangles() { Triangulation<2, 2> tria_quad; Triangulation<2, 2> tria_simplex; - CGALPolygon poly; + MappingQ<2, 2> mapping(1); deallog << "name: hyper_cube converted to simplex mesh" << std::endl; GridGenerator::hyper_cube(tria_quad, 0., 1., false); tria_quad.refine_global(2); GridGenerator::convert_hypercube_to_simplex_mesh(tria_quad, tria_simplex); - dealii_tria_to_cgal_polygon(tria_simplex, poly); + auto poly = dealii_tria_to_cgal_polygon(tria_simplex, mapping); deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() << std::endl; diff --git a/tests/cgal/cgal_polygon_03.cc b/tests/cgal/cgal_polygon_03.cc index cbf2c2be9b..39e0fad5bb 100644 --- a/tests/cgal/cgal_polygon_03.cc +++ b/tests/cgal/cgal_polygon_03.cc @@ -26,9 +26,8 @@ using namespace CGALWrappers; -using K = CGAL::Exact_predicates_exact_constructions_kernel; -using CGALPoint2 = CGAL::Point_2; -using CGALPolygon = CGAL::Polygon_2; +using K = CGAL::Exact_predicates_exact_constructions_kernel; + using CGALPolygonWithHoles = CGAL::Polygon_with_holes_2; void @@ -54,11 +53,9 @@ write_volumes(std::vector &poly_vec) void test(unsigned int refinement_1, unsigned int refinement_2) { - Triangulation<2, 2> tria_1; - Triangulation<2, 2> tria_2; - CGALPolygon poly_1; - CGALPolygon poly_2; - std::vector poly_out_vec; + Triangulation<2, 2> tria_1; + Triangulation<2, 2> tria_2; + MappingQ<2, 2> mapping(1); std::vector> names_and_args; @@ -79,7 +76,7 @@ test(unsigned int refinement_1, unsigned int refinement_2) names_and_args[0].first, names_and_args[0].second); tria_1.refine_global(refinement_1); - dealii_tria_to_cgal_polygon(tria_1, poly_1); + auto poly_1 = dealii_tria_to_cgal_polygon(tria_1, mapping); for (const auto &info_pair : names_and_args) { @@ -88,29 +85,28 @@ test(unsigned int refinement_1, unsigned int refinement_2) deallog << "name: " << name << std::endl; GridGenerator::generate_from_name_and_arguments(tria_2, name, args); tria_2.refine_global(refinement_2); - dealii_tria_to_cgal_polygon(tria_2, poly_2); + auto poly_2 = dealii_tria_to_cgal_polygon(tria_2, mapping); - compute_boolean_operation(poly_1, - poly_2, - BooleanOperation::compute_intersection, - poly_out_vec); + auto poly_out_vec = + compute_boolean_operation(poly_1, + poly_2, + BooleanOperation::compute_intersection); deallog << "Intersection: " << poly_out_vec.size() << " polygons" << std::endl; write_volumes(poly_out_vec); - compute_boolean_operation(poly_1, - poly_2, - BooleanOperation::compute_difference, - poly_out_vec); + poly_out_vec = + compute_boolean_operation(poly_1, + poly_2, + BooleanOperation::compute_difference); deallog << "Difference: " << poly_out_vec.size() << " polygons" << std::endl; write_volumes(poly_out_vec); - compute_boolean_operation(poly_1, - poly_2, - BooleanOperation::compute_union, - poly_out_vec); + poly_out_vec = compute_boolean_operation(poly_1, + poly_2, + BooleanOperation::compute_union); deallog << "Union: " << poly_out_vec.size() << " polygons" << std::endl; write_volumes(poly_out_vec); -- 2.39.5