From: Sascha Hofstetter Date: Mon, 30 Jun 2025 19:40:24 +0000 (+0200) Subject: GridTools extract_ordered_boundary_vertices and usage for polygons X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b453031d58fedfa2fbbd43e5b0da786a8136da87;p=dealii.git GridTools extract_ordered_boundary_vertices and usage for polygons tests for extract_ordered_boundary_vertices Added comments and changed polygon.cc Instantiation removed since not needed Doxygen and Typos --- diff --git a/include/deal.II/cgal/polygon.h b/include/deal.II/cgal/polygon.h index 030e0ddb7b..fe38bbe936 100644 --- a/include/deal.II/cgal/polygon.h +++ b/include/deal.II/cgal/polygon.h @@ -65,12 +65,30 @@ namespace CGALWrappers * @param[in] mapping The mapping used to map the vertices of cells */ template - CGAL::Polygon_2 + CGAL::Polygon_with_holes_2 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 + CGAL::Polygon_with_holes_2 + polygon_to_polygon_with_holes( + const CGAL::Polygon_2 &boundary_outside, + const std::vector> &boundary_holes = {}); + + + /** * Perform a Regularized Boolean set-operation on two CGAL::Polygon_2. * @@ -92,9 +110,10 @@ namespace CGALWrappers */ template std::vector> - compute_boolean_operation(const CGAL::Polygon_2 &polygon_1, - const CGAL::Polygon_2 &polygon_2, - const BooleanOperation &boolean_operation); + compute_boolean_operation( + const CGAL::Polygon_with_holes_2 &polygon_1, + const CGAL::Polygon_with_holes_2 &polygon_2, + const BooleanOperation &boolean_operation); } // namespace CGALWrappers DEAL_II_NAMESPACE_CLOSE diff --git a/include/deal.II/grid/grid_tools_topology.h b/include/deal.II/grid/grid_tools_topology.h index 51fe824226..a47534ef11 100644 --- a/include/deal.II/grid/grid_tools_topology.h +++ b/include/deal.II/grid/grid_tools_topology.h @@ -184,6 +184,36 @@ namespace GridTools std::map> get_all_vertices_at_boundary(const Triangulation &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`. + * + * 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 + std::vector>>> + extract_ordered_boundary_vertices( + const Triangulation &tria, + const Mapping &mapping = + (ReferenceCells::get_hypercube() +#ifndef _MSC_VER + .template get_default_linear_mapping() +#else + .ReferenceCell::get_default_linear_mapping() +#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 diff --git a/source/cgal/polygon.cc b/source/cgal/polygon.cc index 6e5cfeabe6..827cb03c4a 100644 --- a/source/cgal/polygon.cc +++ b/source/cgal/polygon.cc @@ -18,7 +18,8 @@ #include -#include +#include +#include //should be ok to delete is in grid tools #ifdef DEAL_II_WITH_CGAL @@ -76,111 +77,64 @@ namespace CGALWrappers template - CGAL::Polygon_2 + CGAL::Polygon_with_holes_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, - // counterclockwise second vertex index on second position. - std::map face_vertex_indices; - std::map> 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 outer_boundary; + std::vector> 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 current_polygon; + for (const auto &vertices : boundary) + { + current_polygon.push_back( + dealii_point_to_cgal_point, 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(outer_boundary, + holes.begin(), + holes.end()); + } - CGAL::Polygon_2 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, 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 + CGAL::Polygon_with_holes_2 + polygon_to_polygon_with_holes( + const CGAL::Polygon_2 &boundary_outside, + const std::vector> &boundary_holes) + { + return CGAL::Polygon_with_holes_2(boundary_outside, + boundary_holes.begin(), + boundary_holes.end()); } template std::vector> - compute_boolean_operation(const CGAL::Polygon_2 &polygon_1, - const CGAL::Polygon_2 &polygon_2, - const BooleanOperation &boolean_operation) + compute_boolean_operation( + const CGAL::Polygon_with_holes_2 &polygon_1, + const CGAL::Polygon_with_holes_2 &polygon_2, + const BooleanOperation &boolean_operation) { Assert(!(boolean_operation == BooleanOperation::compute_corefinement), ExcMessage("Corefinement has no usecase for 2D polygons")); diff --git a/source/cgal/polygon.inst.in b/source/cgal/polygon.inst.in index 2622f01c67..4928cedb8c 100644 --- a/source/cgal/polygon.inst.in +++ b/source/cgal/polygon.inst.in @@ -20,13 +20,18 @@ for (cgal_kernel : CGAL_KERNELS) cgal_kernel>(const typename Triangulation<2, 2>::cell_iterator &cell, const Mapping<2, 2> &mapping); - template CGAL::Polygon_2 + template CGAL::Polygon_with_holes_2 dealii_tria_to_cgal_polygon(const Triangulation<2, 2> &tria, const Mapping<2, 2> &mapping); + template CGAL::Polygon_with_holes_2 + polygon_to_polygon_with_holes( + const CGAL::Polygon_2 &boundary_outside, + const std::vector> &boundary_holes); + template std::vector> compute_boolean_operation( - const CGAL::Polygon_2 &polygon_1, - const CGAL::Polygon_2 &polygon_2, - const BooleanOperation &boolean_operation); + const CGAL::Polygon_with_holes_2 &polygon_1, + const CGAL::Polygon_with_holes_2 &polygon_2, + const BooleanOperation &boolean_operation); } diff --git a/source/grid/grid_tools_topology.cc b/source/grid/grid_tools_topology.cc index 27b97c280e..de422ad36f 100644 --- a/source/grid/grid_tools_topology.cc +++ b/source/grid/grid_tools_topology.cc @@ -1675,6 +1675,113 @@ namespace GridTools + template <> + std::vector>>> + 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 face_vertex_indices; + std::map> 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>>> boundaries; + std::vector>> 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 void remove_hanging_nodes(Triangulation &tria, diff --git a/tests/cgal/cgal_polygon_02.cc b/tests/cgal/cgal_polygon_02.cc index 18c09a20a2..f8602666a2 100644 --- a/tests/cgal/cgal_polygon_02.cc +++ b/tests/cgal/cgal_polygon_02.cc @@ -37,12 +37,12 @@ test_quadrilaterals() 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"}}; + {"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) @@ -52,7 +52,8 @@ test_quadrilaterals() 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(tria_in, mapping); + auto poly_w_h = dealii_tria_to_cgal_polygon(tria_in, mapping); + auto poly = poly_w_h.outer_boundary(); deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() << std::endl; @@ -61,6 +62,16 @@ test_quadrilaterals() 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++) @@ -92,7 +103,8 @@ test_triangles() 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(tria_simplex, mapping); + auto poly_w_h = dealii_tria_to_cgal_polygon(tria_simplex, mapping); + auto poly = poly_w_h.outer_boundary(); deallog << "Simple polygon: " << std::boolalpha << poly.is_simple() << std::endl; diff --git a/tests/cgal/cgal_polygon_02.output b/tests/cgal/cgal_polygon_02.output index c20f04614c..4b38b8b5dc 100644 --- a/tests/cgal/cgal_polygon_02.output +++ b/tests/cgal/cgal_polygon_02.output @@ -19,6 +19,31 @@ DEAL::Simple polygon: true 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 diff --git a/tests/grid/grid_tools_extract_ordered_boundary_vertices.cc b/tests/grid/grid_tools_extract_ordered_boundary_vertices.cc new file mode 100644 index 0000000000..181f39b465 --- /dev/null +++ b/tests/grid/grid_tools_extract_ordered_boundary_vertices.cc @@ -0,0 +1,75 @@ +// ------------------------------------------------------------------------ +// +// 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 + +#include +#include +#include +#include + +#include "../tests.h" + + +void +test() +{ + deallog << "dim: " << 2 << ", spacedim: " << 2 << std::endl; + + Triangulation<2, 2> tria_in; + MappingQ<2, 2> mapping(1); + + std::vector> 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; +} diff --git a/tests/grid/grid_tools_extract_ordered_boundary_vertices.output b/tests/grid/grid_tools_extract_ordered_boundary_vertices.output new file mode 100644 index 0000000000..92014b6079 --- /dev/null +++ b/tests/grid/grid_tools_extract_ordered_boundary_vertices.output @@ -0,0 +1,136 @@ + +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::