From: Sascha Hofstetter Date: Tue, 10 Jun 2025 12:02:42 +0000 (+0200) Subject: cgal polygon header with tests X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cc3e02d9773d8c8ceac753438b82e704759d0c02;p=dealii.git cgal polygon header with tests added licence to tests add test case for boolean improved test cases Template initialization in cc file document and rename arguments make indent add doxygen documentation fix renaming Adding comments and changes requested in PR typo --- diff --git a/include/deal.II/cgal/polygon.h b/include/deal.II/cgal/polygon.h new file mode 100644 index 0000000000..0b2c810e9d --- /dev/null +++ b/include/deal.II/cgal/polygon.h @@ -0,0 +1,116 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2025 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +#ifndef dealii_cgal_polygon_h +#define dealii_cgal_polygon_h + +#include + +#ifdef DEAL_II_WITH_CGAL +# include +# include + +# include +# include + + +DEAL_II_NAMESPACE_OPEN + +namespace CGALWrappers +{ + /** + * Build a CGAL::Polygon_2 from a deal.II cell. + * + * The class Polygon_2 is a wrapper around a container of points that can + * be used to represent polygons. + * The points are added in counterclockwise order to a Polygon_2 + * + * More information on this class is available at + * https://doc.cgal.org/latest/Polygon/index.html + * + * The functions are for 2D triangulations in 2D space. + * Projecting 3D points is possible with CGAL but not implemented. + * + * The generated boundary representation can be used to perform + * regularized boolean set-operations. + * + * @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 + dealii_cell_to_cgal_polygon( + const typename Triangulation<2, 2>::cell_iterator &cell, + const Mapping<2, 2> &mapping, + CGAL::Polygon_2 &polygon); + + + + /** + * Convert a deal.II triangulation to a CGAL::Polygon_2. + * + * Triangulations that have holes are not supported. + * + * @param[in] tria The input deal.II triangulation + * @param[out] polygon The output CGAL::Polygon_2 + */ + template + void + dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria, + CGAL::Polygon_2 &polygon); + + + + /** + * Perform a Regularized Boolean set-operation on two CGAL::Polygon_2. + * + * More information on this class is available at + * https://doc.cgal.org/latest/Boolean_set_operations_2/index.html + * + * The output is generally a vector of CGAL::Polygon_2_with_holes. + * + * For the union operation the vector will always have length one. + * + * For the difference operation the second polygon is subtracted + * from the first one. + * + * Corefinement is not supported as boolean operation. + * + * @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); +} // namespace CGALWrappers + +DEAL_II_NAMESPACE_CLOSE + +#else + +// Make sure the scripts that create the C++20 module input files have +// something to latch on if the preprocessor #ifdef above would +// otherwise lead to an empty content of the file. +DEAL_II_NAMESPACE_OPEN +DEAL_II_NAMESPACE_CLOSE + +#endif +#endif diff --git a/source/cgal/CMakeLists.txt b/source/cgal/CMakeLists.txt index af68ec1ad6..3503ecc161 100644 --- a/source/cgal/CMakeLists.txt +++ b/source/cgal/CMakeLists.txt @@ -31,12 +31,14 @@ if(DEAL_II_WITH_CGAL) ${_src} surface_mesh.cc intersections.cc + polygon.cc ) set(_inst ${_inst} surface_mesh.inst.in intersections.inst.in + polygon.inst.in ) endif() diff --git a/source/cgal/polygon.cc b/source/cgal/polygon.cc new file mode 100644 index 0000000000..be023257da --- /dev/null +++ b/source/cgal/polygon.cc @@ -0,0 +1,194 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2025 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +#include + +#include + +#include + +#include + +#ifdef DEAL_II_WITH_CGAL + +# include + +DEAL_II_NAMESPACE_OPEN + +# ifndef DOXYGEN +// Template implementations +namespace CGALWrappers +{ + template + void + dealii_cell_to_cgal_polygon( + const typename Triangulation<2, 2>::cell_iterator &cell, + const Mapping<2, 2> &mapping, + CGAL::Polygon_2 &polygon) + { + const auto &vertices = mapping.get_vertices(cell); + polygon.clear(); + if (cell->reference_cell() == ReferenceCells::Triangle) + { + polygon.push_back( + CGALWrappers::dealii_point_to_cgal_point, + 2>(vertices[0])); + polygon.push_back( + CGALWrappers::dealii_point_to_cgal_point, + 2>(vertices[1])); + polygon.push_back( + CGALWrappers::dealii_point_to_cgal_point, + 2>(vertices[2])); + } + else if (cell->reference_cell() == ReferenceCells::Quadrilateral) + { + polygon.push_back( + CGALWrappers::dealii_point_to_cgal_point, + 2>(vertices[0])); + polygon.push_back( + CGALWrappers::dealii_point_to_cgal_point, + 2>(vertices[1])); + polygon.push_back( + CGALWrappers::dealii_point_to_cgal_point, + 2>(vertices[3])); + polygon.push_back( + CGALWrappers::dealii_point_to_cgal_point, + 2>(vertices[2])); + } + else + { + DEAL_II_ASSERT_UNREACHABLE(); + } + } + + + + template + void + dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria, + CGAL::Polygon_2 &polygon) + { + // This map holds the two vertex indices of each face. + // Counterclockwise first vertex index on first position, + // counterclockwise second vertex index on second position. + std::map face_vertex_indices; + + // Iterate over all active cells at the boundary + for (const auto &cell : tria.active_cell_iterators()) + { + for (const unsigned int i : cell->face_indices()) + { + const typename Triangulation<2, 2>::face_iterator &face = + cell->face(i); + // Ensure that vertex indices of the face are in + // counterclockwise order inserted in the map. + if (face->at_boundary()) + { + if (cell->reference_cell() == ReferenceCells::Quadrilateral && + (i == 0 || i == 3)) + { + face_vertex_indices[face->vertex_index(1)] = + face->vertex_index(0); + } + else + { + face_vertex_indices[face->vertex_index(0)] = + face->vertex_index(1); + } + } + } + } + + const auto &vertices = tria.get_vertices(); + polygon.clear(); + + // Vertex to start counterclockwise insertion + const unsigned int start_index = face_vertex_indices.begin()->first; + unsigned int current_index = start_index; + + // As long as still entries in the map, use last vertex index to + // find next vertex index in counterclockwise order + while (face_vertex_indices.size() > 0) + { + polygon.push_back( + dealii_point_to_cgal_point, 2>( + vertices[current_index])); + + const auto it = face_vertex_indices.find(current_index); + + // If the boundary is one closed loop, the next vertex index + // must exist as key until the map is empty. + Assert(it != face_vertex_indices.end(), + ExcMessage("Triangulation might contain holes")); + + current_index = it->second; + face_vertex_indices.erase(it); + + // Ensure that last vertex index is the start index + Assert( + !(face_vertex_indices.size() == 0 && current_index != start_index), + ExcMessage( + "This should not occur, reasons might be a non closed boundary or a bug in this function")); + } + } + + + + 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) + { + Assert(!(boolean_operation == BooleanOperation::compute_corefinement), + ExcMessage("Corefinement has no usecase for 2D polygons")); + + polygon_out.clear(); + + if (boolean_operation == BooleanOperation::compute_intersection) + { + CGAL::intersection(polygon_1, + polygon_2, + std::back_inserter(polygon_out)); + } + else if (boolean_operation == BooleanOperation::compute_difference) + { + CGAL::difference(polygon_1, polygon_2, std::back_inserter(polygon_out)); + } + else if (boolean_operation == BooleanOperation::compute_union) + { + polygon_out.resize(1); + CGAL::join(polygon_1, polygon_2, polygon_out[0]); + } + } + +// Explicit instantiations. +// +// We don't build the instantiations.inst file if deal.II isn't +// configured with CGAL, but doxygen doesn't know that and tries to +// find that file anyway for parsing -- which then of course it fails +// on. So exclude the following from doxygen consideration. +# ifndef DOXYGEN +# include "cgal/polygon.inst" +# endif + +} // namespace CGALWrappers +# endif + + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/cgal/polygon.inst.in b/source/cgal/polygon.inst.in new file mode 100644 index 0000000000..9a62764b36 --- /dev/null +++ b/source/cgal/polygon.inst.in @@ -0,0 +1,32 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2025 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + + + +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 void dealii_tria_to_cgal_polygon( + const Triangulation<2, 2> &tria, CGAL::Polygon_2 &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); + } diff --git a/tests/cgal/cgal_polygon_01.cc b/tests/cgal/cgal_polygon_01.cc new file mode 100644 index 0000000000..0a23fbe1e6 --- /dev/null +++ b/tests/cgal/cgal_polygon_01.cc @@ -0,0 +1,108 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2025 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Convert a deal.II reference cell to a cgal polygon_2. +#include + +#include + +#include + +#include +#include + +#include "../tests.h" + +using namespace CGALWrappers; + +using K = CGAL::Exact_predicates_exact_constructions_kernel; +using CGALPoint2 = CGAL::Point_2; +using CGALPolygon = CGAL::Polygon_2; + +void +test_quadrilateral() +{ + using namespace ReferenceCells; + ReferenceCell r_cell = Quadrilateral; + Triangulation<2, 2> tria; + + const auto mapping = r_cell.template get_default_mapping<2, 2>(1); + + GridGenerator::reference_cell(tria, r_cell); + + CGALPolygon poly; + + const auto cell = tria.begin_active(); + dealii_cell_to_cgal_polygon(cell, *mapping, poly); + + deallog << "Quadrilateral: " << std::endl; + deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() + << std::endl; + deallog << "Counterclockwise oriented polygon: " << std::boolalpha + << poly.is_counterclockwise_oriented() << std::endl; + deallog << "deal vertices: " << cell->n_vertices() << ", cgal vertices " + << poly.size() << std::endl; + deallog << "deal measure: " << GridTools::volume(tria) + << ", cgal measure: " << poly.area() << std::endl; + + for (size_t i = 0; i < poly.size(); ++i) + { + deallog << "Polygon vertex " << i << ": " << poly[i] << std::endl; + } + deallog << std::endl; +} + + + +void +test_triangle() +{ + using namespace ReferenceCells; + ReferenceCell r_cell = Triangle; + Triangulation<2, 2> tria; + + const auto mapping = r_cell.template get_default_mapping<2, 2>(1); + + GridGenerator::reference_cell(tria, r_cell); + + CGALPolygon poly; + + const auto cell = tria.begin_active(); + dealii_cell_to_cgal_polygon(cell, *mapping, poly); + + deallog << "Triangle: " << std::endl; + deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() + << std::endl; + deallog << "Counterclockwise oriented polygon: " << std::boolalpha + << poly.is_counterclockwise_oriented() << std::endl; + deallog << "deal vertices: " << cell->n_vertices() << ", cgal vertices " + << poly.size() << std::endl; + deallog << "deal measure: " << GridTools::volume(tria) + << ", cgal measure: " << poly.area() << std::endl; + + for (size_t i = 0; i < poly.size(); ++i) + { + deallog << "Polygon vertex " << i << ": " << poly[i] << std::endl; + } +} + + + +int +main() +{ + initlog(); + test_quadrilateral(); + test_triangle(); +} diff --git a/tests/cgal/cgal_polygon_01.output b/tests/cgal/cgal_polygon_01.output new file mode 100644 index 0000000000..989b3eae4d --- /dev/null +++ b/tests/cgal/cgal_polygon_01.output @@ -0,0 +1,19 @@ + +DEAL::Quadrilateral: +DEAL::Simple polygon: true +DEAL::Counterclockwise oriented polygon: true +DEAL::deal vertices: 4, cgal vertices 4 +DEAL::deal measure: 1.00000, cgal measure: 1.00000 +DEAL::Polygon vertex 0: 0.00000 0.00000 +DEAL::Polygon vertex 1: 1.00000 0.00000 +DEAL::Polygon vertex 2: 1.00000 1.00000 +DEAL::Polygon vertex 3: 0.00000 1.00000 +DEAL:: +DEAL::Triangle: +DEAL::Simple polygon: true +DEAL::Counterclockwise oriented polygon: true +DEAL::deal vertices: 3, cgal vertices 3 +DEAL::deal measure: 0.500000, cgal measure: 0.500000 +DEAL::Polygon vertex 0: 0.00000 0.00000 +DEAL::Polygon vertex 1: 1.00000 0.00000 +DEAL::Polygon vertex 2: 0.00000 1.00000 diff --git a/tests/cgal/cgal_polygon_02.cc b/tests/cgal/cgal_polygon_02.cc new file mode 100644 index 0000000000..a3792ab9e9 --- /dev/null +++ b/tests/cgal/cgal_polygon_02.cc @@ -0,0 +1,131 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2025 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Convert a deal.II triangulation to a CGAL polygon_2. +#include + +#include + +#include + +#include +#include + +#include "../tests.h" + +using namespace CGALWrappers; + +using K = CGAL::Exact_predicates_exact_constructions_kernel; +using CGALPoint2 = CGAL::Point_2; +using CGALPolygon = CGAL::Polygon_2; + +void +test_quadrilaterals() +{ + Triangulation<2, 2> tria_in; + CGALPolygon poly; + + std::vector> names_and_args; + + // only for meshes that have no holes + // extension to CGAL::Polygon_with_holes possible + names_and_args = {{"hyper_cube", "0.0 : 1.0 : false"}, + {"hyper_ball_balanced", "0.,0. : 1. "}, + {"hyper_L", "0.0 : 1.0 : false"}, + {"simplex", "0.0, 0.0 ; 1.0 , 0.0 ; 0.0, 1.0"}}; + + + for (const auto &info_pair : names_and_args) + { + auto name = info_pair.first; + auto args = info_pair.second; + 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); + + deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() + << std::endl; + deallog << "Counterclockwise oriented polygon: " << std::boolalpha + << poly.is_counterclockwise_oriented() << std::endl; + deallog << "deal measure: " << GridTools::volume(tria_in) + << ", cgal measure: " << poly.area() << std::endl; + + // every vertex of the tria should be inside or on the boundary + const auto &vertices = tria_in.get_vertices(); + for (unsigned int i = 0; i < vertices.size(); i++) + { + auto bounded_side = poly.bounded_side( + dealii_point_to_cgal_point(vertices[i])); + Assert( + bounded_side == CGAL::ON_BOUNDED_SIDE || + bounded_side == CGAL::ON_BOUNDARY, + ExcMessage( + "The polygon area doesn't contain all triangulation points")); + } + + tria_in.clear(); + poly.clear(); + deallog << std::endl; + } +} + + + +void +test_triangles() +{ + Triangulation<2, 2> tria_quad; + Triangulation<2, 2> tria_simplex; + CGALPolygon poly; + + 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); + + deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() + << std::endl; + deallog << "Counterclockwise oriented polygon: " << std::boolalpha + << poly.is_counterclockwise_oriented() << std::endl; + deallog << "deal measure: " << GridTools::volume(tria_simplex) + << ", cgal measure: " << poly.area() << std::endl; + + // every vertex of the tria should be inside or on the boundary + const auto &vertices = tria_simplex.get_vertices(); + for (unsigned int i = 0; i < vertices.size(); i++) + { + auto bounded_side = poly.bounded_side( + dealii_point_to_cgal_point(vertices[i])); + Assert(bounded_side == CGAL::ON_BOUNDED_SIDE || + bounded_side == CGAL::ON_BOUNDARY, + ExcMessage( + "The polygon area doesn't contain all triangulation points")); + } + + tria_quad.clear(); + tria_simplex.clear(); + poly.clear(); +} + + + +int +main() +{ + initlog(); + test_quadrilaterals(); + test_triangles(); +} diff --git a/tests/cgal/cgal_polygon_02.output b/tests/cgal/cgal_polygon_02.output new file mode 100644 index 0000000000..c20f04614c --- /dev/null +++ b/tests/cgal/cgal_polygon_02.output @@ -0,0 +1,25 @@ + +DEAL::name: hyper_cube +DEAL::Simple polygon: true +DEAL::Counterclockwise oriented polygon: true +DEAL::deal measure: 1.00000, cgal measure: 1.00000 +DEAL:: +DEAL::name: hyper_ball_balanced +DEAL::Simple polygon: true +DEAL::Counterclockwise oriented polygon: true +DEAL::deal measure: 3.12145, cgal measure: 3.12145 +DEAL:: +DEAL::name: hyper_L +DEAL::Simple polygon: true +DEAL::Counterclockwise oriented polygon: true +DEAL::deal measure: 0.750000, cgal measure: 0.750000 +DEAL:: +DEAL::name: simplex +DEAL::Simple polygon: true +DEAL::Counterclockwise oriented polygon: true +DEAL::deal measure: 0.500000, cgal measure: 0.500000 +DEAL:: +DEAL::name: hyper_cube converted to simplex mesh +DEAL::Simple polygon: true +DEAL::Counterclockwise oriented polygon: true +DEAL::deal measure: 1.00000, cgal measure: 1.00000 diff --git a/tests/cgal/cgal_polygon_03.cc b/tests/cgal/cgal_polygon_03.cc new file mode 100644 index 0000000000..cbf2c2be9b --- /dev/null +++ b/tests/cgal/cgal_polygon_03.cc @@ -0,0 +1,130 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2025 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + +// Boolean operations on polygons +#include + +#include + +#include + +#include +#include + +#include "../tests.h" + +using namespace CGALWrappers; + +using K = CGAL::Exact_predicates_exact_constructions_kernel; +using CGALPoint2 = CGAL::Point_2; +using CGALPolygon = CGAL::Polygon_2; +using CGALPolygonWithHoles = CGAL::Polygon_with_holes_2; + +void +write_volumes(std::vector &poly_vec) +{ + deallog << "With volumes: "; + for (const auto &polygon : poly_vec) + { + if (polygon.has_holes()) + { + for (const auto &hole : polygon.holes()) + { + deallog << hole.area() << " "; + } + } + deallog << polygon.outer_boundary().area() << " "; + } + deallog << std::endl; +} + + + +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; + + std::vector> names_and_args; + + // only for operations on meshes that have no holes + // extension to CGAL::Polygon_with_holes possible + names_and_args = {{"hyper_cube", "0.0 : 1.0 : false"}, + {"hyper_cube", "-0.1 : 0.9 : false"}, + {"hyper_rectangle", "0.0, -0.1 : 1.0 , 0.5 : false"}, + {"hyper_rectangle", "-0.1, -0.1 : 0.5 , 1.1 : false"}, + {"hyper_rectangle", "-0.1, 0.5 : 0.0 , 1.1 : false"}, + {"simplex", "0.0, 0.0 ; 1.0 , 0.0 ; 0.0, 1.0"}, + {"simplex", "-0.1, -0.1 ; 0.9 , -0.1 ; -0.1, 0.9"}, + {"simplex", "-0.5, -0.5 ; 0.5 , 0.5 ; -0.5, 1.5"}, + {"hyper_ball_balanced", "0.0,0.0 : 1.0"}, + {"hyper_ball_balanced", "-0.1 ,1.0 : 0.5"}}; + + GridGenerator::generate_from_name_and_arguments(tria_1, + 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); + + for (const auto &info_pair : names_and_args) + { + auto name = info_pair.first; + auto args = info_pair.second; + 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); + + compute_boolean_operation(poly_1, + poly_2, + BooleanOperation::compute_intersection, + poly_out_vec); + 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); + 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); + deallog << "Union: " << poly_out_vec.size() << " polygons" << std::endl; + write_volumes(poly_out_vec); + + tria_2.clear(); + poly_2.clear(); + deallog << std::endl; + } +} + + +int +main() +{ + initlog(); + test(0, 0); + test(2, 2); +} diff --git a/tests/cgal/cgal_polygon_03.output b/tests/cgal/cgal_polygon_03.output new file mode 100644 index 0000000000..2de58bb072 --- /dev/null +++ b/tests/cgal/cgal_polygon_03.output @@ -0,0 +1,161 @@ + +DEAL::name: hyper_cube +DEAL::Intersection: 1 polygons +DEAL::With volumes: 1.00000 +DEAL::Difference: 0 polygons +DEAL::With volumes: +DEAL::Union: 1 polygons +DEAL::With volumes: 1.00000 +DEAL:: +DEAL::name: hyper_cube +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.810000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.190000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.19000 +DEAL:: +DEAL::name: hyper_rectangle +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.10000 +DEAL:: +DEAL::name: hyper_rectangle +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.22000 +DEAL:: +DEAL::name: hyper_rectangle +DEAL::Intersection: 0 polygons +DEAL::With volumes: +DEAL::Difference: 1 polygons +DEAL::With volumes: 1.00000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.06000 +DEAL:: +DEAL::name: simplex +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.00000 +DEAL:: +DEAL::name: simplex +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.320000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.680000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.18000 +DEAL:: +DEAL::name: simplex +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.250000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.750000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.75000 +DEAL:: +DEAL::name: hyper_ball_balanced +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.707107 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.292893 +DEAL::Union: 1 polygons +DEAL::With volumes: 3.12132 +DEAL:: +DEAL::name: hyper_ball_balanced +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.128848 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.871152 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.57826 +DEAL:: +DEAL::name: hyper_cube +DEAL::Intersection: 1 polygons +DEAL::With volumes: 1.00000 +DEAL::Difference: 0 polygons +DEAL::With volumes: +DEAL::Union: 1 polygons +DEAL::With volumes: 1.00000 +DEAL:: +DEAL::name: hyper_cube +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.810000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.190000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.19000 +DEAL:: +DEAL::name: hyper_rectangle +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.10000 +DEAL:: +DEAL::name: hyper_rectangle +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.22000 +DEAL:: +DEAL::name: hyper_rectangle +DEAL::Intersection: 0 polygons +DEAL::With volumes: +DEAL::Difference: 1 polygons +DEAL::With volumes: 1.00000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.06000 +DEAL:: +DEAL::name: simplex +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.500000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.00000 +DEAL:: +DEAL::name: simplex +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.320000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.680000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.18000 +DEAL:: +DEAL::name: simplex +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.250000 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.750000 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.75000 +DEAL:: +DEAL::name: hyper_ball_balanced +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.780361 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.219639 +DEAL::Union: 1 polygons +DEAL::With volumes: 3.34108 +DEAL:: +DEAL::name: hyper_ball_balanced +DEAL::Intersection: 1 polygons +DEAL::With volumes: 0.145583 +DEAL::Difference: 1 polygons +DEAL::With volumes: 0.854417 +DEAL::Union: 1 polygons +DEAL::With volumes: 1.63478 +DEAL::