]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix create_construction_data_from_triangulation for periodic boundaries 8982/head
authorPeter Munch <peterrmuench@gmail.com>
Tue, 29 Oct 2019 13:11:34 +0000 (14:11 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Sun, 3 Nov 2019 21:47:14 +0000 (22:47 +0100)
include/deal.II/grid/grid_tools.h
source/distributed/fully_distributed_tria_util.cc
source/distributed/tria_base.cc
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
source/grid/grid_tools_dof_handlers.cc
tests/fullydistributed_grids/copy_serial_tria_04.cc
tests/fullydistributed_grids/copy_serial_tria_04.mpirun=1.output
tests/fullydistributed_grids/copy_serial_tria_04.mpirun=4.output [new file with mode: 0644]
tests/fullydistributed_grids/copy_serial_tria_04.mpirun=5.output

index be004daf932d2ccf1b5c57193c1a182e25ddb64e..bbb4808f52e412e6e9e915a4844cabfcb0271b19 100644 (file)
@@ -2940,6 +2940,32 @@ namespace GridTools
     const std::vector<BoundingBox<spacedim>> &local_description,
     MPI_Comm                                  mpi_communicator);
 
+  /**
+   * Collect for a given triangulation all locally relevant vertices that
+   * coincide due to periodicity.
+   *
+   * Coinciding vertices are put into a group, e.g.: [1, 25, 51], which is
+   * labeled by an arbitrary element from it, e.g.: "1". All coinciding vertices
+   * store the label to its group, so that they can quickly access all the
+   * coinciding vertices in that group: e.g.: 51 ->  "1" -> [1, 25, 51]
+   *
+   * @param[in] tria Triangulation.
+   * @param[out] coinciding_vertex_groups A map of equivalence classes (of
+   *             coinciding vertices) labeled by an arbitrary element from them.
+   *             Vertices not coinciding are ignored.
+   * @param[out] vertex_to_coinciding_vertex_group Map of a vertex to the label
+   *             of a group of coinciding vertices. Vertices not contained in
+   *             this vector are not coinciding with any other vertex.
+   *
+   * @author Peter Munch, 2019.
+   */
+  template <int dim, int spacedim>
+  void
+  collect_coinciding_vertices(
+    const Triangulation<dim, spacedim> &               tria,
+    std::map<unsigned int, std::vector<unsigned int>> &coinciding_vertex_groups,
+    std::map<unsigned int, unsigned int> &vertex_to_coinciding_vertex_group);
+
   /**
    * A structure that allows the transfer of cell data of type @p T from one processor
    * to another. It corresponds to a packed buffer that stores a vector of
index 9d63fcd65c6eb9868d4248a42541127e4f93b544..d7acbeaa7a38d77c11059a90554edf0c69104687 100644 (file)
@@ -154,34 +154,34 @@ namespace parallel
         // store the communicator
         construction_data.comm = comm;
 
+
+        std::map<unsigned int, std::vector<unsigned int>>
+                                             coinciding_vertex_groups;
+        std::map<unsigned int, unsigned int> vertex_to_coinciding_vertex_group;
+
+        GridTools::collect_coinciding_vertices(
+          tria, coinciding_vertex_groups, vertex_to_coinciding_vertex_group);
+
         // helper function, which collects all vertices belonging to a cell
         // (also taking into account periodicity)
         auto add_vertices_of_cell_to_vertices_owned_by_locally_owned_cells =
-          [](TriaIterator<CellAccessor<dim, spacedim>> &cell,
-             std::vector<bool> &vertices_owned_by_locally_owned_cells) {
-            // add vertices belonging to a periodic neighbor
-            for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
-              if (cell->has_periodic_neighbor(f))
-                {
-                  auto face_t = cell->face(f);
-                  auto face_n = cell->periodic_neighbor(f)->face(
-                    cell->periodic_neighbor_face_no(f));
-                  for (unsigned int v = 0;
-                       v < GeometryInfo<dim>::vertices_per_face;
-                       ++v)
-                    {
-                      vertices_owned_by_locally_owned_cells
-                        [face_t->vertex_index(v)] = true;
-                      vertices_owned_by_locally_owned_cells
-                        [face_n->vertex_index(v)] = true;
-                    }
-                }
-
+          [coinciding_vertex_groups, vertex_to_coinciding_vertex_group](
+            TriaIterator<CellAccessor<dim, spacedim>> &cell,
+            std::vector<bool> &vertices_owned_by_locally_owned_cells) {
             // add local vertices
             for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell;
                  ++v)
-              vertices_owned_by_locally_owned_cells[cell->vertex_index(v)] =
-                true;
+              {
+                vertices_owned_by_locally_owned_cells[cell->vertex_index(v)] =
+                  true;
+                const auto coinciding_vertex_group =
+                  vertex_to_coinciding_vertex_group.find(cell->vertex_index(v));
+                if (coinciding_vertex_group !=
+                    vertex_to_coinciding_vertex_group.end())
+                  for (const auto &co_vertex : coinciding_vertex_groups.at(
+                         coinciding_vertex_group->second))
+                    vertices_owned_by_locally_owned_cells[co_vertex] = true;
+              }
           };
 
         // check if multilevel hierarchy should be constructed
index bad673e21e81686fecf05d5f40bb49a9158887c5..dda52d3e1321d73b71040e8bbafc18579d91ac7a 100644 (file)
@@ -364,115 +364,9 @@ namespace parallel
     std::map<unsigned int, std::vector<unsigned int>> coinciding_vertex_groups;
     std::map<unsigned int, unsigned int> vertex_to_coinciding_vertex_group;
 
-    {
-      static const int lookup_table_2d[2][2] =
-        //           flip:
-        {
-          {0, 1}, // false
-          {1, 0}  // true
-        };
-
-      static const int lookup_table_3d[2][2][2][4] =
-        //                   orientation flip  rotation
-        {{{
-            {0, 2, 1, 3}, // false       false false
-            {2, 3, 0, 1}  // false       false true
-          },
-          {
-            {3, 1, 2, 0}, // false       true  false
-            {1, 0, 3, 2}  // false       true  true
-          }},
-         {{
-            {0, 1, 2, 3}, // true        false false
-            {1, 3, 0, 2}  // true        false true
-          },
-          {
-            {3, 2, 1, 0}, // true        true  false
-            {2, 0, 3, 1}  // true        true  true
-          }}};
-
-      // loop over all periodic face pairs
-      for (const auto &pair : this->get_periodic_face_map())
-        {
-          if (pair.first.first->level() != pair.second.first.first->level())
-            continue;
-
-          const auto face_a = pair.first.first->face(pair.first.second);
-          const auto face_b =
-            pair.second.first.first->face(pair.second.first.second);
-          const auto mask = pair.second.second;
-
-          // loop over all vertices on face
-          for (unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_face;
-               ++i)
-            {
-              const bool face_orientation = mask[0];
-              const bool face_flip        = mask[1];
-              const bool face_rotation    = mask[2];
-
-              // find the right local vertex index for the second face
-              unsigned int j = 0;
-              switch (dim)
-                {
-                  case 1:
-                    j = i;
-                    break;
-                  case 2:
-                    j = lookup_table_2d[face_flip][i];
-                    break;
-                  case 3:
-                    j = lookup_table_3d[face_orientation][face_flip]
-                                       [face_rotation][i];
-                    break;
-                  default:
-                    AssertThrow(false, ExcNotImplemented());
-                }
-
-              // get vertex indices and store in map
-              const auto   vertex_a = face_a->vertex_index(i);
-              const auto   vertex_b = face_b->vertex_index(j);
-              unsigned int temp     = std::min(vertex_a, vertex_b);
-
-              auto it_a = vertex_to_coinciding_vertex_group.find(vertex_a);
-              if (it_a != vertex_to_coinciding_vertex_group.end())
-                temp = std::min(temp, it_a->second);
-
-              auto it_b = vertex_to_coinciding_vertex_group.find(vertex_b);
-              if (it_b != vertex_to_coinciding_vertex_group.end())
-                temp = std::min(temp, it_b->second);
-
-              if (it_a != vertex_to_coinciding_vertex_group.end())
-                it_a->second = temp;
-              else
-                vertex_to_coinciding_vertex_group[vertex_a] = temp;
-
-              if (it_b != vertex_to_coinciding_vertex_group.end())
-                it_b->second = temp;
-              else
-                vertex_to_coinciding_vertex_group[vertex_b] = temp;
-            }
-        }
-
-      // 1b) compress map: let vertices point to the coinciding vertex with
-      //     the smallest index
-      for (auto &p : vertex_to_coinciding_vertex_group)
-        {
-          if (p.first == p.second)
-            continue;
-          unsigned int temp = p.second;
-          while (temp != vertex_to_coinciding_vertex_group[temp])
-            temp = vertex_to_coinciding_vertex_group[temp];
-          p.second = temp;
-        }
-
-      // 1c) create a map: smallest index of coinciding index -> all
-      // coinciding indices
-      for (auto p : vertex_to_coinciding_vertex_group)
-        coinciding_vertex_groups[p.second] = {};
-
-      for (auto p : vertex_to_coinciding_vertex_group)
-        coinciding_vertex_groups[p.second].push_back(p.first);
-    }
+    GridTools::collect_coinciding_vertices(*this,
+                                           coinciding_vertex_groups,
+                                           vertex_to_coinciding_vertex_group);
 
     // 2) collect vertices belonging to local cells
     std::vector<bool> vertex_of_own_cell(this->n_vertices(), false);
index efad61ecbd63cb42b6e80bbc1c44c8d96bacfb21..c0cb7035df2e5937404fbf3e6615f7f2354287d3 100644 (file)
@@ -5651,6 +5651,127 @@ namespace GridTools
 
 
 
+  template <int dim, int spacedim>
+  void
+  collect_coinciding_vertices(
+    const Triangulation<dim, spacedim> &               tria,
+    std::map<unsigned int, std::vector<unsigned int>> &coinciding_vertex_groups,
+    std::map<unsigned int, unsigned int> &vertex_to_coinciding_vertex_group)
+  {
+    // 1) determine for each vertex a vertex it concides with and
+    //    put it into a map
+    {
+      static const int lookup_table_2d[2][2] =
+        //           flip:
+        {
+          {0, 1}, // false
+          {1, 0}  // true
+        };
+
+      static const int lookup_table_3d[2][2][2][4] =
+        //                   orientation flip  rotation
+        {{{
+            {0, 2, 1, 3}, // false       false false
+            {2, 3, 0, 1}  // false       false true
+          },
+          {
+            {3, 1, 2, 0}, // false       true  false
+            {1, 0, 3, 2}  // false       true  true
+          }},
+         {{
+            {0, 1, 2, 3}, // true        false false
+            {1, 3, 0, 2}  // true        false true
+          },
+          {
+            {3, 2, 1, 0}, // true        true  false
+            {2, 0, 3, 1}  // true        true  true
+          }}};
+
+      // loop over all periodic face pairs
+      for (const auto &pair : tria.get_periodic_face_map())
+        {
+          if (pair.first.first->level() != pair.second.first.first->level())
+            continue;
+
+          const auto face_a = pair.first.first->face(pair.first.second);
+          const auto face_b =
+            pair.second.first.first->face(pair.second.first.second);
+          const auto mask = pair.second.second;
+
+          // loop over all vertices on face
+          for (unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_face;
+               ++i)
+            {
+              const bool face_orientation = mask[0];
+              const bool face_flip        = mask[1];
+              const bool face_rotation    = mask[2];
+
+              // find the right local vertex index for the second face
+              unsigned int j = 0;
+              switch (dim)
+                {
+                  case 1:
+                    j = i;
+                    break;
+                  case 2:
+                    j = lookup_table_2d[face_flip][i];
+                    break;
+                  case 3:
+                    j = lookup_table_3d[face_orientation][face_flip]
+                                       [face_rotation][i];
+                    break;
+                  default:
+                    AssertThrow(false, ExcNotImplemented());
+                }
+
+              // get vertex indices and store in map
+              const auto   vertex_a = face_a->vertex_index(i);
+              const auto   vertex_b = face_b->vertex_index(j);
+              unsigned int temp     = std::min(vertex_a, vertex_b);
+
+              auto it_a = vertex_to_coinciding_vertex_group.find(vertex_a);
+              if (it_a != vertex_to_coinciding_vertex_group.end())
+                temp = std::min(temp, it_a->second);
+
+              auto it_b = vertex_to_coinciding_vertex_group.find(vertex_b);
+              if (it_b != vertex_to_coinciding_vertex_group.end())
+                temp = std::min(temp, it_b->second);
+
+              if (it_a != vertex_to_coinciding_vertex_group.end())
+                it_a->second = temp;
+              else
+                vertex_to_coinciding_vertex_group[vertex_a] = temp;
+
+              if (it_b != vertex_to_coinciding_vertex_group.end())
+                it_b->second = temp;
+              else
+                vertex_to_coinciding_vertex_group[vertex_b] = temp;
+            }
+        }
+
+      // 2) compress map: let vertices point to the coinciding vertex with
+      //    the smallest index
+      for (auto &p : vertex_to_coinciding_vertex_group)
+        {
+          if (p.first == p.second)
+            continue;
+          unsigned int temp = p.second;
+          while (temp != vertex_to_coinciding_vertex_group[temp])
+            temp = vertex_to_coinciding_vertex_group[temp];
+          p.second = temp;
+        }
+
+      // 3) create a map: smallest index of coinciding index -> all
+      //    coinciding indices
+      for (auto p : vertex_to_coinciding_vertex_group)
+        coinciding_vertex_groups[p.second] = {};
+
+      for (auto p : vertex_to_coinciding_vertex_group)
+        coinciding_vertex_groups[p.second].push_back(p.first);
+    }
+  }
+
+
 } /* namespace GridTools */
 
 
index edf8ab5280d61fe17c4eaf43aaaa59e3739151c2..78ff47ec8bf6a83cb39ecdf88cc64c7360ce5c19 100644 (file)
@@ -404,6 +404,12 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
       regularize_corner_cells(
         Triangulation<deal_II_dimension, deal_II_space_dimension> &,
         double);
+
+      template void
+      collect_coinciding_vertices(
+        const Triangulation<deal_II_dimension, deal_II_space_dimension> &,
+        std::map<unsigned int, std::vector<unsigned int>> &,
+        std::map<unsigned int, unsigned int> &);
     \}
 #endif
   }
index 200f295a708e9360fd81e5a5c514619cd0d344d2..5da7755d2030bea813304187a27b0f3c2410d3f8 100644 (file)
@@ -2293,7 +2293,7 @@ namespace GridTools
       (pairs1.size() > 0 ||
        (dynamic_cast<
           const parallel::fullydistributed::Triangulation<dim, space_dim> *>(
-          &pairs1.begin()->first->get_triangulation()) != nullptr)),
+          &mesh.begin()->get_triangulation()) != nullptr)),
       ExcMessage("No new periodic face pairs have been found. "
                  "Are you sure that you've selected the correct boundary "
                  "id's and that the coarsest level mesh is colorized?"));
index 64fd34e9377b6d7bf2542515070d4661ae6d5f54..673bfcde7a2b3dd0ac8df8114c93a38f6c5f4b00 100644 (file)
@@ -56,8 +56,25 @@ test(const int n_refinements, const int n_subdivisions, MPI_Comm comm)
         else if (std::fabs(cell->face(face_number)->center()(0) - right) <
                  1e-12)
           cell->face(face_number)->set_all_boundary_ids(2);
-
-    GridTools::collect_periodic_faces(tria, 1, 2, 0, periodic_faces);
+        else if (dim >= 2 &&
+                 std::fabs(cell->face(face_number)->center()(1) - left) < 1e-12)
+          cell->face(face_number)->set_all_boundary_ids(3);
+        else if (dim >= 2 && std::fabs(cell->face(face_number)->center()(1) -
+                                       right) < 1e-12)
+          cell->face(face_number)->set_all_boundary_ids(4);
+        else if (dim >= 3 &&
+                 std::fabs(cell->face(face_number)->center()(2) - left) < 1e-12)
+          cell->face(face_number)->set_all_boundary_ids(5);
+        else if (dim >= 3 && std::fabs(cell->face(face_number)->center()(2) -
+                                       right) < 1e-12)
+          cell->face(face_number)->set_all_boundary_ids(6);
+
+    if (dim >= 1)
+      GridTools::collect_periodic_faces(tria, 1, 2, 0, periodic_faces);
+    if (dim >= 2)
+      GridTools::collect_periodic_faces(tria, 3, 4, 1, periodic_faces);
+    if (dim >= 3)
+      GridTools::collect_periodic_faces(tria, 5, 6, 2, periodic_faces);
 
     tria.add_periodicity(periodic_faces);
   };
@@ -120,7 +137,7 @@ main(int argc, char *argv[])
   }
 
   {
-    deallog.push("2d");
+    deallog.push("3d");
     const int n_refinements  = 3;
     const int n_subdivisions = 4;
     test<3>(n_refinements, n_subdivisions, comm);
index a1c457fa4f720b5be417b6c26e8ba8a2e020df5f..ed97205217d80fc140660604e4c2382209853614 100644 (file)
@@ -13,10 +13,10 @@ DEAL:0:2d::
 DEAL:0:2d::n_dofs:                             16641
 DEAL:0:2d::n_locally_owned_dofs:               16641
 DEAL:0:2d::
-DEAL:0:2d::n_levels:                  1
-DEAL:0:2d::n_cells:                   32768
-DEAL:0:2d::n_active_cells:            32768
-DEAL:0:2d::
-DEAL:0:2d::n_dofs:                             274625
-DEAL:0:2d::n_locally_owned_dofs:               274625
-DEAL:0:2d::
+DEAL:0:3d::n_levels:                  1
+DEAL:0:3d::n_cells:                   32768
+DEAL:0:3d::n_active_cells:            32768
+DEAL:0:3d::
+DEAL:0:3d::n_dofs:                             274625
+DEAL:0:3d::n_locally_owned_dofs:               274625
+DEAL:0:3d::
diff --git a/tests/fullydistributed_grids/copy_serial_tria_04.mpirun=4.output b/tests/fullydistributed_grids/copy_serial_tria_04.mpirun=4.output
new file mode 100644 (file)
index 0000000..0046a78
--- /dev/null
@@ -0,0 +1,91 @@
+
+DEAL:0:1d::n_levels:                  1
+DEAL:0:1d::n_cells:                   34
+DEAL:0:1d::n_active_cells:            34
+DEAL:0:1d::
+DEAL:0:1d::n_dofs:                             257
+DEAL:0:1d::n_locally_owned_dofs:               65
+DEAL:0:1d::
+DEAL:0:2d::n_levels:                  1
+DEAL:0:2d::n_cells:                   1156
+DEAL:0:2d::n_active_cells:            1156
+DEAL:0:2d::
+DEAL:0:2d::n_dofs:                             16641
+DEAL:0:2d::n_locally_owned_dofs:               4225
+DEAL:0:2d::
+DEAL:0:3d::n_levels:                  1
+DEAL:0:3d::n_cells:                   9792
+DEAL:0:3d::n_active_cells:            9792
+DEAL:0:3d::
+DEAL:0:3d::n_dofs:                             274625
+DEAL:0:3d::n_locally_owned_dofs:               70785
+DEAL:0:3d::
+
+DEAL:1:1d::n_levels:                  1
+DEAL:1:1d::n_cells:                   34
+DEAL:1:1d::n_active_cells:            34
+DEAL:1:1d::
+DEAL:1:1d::n_dofs:                             257
+DEAL:1:1d::n_locally_owned_dofs:               64
+DEAL:1:1d::
+DEAL:1:2d::n_levels:                  1
+DEAL:1:2d::n_cells:                   1156
+DEAL:1:2d::n_active_cells:            1156
+DEAL:1:2d::
+DEAL:1:2d::n_dofs:                             16641
+DEAL:1:2d::n_locally_owned_dofs:               4160
+DEAL:1:2d::
+DEAL:1:3d::n_levels:                  1
+DEAL:1:3d::n_cells:                   9792
+DEAL:1:3d::n_active_cells:            9792
+DEAL:1:3d::
+DEAL:1:3d::n_dofs:                             274625
+DEAL:1:3d::n_locally_owned_dofs:               68640
+DEAL:1:3d::
+
+
+DEAL:2:1d::n_levels:                  1
+DEAL:2:1d::n_cells:                   34
+DEAL:2:1d::n_active_cells:            34
+DEAL:2:1d::
+DEAL:2:1d::n_dofs:                             257
+DEAL:2:1d::n_locally_owned_dofs:               64
+DEAL:2:1d::
+DEAL:2:2d::n_levels:                  1
+DEAL:2:2d::n_cells:                   1156
+DEAL:2:2d::n_active_cells:            1156
+DEAL:2:2d::
+DEAL:2:2d::n_dofs:                             16641
+DEAL:2:2d::n_locally_owned_dofs:               4160
+DEAL:2:2d::
+DEAL:2:3d::n_levels:                  1
+DEAL:2:3d::n_cells:                   9792
+DEAL:2:3d::n_active_cells:            9792
+DEAL:2:3d::
+DEAL:2:3d::n_dofs:                             274625
+DEAL:2:3d::n_locally_owned_dofs:               68640
+DEAL:2:3d::
+
+
+DEAL:3:1d::n_levels:                  1
+DEAL:3:1d::n_cells:                   34
+DEAL:3:1d::n_active_cells:            34
+DEAL:3:1d::
+DEAL:3:1d::n_dofs:                             257
+DEAL:3:1d::n_locally_owned_dofs:               64
+DEAL:3:1d::
+DEAL:3:2d::n_levels:                  1
+DEAL:3:2d::n_cells:                   1156
+DEAL:3:2d::n_active_cells:            1156
+DEAL:3:2d::
+DEAL:3:2d::n_dofs:                             16641
+DEAL:3:2d::n_locally_owned_dofs:               4096
+DEAL:3:2d::
+DEAL:3:3d::n_levels:                  1
+DEAL:3:3d::n_cells:                   9792
+DEAL:3:3d::n_active_cells:            9792
+DEAL:3:3d::
+DEAL:3:3d::n_dofs:                             274625
+DEAL:3:3d::n_locally_owned_dofs:               66560
+DEAL:3:3d::
+
index 8c70249a7b7bae4dfa7766b3d82595c26af4b5b4..f963884650d58b192d29362f01bf29d04ff2534f 100644 (file)
@@ -7,19 +7,19 @@ DEAL:0:1d::n_dofs:                             257
 DEAL:0:1d::n_locally_owned_dofs:               53
 DEAL:0:1d::
 DEAL:0:2d::n_levels:                  1
-DEAL:0:2d::n_cells:                   918
-DEAL:0:2d::n_active_cells:            918
+DEAL:0:2d::n_cells:                   952
+DEAL:0:2d::n_active_cells:            952
 DEAL:0:2d::
 DEAL:0:2d::n_dofs:                             16641
 DEAL:0:2d::n_locally_owned_dofs:               3409
 DEAL:0:2d::
-DEAL:0:2d::n_levels:                  1
-DEAL:0:2d::n_cells:                   7816
-DEAL:0:2d::n_active_cells:            7816
-DEAL:0:2d::
-DEAL:0:2d::n_dofs:                             274625
-DEAL:0:2d::n_locally_owned_dofs:               57409
-DEAL:0:2d::
+DEAL:0:3d::n_levels:                  1
+DEAL:0:3d::n_cells:                   8312
+DEAL:0:3d::n_active_cells:            8312
+DEAL:0:3d::
+DEAL:0:3d::n_dofs:                             274625
+DEAL:0:3d::n_locally_owned_dofs:               57409
+DEAL:0:3d::
 
 DEAL:1:1d::n_levels:                  1
 DEAL:1:1d::n_cells:                   28
@@ -29,19 +29,19 @@ DEAL:1:1d::n_dofs:                             257
 DEAL:1:1d::n_locally_owned_dofs:               52
 DEAL:1:1d::
 DEAL:1:2d::n_levels:                  1
-DEAL:1:2d::n_cells:                   950
-DEAL:1:2d::n_active_cells:            950
+DEAL:1:2d::n_cells:                   984
+DEAL:1:2d::n_active_cells:            984
 DEAL:1:2d::
 DEAL:1:2d::n_dofs:                             16641
 DEAL:1:2d::n_locally_owned_dofs:               3344
 DEAL:1:2d::
-DEAL:1:2d::n_levels:                  1
-DEAL:1:2d::n_cells:                   8842
-DEAL:1:2d::n_active_cells:            8842
-DEAL:1:2d::
-DEAL:1:2d::n_dofs:                             274625
-DEAL:1:2d::n_locally_owned_dofs:               55264
-DEAL:1:2d::
+DEAL:1:3d::n_levels:                  1
+DEAL:1:3d::n_cells:                   9208
+DEAL:1:3d::n_active_cells:            9208
+DEAL:1:3d::
+DEAL:1:3d::n_dofs:                             274625
+DEAL:1:3d::n_locally_owned_dofs:               55264
+DEAL:1:3d::
 
 
 DEAL:2:1d::n_levels:                  1
@@ -58,13 +58,13 @@ DEAL:2:2d::
 DEAL:2:2d::n_dofs:                             16641
 DEAL:2:2d::n_locally_owned_dofs:               3296
 DEAL:2:2d::
-DEAL:2:2d::n_levels:                  1
-DEAL:2:2d::n_cells:                   8916
-DEAL:2:2d::n_active_cells:            8916
-DEAL:2:2d::
-DEAL:2:2d::n_dofs:                             274625
-DEAL:2:2d::n_locally_owned_dofs:               54944
-DEAL:2:2d::
+DEAL:2:3d::n_levels:                  1
+DEAL:2:3d::n_cells:                   9512
+DEAL:2:3d::n_active_cells:            9512
+DEAL:2:3d::
+DEAL:2:3d::n_dofs:                             274625
+DEAL:2:3d::n_locally_owned_dofs:               54944
+DEAL:2:3d::
 
 
 DEAL:3:1d::n_levels:                  1
@@ -75,19 +75,19 @@ DEAL:3:1d::n_dofs:                             257
 DEAL:3:1d::n_locally_owned_dofs:               52
 DEAL:3:1d::
 DEAL:3:2d::n_levels:                  1
-DEAL:3:2d::n_cells:                   950
-DEAL:3:2d::n_active_cells:            950
+DEAL:3:2d::n_cells:                   984
+DEAL:3:2d::n_active_cells:            984
 DEAL:3:2d::
 DEAL:3:2d::n_dofs:                             16641
 DEAL:3:2d::n_locally_owned_dofs:               3312
 DEAL:3:2d::
-DEAL:3:2d::n_levels:                  1
-DEAL:3:2d::n_cells:                   8842
-DEAL:3:2d::n_active_cells:            8842
-DEAL:3:2d::
-DEAL:3:2d::n_dofs:                             274625
-DEAL:3:2d::n_locally_owned_dofs:               54080
-DEAL:3:2d::
+DEAL:3:3d::n_levels:                  1
+DEAL:3:3d::n_cells:                   9208
+DEAL:3:3d::n_active_cells:            9208
+DEAL:3:3d::
+DEAL:3:3d::n_dofs:                             274625
+DEAL:3:3d::n_locally_owned_dofs:               54080
+DEAL:3:3d::
 
 
 DEAL:4:1d::n_levels:                  1
@@ -98,17 +98,17 @@ DEAL:4:1d::n_dofs:                             257
 DEAL:4:1d::n_locally_owned_dofs:               52
 DEAL:4:1d::
 DEAL:4:2d::n_levels:                  1
-DEAL:4:2d::n_cells:                   918
-DEAL:4:2d::n_active_cells:            918
+DEAL:4:2d::n_cells:                   952
+DEAL:4:2d::n_active_cells:            952
 DEAL:4:2d::
 DEAL:4:2d::n_dofs:                             16641
 DEAL:4:2d::n_locally_owned_dofs:               3280
 DEAL:4:2d::
-DEAL:4:2d::n_levels:                  1
-DEAL:4:2d::n_cells:                   7816
-DEAL:4:2d::n_active_cells:            7816
-DEAL:4:2d::
-DEAL:4:2d::n_dofs:                             274625
-DEAL:4:2d::n_locally_owned_dofs:               52928
-DEAL:4:2d::
+DEAL:4:3d::n_levels:                  1
+DEAL:4:3d::n_cells:                   8312
+DEAL:4:3d::n_active_cells:            8312
+DEAL:4:3d::
+DEAL:4:3d::n_dofs:                             274625
+DEAL:4:3d::n_locally_owned_dofs:               52928
+DEAL:4:3d::
 

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.