* @param[in] mapping The mapping used to map the vertices of cells
*/
template <typename KernelType>
- CGAL::Polygon_2<KernelType>
+ CGAL::Polygon_with_holes_2<KernelType>
dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria,
const Mapping<2, 2> &mapping);
+ /**
+ * Constructs a Polygon_with_holes_2 from the input Polygon_2.
+ * Further polygons for holes are optional.
+ *
+ * Polygon_with_holes_2 has function .outer_boundary() as well as
+ * .holes() or .holes_begin() and .holes_end()
+ *
+ * @param[in] boundary_outside Polygon for the outer boundary
+ * @param[in] boundary_holes Polygons for the holes
+ */
+ template <typename KernelType>
+ CGAL::Polygon_with_holes_2<KernelType>
+ polygon_to_polygon_with_holes(
+ const CGAL::Polygon_2<KernelType> &boundary_outside,
+ const std::vector<CGAL::Polygon_2<KernelType>> &boundary_holes = {});
+
+
+
/**
* Perform a Regularized Boolean set-operation on two CGAL::Polygon_2.
*
*/
template <typename KernelType>
std::vector<CGAL::Polygon_with_holes_2<KernelType>>
- compute_boolean_operation(const CGAL::Polygon_2<KernelType> &polygon_1,
- const CGAL::Polygon_2<KernelType> &polygon_2,
- const BooleanOperation &boolean_operation);
+ compute_boolean_operation(
+ const CGAL::Polygon_with_holes_2<KernelType> &polygon_1,
+ const CGAL::Polygon_with_holes_2<KernelType> &polygon_2,
+ const BooleanOperation &boolean_operation);
} // namespace CGALWrappers
DEAL_II_NAMESPACE_CLOSE
std::map<unsigned int, Point<spacedim>>
get_all_vertices_at_boundary(const Triangulation<dim, spacedim> &tria);
+ /**
+ * Return a std::vector that contains a map for each closed boundary.
+ * Each closed boundary is described by a std:vector of pairs
+ * `vertex index , Point<spacedim>`.
+ *
+ * The vertices are in counter-clockwise ordering for the outer boundary.
+ * The vertices are in clockwise ordering for inner boundaries (holes).
+ *
+ * It is generally not guaranteed that the first entry of the vector is the
+ * outer boundary. However, when cell 0 is located on the outer boundary it is
+ * the case.
+ *
+ * The mapping argument enables the use of e.g. MappingQEulerian.
+ *
+ * @param[in] tria The Triangulation object.
+ * @param[in] mapping The mapping used to map the vertices of the cell
+ */
+ template <int dim, int spacedim>
+ std::vector<std::vector<std::pair<unsigned int, Point<spacedim>>>>
+ extract_ordered_boundary_vertices(
+ const Triangulation<dim, spacedim> &tria,
+ const Mapping<dim, spacedim> &mapping =
+ (ReferenceCells::get_hypercube<dim>()
+#ifndef _MSC_VER
+ .template get_default_linear_mapping<dim, spacedim>()
+#else
+ .ReferenceCell::get_default_linear_mapping<dim, spacedim>()
+#endif
+ ));
+
/**
* Remove hanging nodes from a grid. If the @p isotropic parameter is set
* to @p false (default) this function detects cells with hanging nodes and
#include <deal.II/fe/mapping.h>
-#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_tools_topology.h>
+#include <deal.II/grid/tria.h> //should be ok to delete is in grid tools
#ifdef DEAL_II_WITH_CGAL
template <typename KernelType>
- CGAL::Polygon_2<KernelType>
+ CGAL::Polygon_with_holes_2<KernelType>
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,
- // counterclockwise second vertex index on second position.
- std::map<unsigned int, unsigned int> face_vertex_indices;
- std::map<unsigned int, Point<2>> vertex_to_point;
-
- // Iterate over all active cells at the boundary
- for (const auto &cell : tria.active_cell_iterators())
+ // Outer boundary vertices counter clockwise
+ // and hole vertices clockwise ordered
+ auto boundaries =
+ GridTools::extract_ordered_boundary_vertices(tria, mapping);
+
+ CGAL::Polygon_2<KernelType> outer_boundary;
+ std::vector<CGAL::Polygon_2<KernelType>> holes;
+ holes.reserve(boundaries.size() - 1);
+
+ for (const auto &boundary : boundaries)
{
- for (const unsigned int f : cell->face_indices())
+ // Add vertices of current boundary to polygon
+ CGAL::Polygon_2<KernelType> current_polygon;
+ for (const auto &vertices : boundary)
+ {
+ current_polygon.push_back(
+ dealii_point_to_cgal_point<CGAL::Point_2<KernelType>, 2>(
+ vertices.second));
+ }
+ // Decide if outer boundary or hole
+ if (current_polygon.is_counterclockwise_oriented())
{
- if (cell->face(f)->at_boundary())
- {
- // get mapped vertices of the cell
- const auto v_mapped = mapping.get_vertices(cell);
- const unsigned int v0 = cell->face(f)->vertex_index(0);
- const unsigned int v1 = cell->face(f)->vertex_index(1);
-
- if (cell->reference_cell() == ReferenceCells::Triangle)
- {
- // add indices and first mapped vertex of the face
- vertex_to_point[v0] = v_mapped[f];
- face_vertex_indices[v0] = v1;
- }
- else if (cell->reference_cell() ==
- ReferenceCells::Quadrilateral)
- {
- // Ensure that vertex indices of the face are in
- // counterclockwise order inserted in the map.
- if (f == 0 || f == 3)
- {
- // add indices and first mapped vertex of the face
- vertex_to_point[v1] =
- v_mapped[GeometryInfo<2>::face_to_cell_vertices(f,
- 1)];
- face_vertex_indices[v1] = v0;
- }
- else
- {
- // add indices and first mapped vertex of the face
- vertex_to_point[v0] =
- v_mapped[GeometryInfo<2>::face_to_cell_vertices(f,
- 0)];
- face_vertex_indices[v0] = v1;
- }
- }
- else
- {
- DEAL_II_ASSERT_UNREACHABLE();
- }
- }
+ outer_boundary = current_polygon;
+ }
+ else
+ {
+ holes.push_back(current_polygon);
}
}
+ return CGAL::Polygon_with_holes_2<KernelType>(outer_boundary,
+ holes.begin(),
+ holes.end());
+ }
- CGAL::Polygon_2<KernelType> polygon;
-
- // 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)
- {
- const auto vertex_it = vertex_to_point.find(current_index);
- Assert(vertex_it != vertex_to_point.end(),
- ExcMessage("This should not occur, please report bug"));
- polygon.push_back(
- dealii_point_to_cgal_point<CGAL::Point_2<KernelType>, 2>(
- vertex_it->second));
- vertex_to_point.erase(vertex_it);
-
- 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.
- // This condition can be used to extend the code for
- // trianglations with holes
- 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"));
- }
-
- return polygon;
+ template <typename KernelType>
+ CGAL::Polygon_with_holes_2<KernelType>
+ polygon_to_polygon_with_holes(
+ const CGAL::Polygon_2<KernelType> &boundary_outside,
+ const std::vector<CGAL::Polygon_2<KernelType>> &boundary_holes)
+ {
+ return CGAL::Polygon_with_holes_2<KernelType>(boundary_outside,
+ boundary_holes.begin(),
+ boundary_holes.end());
}
template <typename KernelType>
std::vector<CGAL::Polygon_with_holes_2<KernelType>>
- compute_boolean_operation(const CGAL::Polygon_2<KernelType> &polygon_1,
- const CGAL::Polygon_2<KernelType> &polygon_2,
- const BooleanOperation &boolean_operation)
+ compute_boolean_operation(
+ const CGAL::Polygon_with_holes_2<KernelType> &polygon_1,
+ const CGAL::Polygon_with_holes_2<KernelType> &polygon_2,
+ const BooleanOperation &boolean_operation)
{
Assert(!(boolean_operation == BooleanOperation::compute_corefinement),
ExcMessage("Corefinement has no usecase for 2D polygons"));
cgal_kernel>(const typename Triangulation<2, 2>::cell_iterator &cell,
const Mapping<2, 2> &mapping);
- template CGAL::Polygon_2<cgal_kernel>
+ template CGAL::Polygon_with_holes_2<cgal_kernel>
dealii_tria_to_cgal_polygon<cgal_kernel>(const Triangulation<2, 2> &tria,
const Mapping<2, 2> &mapping);
+ template CGAL::Polygon_with_holes_2<cgal_kernel>
+ polygon_to_polygon_with_holes(
+ const CGAL::Polygon_2<cgal_kernel> &boundary_outside,
+ const std::vector<CGAL::Polygon_2<cgal_kernel>> &boundary_holes);
+
template std::vector<CGAL::Polygon_with_holes_2<cgal_kernel>>
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);
+ const CGAL::Polygon_with_holes_2<cgal_kernel> &polygon_1,
+ const CGAL::Polygon_with_holes_2<cgal_kernel> &polygon_2,
+ const BooleanOperation &boolean_operation);
}
+ template <>
+ std::vector<std::vector<std::pair<unsigned int, Point<2>>>>
+ extract_ordered_boundary_vertices<2, 2>(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,
+ // counterclockwise second vertex index on second position.
+ std::map<unsigned int, unsigned int> face_vertex_indices;
+ std::map<unsigned int, Point<2>> vertex_to_point;
+
+ // Iterate over all active cells at the boundary
+ for (const auto &cell : tria.active_cell_iterators())
+ {
+ for (const unsigned int f : cell->face_indices())
+ {
+ if (cell->face(f)->at_boundary())
+ {
+ // get mapped vertices of the cell
+ const auto v_mapped = mapping.get_vertices(cell);
+ const unsigned int v0 = cell->face(f)->vertex_index(0);
+ const unsigned int v1 = cell->face(f)->vertex_index(1);
+
+ if (cell->reference_cell() == ReferenceCells::Triangle)
+ {
+ // add indices and first mapped vertex of the face
+ vertex_to_point[v0] = v_mapped[f];
+ face_vertex_indices[v0] = v1;
+ }
+ else if (cell->reference_cell() ==
+ ReferenceCells::Quadrilateral)
+ {
+ // Ensure that vertex indices of the face are in
+ // counterclockwise order inserted in the map.
+ if (f == 0 || f == 3)
+ {
+ // add indices and first mapped vertex of the face
+ vertex_to_point[v1] =
+ v_mapped[GeometryInfo<2>::face_to_cell_vertices(f,
+ 1)];
+ face_vertex_indices[v1] = v0;
+ }
+ else
+ {
+ // add indices and first mapped vertex of the face
+ vertex_to_point[v0] =
+ v_mapped[GeometryInfo<2>::face_to_cell_vertices(f,
+ 0)];
+ face_vertex_indices[v0] = v1;
+ }
+ }
+ else
+ {
+ DEAL_II_ASSERT_UNREACHABLE();
+ }
+ }
+ }
+ }
+
+ std::vector<std::vector<std::pair<unsigned int, Point<2>>>> boundaries;
+ std::vector<std::pair<unsigned int, Point<2>>> current_boundary;
+
+ // Vertex to start counterclockwise insertion
+ 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
+ while (face_vertex_indices.size() > 0)
+ {
+ const auto vertex_it = vertex_to_point.find(current_index);
+ Assert(vertex_it != vertex_to_point.end(),
+ ExcMessage("This should not occur, please report bug"));
+ current_boundary.emplace_back(vertex_it->first, vertex_it->second);
+ vertex_to_point.erase(vertex_it);
+
+ 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);
+
+ // traversed one closed boundary loop
+ if (current_index == start_index)
+ {
+ boundaries.push_back(current_boundary);
+ current_boundary.clear();
+
+ if (face_vertex_indices.size() == 0)
+ {
+ break;
+ }
+
+ // Take arbitrary remaining vertex as new start
+ // for next boundary loop
+ start_index = face_vertex_indices.begin()->first;
+ current_index = start_index;
+ }
+ }
+ return boundaries;
+ }
+
+
+
template <int dim, int spacedim>
void
remove_hanging_nodes(Triangulation<dim, spacedim> &tria,
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"}};
+ {"simplex", "0.0, 0.0 ; 1.0 , 0.0 ; 0.0, 1.0"},
+ {"channel_with_cylinder", "0.03 : 2 : 2.0 : false"},
+ {"cheese", "2 , 2"}};
for (const auto &info_pair : names_and_args)
deallog << "name: " << name << std::endl;
GridGenerator::generate_from_name_and_arguments(tria_in, name, args);
tria_in.refine_global(2);
- auto poly = dealii_tria_to_cgal_polygon<K>(tria_in, mapping);
+ auto poly_w_h = dealii_tria_to_cgal_polygon<K>(tria_in, mapping);
+ auto poly = poly_w_h.outer_boundary();
deallog << "Simple polygon: " << std::boolalpha << poly.is_simple()
<< std::endl;
deallog << "deal measure: " << GridTools::volume(tria_in)
<< ", cgal measure: " << poly.area() << std::endl;
+ for (const auto &hole : poly_w_h.holes())
+ {
+ deallog << "Hole is simple polygon: " << std::boolalpha
+ << hole.is_simple() << std::endl;
+ deallog << "Hole is counterclockwise oriented polygon: "
+ << std::boolalpha << hole.is_counterclockwise_oriented()
+ << std::endl;
+ deallog << "Hole has measure: " << hole.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++)
GridGenerator::hyper_cube(tria_quad, 0., 1., false);
tria_quad.refine_global(2);
GridGenerator::convert_hypercube_to_simplex_mesh(tria_quad, tria_simplex);
- auto poly = dealii_tria_to_cgal_polygon<K>(tria_simplex, mapping);
+ auto poly_w_h = dealii_tria_to_cgal_polygon<K>(tria_simplex, mapping);
+ auto poly = poly_w_h.outer_boundary();
deallog << "Simple polygon: " << std::boolalpha << poly.is_simple()
<< std::endl;
DEAL::Counterclockwise oriented polygon: true
DEAL::deal measure: 0.500000, cgal measure: 0.500000
DEAL::
+DEAL::name: channel_with_cylinder
+DEAL::Simple polygon: true
+DEAL::Counterclockwise oriented polygon: true
+DEAL::deal measure: 0.894196, cgal measure: 0.902000
+DEAL::Hole is simple polygon: true
+DEAL::Hole is counterclockwise oriented polygon: false
+DEAL::Hole has measure: -0.00780361
+DEAL::
+DEAL::name: cheese
+DEAL::Simple polygon: true
+DEAL::Counterclockwise oriented polygon: true
+DEAL::deal measure: 21.0000, cgal measure: 25.0000
+DEAL::Hole is simple polygon: true
+DEAL::Hole is counterclockwise oriented polygon: false
+DEAL::Hole has measure: -1.00000
+DEAL::Hole is simple polygon: true
+DEAL::Hole is counterclockwise oriented polygon: false
+DEAL::Hole has measure: -1.00000
+DEAL::Hole is simple polygon: true
+DEAL::Hole is counterclockwise oriented polygon: false
+DEAL::Hole has measure: -1.00000
+DEAL::Hole is simple polygon: true
+DEAL::Hole is counterclockwise oriented polygon: false
+DEAL::Hole has measure: -1.00000
+DEAL::
DEAL::name: hyper_cube converted to simplex mesh
DEAL::Simple polygon: true
DEAL::Counterclockwise oriented polygon: true
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2005 - 2018 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.
+//
+// ------------------------------------------------------------------------
+
+// Check extract used_vertices.
+
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/tria.h>
+
+#include "../tests.h"
+
+
+void
+test()
+{
+ deallog << "dim: " << 2 << ", spacedim: " << 2 << std::endl;
+
+ Triangulation<2, 2> tria_in;
+ MappingQ<2, 2> mapping(1);
+
+ std::vector<std::pair<std::string, std::string>> names_and_args;
+
+ names_and_args = {{"hyper_cube", "0.0 : 1.0 : false"},
+ {"hyper_ball_balanced", "0.,0. : 1. "},
+ {"simplex", "0.0, 0.0 ; 1.0 , 0.0 ; 0.0, 1.0"},
+ {"channel_with_cylinder", "0.03 : 2 : 2.0 : false"},
+ {"cheese", "2 , 2"}};
+
+ 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);
+
+ auto boundaries =
+ GridTools::extract_ordered_boundary_vertices(tria_in, mapping);
+
+ unsigned int count = 0;
+ for (const auto &boundary : boundaries)
+ {
+ deallog << "Closed boundary number: " << count++ << std::endl;
+ for (const auto &v : boundary)
+ {
+ deallog << "Vertex: " << v.first << ": " << v.second << std::endl;
+ }
+ }
+
+ tria_in.clear();
+ deallog << std::endl;
+ }
+};
+
+
+int
+main()
+{
+ initlog();
+ test();
+ return 0;
+}
--- /dev/null
+
+DEAL::dim: 2, spacedim: 2
+DEAL::Name: hyper_cube
+DEAL::Closed boundary number: 0
+DEAL::Vertex: 0: 0.00000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000
+DEAL::Vertex: 3: 1.00000 1.00000
+DEAL::Vertex: 2: 0.00000 1.00000
+DEAL::
+DEAL::Name: hyper_ball_balanced
+DEAL::Closed boundary number: 0
+DEAL::Vertex: 1: 1.00000 0.00000
+DEAL::Vertex: 6: 0.707107 0.707107
+DEAL::Vertex: 5: 0.00000 1.00000
+DEAL::Vertex: 10: -0.707107 0.707107
+DEAL::Vertex: 9: -1.00000 0.00000
+DEAL::Vertex: 14: -0.707107 -0.707107
+DEAL::Vertex: 13: 0.00000 -1.00000
+DEAL::Vertex: 16: 0.707107 -0.707107
+DEAL::
+DEAL::Name: simplex
+DEAL::Closed boundary number: 0
+DEAL::Vertex: 0: 0.00000 0.00000
+DEAL::Vertex: 3: 0.500000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000
+DEAL::Vertex: 4: 0.500000 0.500000
+DEAL::Vertex: 2: 0.00000 1.00000
+DEAL::Vertex: 5: 0.00000 0.500000
+DEAL::
+DEAL::Name: channel_with_cylinder
+DEAL::Closed boundary number: 0
+DEAL::Vertex: 0: 0.00000 0.00000
+DEAL::Vertex: 1: 0.100000 0.00000
+DEAL::Vertex: 2: 0.200000 0.00000
+DEAL::Vertex: 3: 0.300000 0.00000
+DEAL::Vertex: 4: 0.400000 0.00000
+DEAL::Vertex: 5: 0.500000 0.00000
+DEAL::Vertex: 6: 0.600000 0.00000
+DEAL::Vertex: 7: 0.700000 0.00000
+DEAL::Vertex: 8: 0.800000 0.00000
+DEAL::Vertex: 9: 0.900000 0.00000
+DEAL::Vertex: 10: 1.00000 0.00000
+DEAL::Vertex: 11: 1.10000 0.00000
+DEAL::Vertex: 12: 1.20000 0.00000
+DEAL::Vertex: 13: 1.30000 0.00000
+DEAL::Vertex: 14: 1.40000 0.00000
+DEAL::Vertex: 15: 1.50000 0.00000
+DEAL::Vertex: 16: 1.60000 0.00000
+DEAL::Vertex: 17: 1.70000 0.00000
+DEAL::Vertex: 18: 1.80000 0.00000
+DEAL::Vertex: 19: 1.90000 0.00000
+DEAL::Vertex: 20: 2.00000 0.00000
+DEAL::Vertex: 21: 2.10000 0.00000
+DEAL::Vertex: 22: 2.20000 0.00000
+DEAL::Vertex: 45: 2.20000 0.102500
+DEAL::Vertex: 67: 2.20000 0.205000
+DEAL::Vertex: 90: 2.20000 0.307500
+DEAL::Vertex: 113: 2.20000 0.410000
+DEAL::Vertex: 112: 2.10000 0.410000
+DEAL::Vertex: 111: 2.00000 0.410000
+DEAL::Vertex: 110: 1.90000 0.410000
+DEAL::Vertex: 109: 1.80000 0.410000
+DEAL::Vertex: 108: 1.70000 0.410000
+DEAL::Vertex: 107: 1.60000 0.410000
+DEAL::Vertex: 106: 1.50000 0.410000
+DEAL::Vertex: 105: 1.40000 0.410000
+DEAL::Vertex: 104: 1.30000 0.410000
+DEAL::Vertex: 103: 1.20000 0.410000
+DEAL::Vertex: 102: 1.10000 0.410000
+DEAL::Vertex: 101: 1.00000 0.410000
+DEAL::Vertex: 100: 0.900000 0.410000
+DEAL::Vertex: 99: 0.800000 0.410000
+DEAL::Vertex: 98: 0.700000 0.410000
+DEAL::Vertex: 97: 0.600000 0.410000
+DEAL::Vertex: 96: 0.500000 0.410000
+DEAL::Vertex: 95: 0.400000 0.410000
+DEAL::Vertex: 94: 0.300000 0.410000
+DEAL::Vertex: 93: 0.200000 0.410000
+DEAL::Vertex: 92: 0.100000 0.410000
+DEAL::Vertex: 91: 0.00000 0.410000
+DEAL::Vertex: 68: 0.00000 0.307500
+DEAL::Vertex: 46: 0.00000 0.205000
+DEAL::Vertex: 23: 0.00000 0.102500
+DEAL::Closed boundary number: 1
+DEAL::Vertex: 130: 0.250000 0.200000
+DEAL::Vertex: 137: 0.235355 0.164645
+DEAL::Vertex: 136: 0.200000 0.150000
+DEAL::Vertex: 135: 0.164645 0.164645
+DEAL::Vertex: 134: 0.150000 0.200000
+DEAL::Vertex: 133: 0.164645 0.235355
+DEAL::Vertex: 132: 0.200000 0.250000
+DEAL::Vertex: 131: 0.235355 0.235355
+DEAL::
+DEAL::Name: cheese
+DEAL::Closed boundary number: 0
+DEAL::Vertex: 0: 0.00000 0.00000
+DEAL::Vertex: 1: 1.00000 0.00000
+DEAL::Vertex: 2: 2.00000 0.00000
+DEAL::Vertex: 3: 3.00000 0.00000
+DEAL::Vertex: 4: 4.00000 0.00000
+DEAL::Vertex: 5: 5.00000 0.00000
+DEAL::Vertex: 11: 5.00000 1.00000
+DEAL::Vertex: 17: 5.00000 2.00000
+DEAL::Vertex: 23: 5.00000 3.00000
+DEAL::Vertex: 29: 5.00000 4.00000
+DEAL::Vertex: 35: 5.00000 5.00000
+DEAL::Vertex: 34: 4.00000 5.00000
+DEAL::Vertex: 33: 3.00000 5.00000
+DEAL::Vertex: 32: 2.00000 5.00000
+DEAL::Vertex: 31: 1.00000 5.00000
+DEAL::Vertex: 30: 0.00000 5.00000
+DEAL::Vertex: 24: 0.00000 4.00000
+DEAL::Vertex: 18: 0.00000 3.00000
+DEAL::Vertex: 12: 0.00000 2.00000
+DEAL::Vertex: 6: 0.00000 1.00000
+DEAL::Closed boundary number: 1
+DEAL::Vertex: 7: 1.00000 1.00000
+DEAL::Vertex: 13: 1.00000 2.00000
+DEAL::Vertex: 14: 2.00000 2.00000
+DEAL::Vertex: 8: 2.00000 1.00000
+DEAL::Closed boundary number: 2
+DEAL::Vertex: 9: 3.00000 1.00000
+DEAL::Vertex: 15: 3.00000 2.00000
+DEAL::Vertex: 16: 4.00000 2.00000
+DEAL::Vertex: 10: 4.00000 1.00000
+DEAL::Closed boundary number: 3
+DEAL::Vertex: 19: 1.00000 3.00000
+DEAL::Vertex: 25: 1.00000 4.00000
+DEAL::Vertex: 26: 2.00000 4.00000
+DEAL::Vertex: 20: 2.00000 3.00000
+DEAL::Closed boundary number: 4
+DEAL::Vertex: 21: 3.00000 3.00000
+DEAL::Vertex: 27: 3.00000 4.00000
+DEAL::Vertex: 28: 4.00000 4.00000
+DEAL::Vertex: 22: 4.00000 3.00000
+DEAL::