--- /dev/null
+// ------------------------------------------------------------------------
+//
+// 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 <deal.II/base/config.h>
+
+#ifdef DEAL_II_WITH_CGAL
+# include <deal.II/cgal/point_conversion.h>
+# include <deal.II/cgal/utilities.h>
+
+# include <CGAL/Polygon_2.h>
+# include <CGAL/Polygon_with_holes_2.h>
+
+
+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 <typename KernelType>
+ void
+ dealii_cell_to_cgal_polygon(
+ const typename Triangulation<2, 2>::cell_iterator &cell,
+ const Mapping<2, 2> &mapping,
+ CGAL::Polygon_2<KernelType> &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 <typename KernelType>
+ void
+ dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria,
+ CGAL::Polygon_2<KernelType> &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 <typename KernelType>
+ void
+ compute_boolean_operation(
+ const CGAL::Polygon_2<KernelType> &polygon_1,
+ const CGAL::Polygon_2<KernelType> &polygon_2,
+ const BooleanOperation &boolean_operation,
+ std::vector<CGAL::Polygon_with_holes_2<KernelType>> &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
${_src}
surface_mesh.cc
intersections.cc
+ polygon.cc
)
set(_inst
${_inst}
surface_mesh.inst.in
intersections.inst.in
+ polygon.inst.in
)
endif()
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// 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 <deal.II/base/config.h>
+
+#include <deal.II/cgal/polygon.h>
+
+#include <deal.II/fe/mapping.h>
+
+#include <deal.II/grid/tria.h>
+
+#ifdef DEAL_II_WITH_CGAL
+
+# include <CGAL/Boolean_set_operations_2.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+# ifndef DOXYGEN
+// Template implementations
+namespace CGALWrappers
+{
+ template <typename KernelType>
+ void
+ dealii_cell_to_cgal_polygon(
+ const typename Triangulation<2, 2>::cell_iterator &cell,
+ const Mapping<2, 2> &mapping,
+ CGAL::Polygon_2<KernelType> &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<CGAL::Point_2<KernelType>,
+ 2>(vertices[0]));
+ polygon.push_back(
+ CGALWrappers::dealii_point_to_cgal_point<CGAL::Point_2<KernelType>,
+ 2>(vertices[1]));
+ polygon.push_back(
+ CGALWrappers::dealii_point_to_cgal_point<CGAL::Point_2<KernelType>,
+ 2>(vertices[2]));
+ }
+ else if (cell->reference_cell() == ReferenceCells::Quadrilateral)
+ {
+ polygon.push_back(
+ CGALWrappers::dealii_point_to_cgal_point<CGAL::Point_2<KernelType>,
+ 2>(vertices[0]));
+ polygon.push_back(
+ CGALWrappers::dealii_point_to_cgal_point<CGAL::Point_2<KernelType>,
+ 2>(vertices[1]));
+ polygon.push_back(
+ CGALWrappers::dealii_point_to_cgal_point<CGAL::Point_2<KernelType>,
+ 2>(vertices[3]));
+ polygon.push_back(
+ CGALWrappers::dealii_point_to_cgal_point<CGAL::Point_2<KernelType>,
+ 2>(vertices[2]));
+ }
+ else
+ {
+ DEAL_II_ASSERT_UNREACHABLE();
+ }
+ }
+
+
+
+ template <typename KernelType>
+ void
+ dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria,
+ CGAL::Polygon_2<KernelType> &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<unsigned int, unsigned int> 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<CGAL::Point_2<KernelType>, 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 <typename KernelType>
+ void
+ compute_boolean_operation(
+ const CGAL::Polygon_2<KernelType> &polygon_1,
+ const CGAL::Polygon_2<KernelType> &polygon_2,
+ const BooleanOperation &boolean_operation,
+ std::vector<CGAL::Polygon_with_holes_2<KernelType>> &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
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// 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<cgal_kernel>(
+ const typename Triangulation<2, 2>::cell_iterator &cell,
+ const Mapping<2, 2> &mapping,
+ CGAL::Polygon_2<cgal_kernel> &polygon);
+
+ template void dealii_tria_to_cgal_polygon<cgal_kernel>(
+ const Triangulation<2, 2> &tria, CGAL::Polygon_2<cgal_kernel> &polygon);
+
+ template void compute_boolean_operation<cgal_kernel>(
+ const CGAL::Polygon_2<cgal_kernel> &polygon_1,
+ const CGAL::Polygon_2<cgal_kernel> &polygon_2,
+ const BooleanOperation &boolean_operation,
+ std::vector<CGAL::Polygon_with_holes_2<cgal_kernel>> &polygon_out);
+ }
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// 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 <deal.II/base/config.h>
+
+#include <deal.II/cgal/polygon.h>
+
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include "../tests.h"
+
+using namespace CGALWrappers;
+
+using K = CGAL::Exact_predicates_exact_constructions_kernel;
+using CGALPoint2 = CGAL::Point_2<K>;
+using CGALPolygon = CGAL::Polygon_2<K>;
+
+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();
+}
--- /dev/null
+
+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
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// 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 <deal.II/base/config.h>
+
+#include <deal.II/cgal/polygon.h>
+
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include "../tests.h"
+
+using namespace CGALWrappers;
+
+using K = CGAL::Exact_predicates_exact_constructions_kernel;
+using CGALPoint2 = CGAL::Point_2<K>;
+using CGALPolygon = CGAL::Polygon_2<K>;
+
+void
+test_quadrilaterals()
+{
+ Triangulation<2, 2> tria_in;
+ CGALPolygon poly;
+
+ std::vector<std::pair<std::string, std::string>> 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<CGALPoint2, 2>(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<CGALPoint2, 2>(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();
+}
--- /dev/null
+
+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
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// 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 <deal.II/base/config.h>
+
+#include <deal.II/cgal/polygon.h>
+
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include "../tests.h"
+
+using namespace CGALWrappers;
+
+using K = CGAL::Exact_predicates_exact_constructions_kernel;
+using CGALPoint2 = CGAL::Point_2<K>;
+using CGALPolygon = CGAL::Polygon_2<K>;
+using CGALPolygonWithHoles = CGAL::Polygon_with_holes_2<K>;
+
+void
+write_volumes(std::vector<CGALPolygonWithHoles> &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<CGALPolygonWithHoles> poly_out_vec;
+
+ std::vector<std::pair<std::string, std::string>> 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);
+}
--- /dev/null
+
+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::