]> https://gitweb.dealii.org/ - dealii.git/commitdiff
GridTools extract_ordered_boundary_vertices and usage for polygons
authorSascha Hofstetter <hofssasc@gmail.com>
Mon, 30 Jun 2025 19:40:24 +0000 (21:40 +0200)
committerSascha Hofstetter <hofssasc@gmail.com>
Wed, 2 Jul 2025 18:32:43 +0000 (20:32 +0200)
tests for extract_ordered_boundary_vertices

Added comments and changed polygon.cc

Instantiation removed since not needed

Doxygen and Typos

include/deal.II/cgal/polygon.h
include/deal.II/grid/grid_tools_topology.h
source/cgal/polygon.cc
source/cgal/polygon.inst.in
source/grid/grid_tools_topology.cc
tests/cgal/cgal_polygon_02.cc
tests/cgal/cgal_polygon_02.output
tests/grid/grid_tools_extract_ordered_boundary_vertices.cc [new file with mode: 0644]
tests/grid/grid_tools_extract_ordered_boundary_vertices.output [new file with mode: 0644]

index 030e0ddb7b2d7dbb95fa0adce88719df3bf156d6..fe38bbe936a4f9129146229f1403a42ef45ce036 100644 (file)
@@ -65,12 +65,30 @@ namespace CGALWrappers
    * @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.
    *
@@ -92,9 +110,10 @@ namespace CGALWrappers
    */
   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
index 51fe8242269d87f964fab13b62226fd0d8d57ea5..a47534ef11680f20118e2736efe3160d5d937600 100644 (file)
@@ -184,6 +184,36 @@ namespace GridTools
   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
index 6e5cfeabe6d304a9513f70240a8de64eec35baa7..827cb03c4ada7ac7185779b591375e730e0bba37 100644 (file)
@@ -18,7 +18,8 @@
 
 #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
 
@@ -76,111 +77,64 @@ namespace CGALWrappers
 
 
   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"));
index 2622f01c67ceb296b93ca602d6ea960124e8fbbc..4928cedb8cf610c484aad9954c8d0ecd88636c68 100644 (file)
@@ -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<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);
   }
index 27b97c280eaa775e6b724c98146b852c01f36921..de422ad36fda0d184f9f8daf2ea9dd70a82e4065 100644 (file)
@@ -1675,6 +1675,113 @@ namespace GridTools
 
 
 
+  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,
index 18c09a20a2e9c6849462bd90653d92f62e42dabc..f8602666a2626933f9a16d640dd15b5702072294 100644 (file)
@@ -37,12 +37,12 @@ test_quadrilaterals()
 
   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)
@@ -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<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;
@@ -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<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;
index c20f04614cfa715faaba852d164123ce3039f193..4b38b8b5dcfc6c1905cad693e6734a4864613f11 100644 (file)
@@ -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 (file)
index 0000000..181f39b
--- /dev/null
@@ -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 <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;
+}
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 (file)
index 0000000..92014b6
--- /dev/null
@@ -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::

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.