]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Treat cells to be h-coarsened as their parent in `hp::Refinement::limit_p_level_diffe... 11943/head
authorMarc Fehling <mafehling.git@gmail.com>
Thu, 18 Mar 2021 20:19:40 +0000 (14:19 -0600)
committerMarc Fehling <mafehling.git@gmail.com>
Mon, 22 Mar 2021 19:39:53 +0000 (13:39 -0600)
22 files changed:
include/deal.II/hp/refinement.h
source/hp/refinement.cc
tests/hp/limit_p_level_difference_01.cc [moved from tests/hp/prepare_coarsening_and_refinement_01.cc with 100% similarity]
tests/hp/limit_p_level_difference_01.output [moved from tests/hp/prepare_coarsening_and_refinement_01.output with 100% similarity]
tests/hp/limit_p_level_difference_02.cc [moved from tests/hp/prepare_coarsening_and_refinement_02.cc with 100% similarity]
tests/hp/limit_p_level_difference_02.output [moved from tests/hp/prepare_coarsening_and_refinement_02.output with 100% similarity]
tests/hp/limit_p_level_difference_03.cc [moved from tests/hp/prepare_coarsening_and_refinement_03.cc with 100% similarity]
tests/hp/limit_p_level_difference_03.output [moved from tests/hp/prepare_coarsening_and_refinement_03.output with 100% similarity]
tests/hp/limit_p_level_difference_04.cc [new file with mode: 0644]
tests/hp/limit_p_level_difference_04.output [new file with mode: 0644]
tests/mpi/limit_p_level_difference_01.cc [moved from tests/mpi/prepare_coarsening_and_refinement_01.cc with 100% similarity]
tests/mpi/limit_p_level_difference_01.with_p4est=true.mpirun=1.output [moved from tests/mpi/prepare_coarsening_and_refinement_01.with_p4est=true.mpirun=1.output with 100% similarity]
tests/mpi/limit_p_level_difference_01.with_p4est=true.mpirun=4.output [moved from tests/mpi/prepare_coarsening_and_refinement_01.with_p4est=true.mpirun=4.output with 100% similarity]
tests/mpi/limit_p_level_difference_02.cc [moved from tests/mpi/prepare_coarsening_and_refinement_02.cc with 100% similarity]
tests/mpi/limit_p_level_difference_02.with_p4est=true.mpirun=1.output [moved from tests/mpi/prepare_coarsening_and_refinement_02.with_p4est=true.mpirun=1.output with 100% similarity]
tests/mpi/limit_p_level_difference_02.with_p4est=true.mpirun=4.output [moved from tests/mpi/prepare_coarsening_and_refinement_02.with_p4est=true.mpirun=4.output with 100% similarity]
tests/mpi/limit_p_level_difference_03.cc [moved from tests/mpi/prepare_coarsening_and_refinement_03.cc with 100% similarity]
tests/mpi/limit_p_level_difference_03.with_p4est=true.mpirun=1.output [moved from tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=1.output with 100% similarity]
tests/mpi/limit_p_level_difference_03.with_p4est=true.mpirun=4.output [moved from tests/mpi/prepare_coarsening_and_refinement_03.with_p4est=true.mpirun=4.output with 100% similarity]
tests/mpi/limit_p_level_difference_04.cc [new file with mode: 0644]
tests/mpi/limit_p_level_difference_04.with_p4est=true.mpirun=1.output [new file with mode: 0644]
tests/mpi/limit_p_level_difference_04.with_p4est=true.mpirun=4.output [new file with mode: 0644]

index b674183cd0810a898c645684116b48bd860cb440..72c35f3f68966499e0750c6a11c106c44cc224dc 100644 (file)
@@ -677,6 +677,13 @@ namespace hp
      * to call this function, nor will it be automatically invoked in any part
      * of the library (contrary to its Triangulation counterpart).
      *
+     * On cells that will be h-coarsened, we enforce the difference criterion as
+     * if it is already a parent cell. That means, we set the level of all
+     * siblings to the highest one among them. In that case, all sibling cells
+     * need to have the h-coarsenening flags set terminally via
+     * Triangulation::prepare_coarsening_and_refinement() beforehand. Otherwise
+     * an assertion will be triggered.
+     *
      * Returns whether any future FE indices have been changed by this function.
      */
     template <int dim, int spacedim>
index edc5cb179488e0f9276e716633688aadc75caf43..def4bb4a8c264c711c841ffb67c79eee0ad29458 100644 (file)
@@ -843,6 +843,99 @@ namespace hp
       // - always raise levels to match criterion, never lower them
       // - exchange level indices on ghost cells
 
+      // Function that updates the level of neighbor to fulfill difference
+      // criterion, and returns whether it was changed.
+      const auto update_neighbor_level =
+        [&future_levels, max_difference](const auto &     neighbor,
+                                         const level_type cell_level) -> bool {
+        Assert(neighbor->is_active(), ExcInternalError());
+        // We only care about locally owned neighbors. If neighbor is a ghost
+        // cell, its future FE index will be updated on the owning process and
+        // communicated at the next loop iteration.
+        if (neighbor->is_locally_owned())
+          {
+            const level_type neighbor_level = static_cast<level_type>(
+              future_levels[neighbor->global_active_cell_index()]);
+
+            // ignore neighbors that are not part of the hierarchy
+            if (neighbor_level == invalid_level)
+              return false;
+
+            if ((cell_level - max_difference) > neighbor_level)
+              {
+                future_levels[neighbor->global_active_cell_index()] =
+                  cell_level - max_difference;
+
+                return true;
+              }
+          }
+
+        return false;
+      };
+
+      // For cells to be h-coarsened, we need to determine a future FE for the
+      // parent cell, which will be the dominated FE among all children
+      // However, if we want to enforce the max_difference criterion on all
+      // cells on the updated mesh, we will need to simulate the updated mesh on
+      // the current mesh.
+      //
+      // As we are working on p-levels, we will set all siblings that will be
+      // coarsened to the highest p-level among them. The parent cell will be
+      // assigned exactly this level in form of the corresponding FE index in
+      // the adaptation process in
+      // Triangulation::execute_coarsening_and_refinement().
+      //
+      // This function takes a cell and sets all its siblings to the highest
+      // p-level among them. Returns whether any future levels have been
+      // changed.
+      const auto prepare_level_for_parent = [&](const auto &neighbor) -> bool {
+        Assert(neighbor->is_active(), ExcInternalError());
+        if (neighbor->coarsen_flag_set() && neighbor->is_locally_owned())
+          {
+            const auto parent = neighbor->parent();
+
+            std::vector<unsigned int> future_levels_children;
+            future_levels_children.reserve(parent->n_children());
+            for (const auto &child : parent->child_iterators())
+              {
+                Assert(child->is_active() && child->coarsen_flag_set(),
+                       (typename dealii::Triangulation<dim, spacedim>::
+                          ExcInconsistentCoarseningFlags()));
+
+                const level_type child_level = static_cast<level_type>(
+                  future_levels[child->global_active_cell_index()]);
+                Assert(child_level != invalid_level,
+                       ExcMessage(
+                         "The FiniteElement on one of the siblings of "
+                         "a cell you are trying to coarsen is not part "
+                         "of the registered p-adaptation hierarchy."));
+                future_levels_children.push_back(child_level);
+              }
+            Assert(!future_levels_children.empty(), ExcInternalError());
+
+            const unsigned int max_level_children =
+              *std::max_element(future_levels_children.begin(),
+                                future_levels_children.end());
+
+            bool children_changed = false;
+            for (const auto &child : parent->child_iterators())
+              // We only care about locally owned children. If child is a ghost
+              // cell, its future FE index will be updated on the owning process
+              // and communicated at the next loop iteration.
+              if (child->is_locally_owned() &&
+                  future_levels[child->global_active_cell_index()] !=
+                    max_level_children)
+                {
+                  future_levels[child->global_active_cell_index()] =
+                    max_level_children;
+                  children_changed = true;
+                }
+            return children_changed;
+          }
+
+        return false;
+      };
+
       bool levels_changed = false;
       bool levels_changed_in_cycle;
       do
@@ -878,64 +971,22 @@ namespace hp
                               const auto neighbor =
                                 cell->neighbor_child_on_subface(f, sf);
 
-                              // We only care about locally owned neighbors. If
-                              // neighbor is a ghost cell, its future FE index
-                              // will be updated on the owning process and
-                              // communicated at the next loop iteration.
-                              if (neighbor->is_locally_owned())
-                                {
-                                  const level_type neighbor_level =
-                                    static_cast<level_type>(
-                                      future_levels
-                                        [neighbor->global_active_cell_index()]);
-
-                                  // ignore neighbors that are not part of the
-                                  // hierarchy
-                                  if (neighbor_level == invalid_level)
-                                    continue;
-
-                                  if ((cell_level - max_difference) >
-                                      neighbor_level)
-                                    {
-                                      future_levels
-                                        [neighbor->global_active_cell_index()] =
-                                          cell_level - max_difference;
-
-                                      levels_changed_in_cycle = true;
-                                    }
-                                }
+                              levels_changed_in_cycle |=
+                                update_neighbor_level(neighbor, cell_level);
+
+                              levels_changed_in_cycle |=
+                                prepare_level_for_parent(neighbor);
                             }
                         }
                       else
                         {
                           const auto neighbor = cell->neighbor(f);
 
-                          // We only care about locally owned neighbors. If
-                          // neighbor is a ghost cell, its future FE index will
-                          // be updated on the owning process and communicated
-                          // at the next loop iteration.
-                          if (neighbor->is_locally_owned())
-                            {
-                              const level_type neighbor_level =
-                                static_cast<level_type>(
-                                  future_levels
-                                    [neighbor->global_active_cell_index()]);
-
-                              // ignore neighbors that are not part of the
-                              // hierarchy
-                              if (neighbor_level == invalid_level)
-                                continue;
-
-                              if ((cell_level - max_difference) >
-                                  neighbor_level)
-                                {
-                                  future_levels
-                                    [neighbor->global_active_cell_index()] =
-                                      cell_level - max_difference;
-
-                                  levels_changed_in_cycle = true;
-                                }
-                            }
+                          levels_changed_in_cycle |=
+                            update_neighbor_level(neighbor, cell_level);
+
+                          levels_changed_in_cycle |=
+                            prepare_level_for_parent(neighbor);
                         }
                     }
               }
diff --git a/tests/hp/limit_p_level_difference_04.cc b/tests/hp/limit_p_level_difference_04.cc
new file mode 100644 (file)
index 0000000..82058de
--- /dev/null
@@ -0,0 +1,116 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// verify restrictions on level differences imposed by
+// DoFHandler::prepare_coarsening_and_refinement()
+// on h-coarsened grids
+
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/hp/refinement.h>
+
+#include "../tests.h"
+
+#include "../test_grids.h"
+
+
+template <int dim>
+void
+test(const unsigned int max_difference)
+{
+  Assert(max_difference > 0, ExcInternalError());
+
+  // setup FE collection
+  const unsigned int fes_size = 2 * max_difference + 2;
+
+  hp::FECollection<dim> fes;
+  while (fes.size() < fes_size)
+    fes.push_back(FE_Q<dim>(1));
+
+  const unsigned int contains_fe_index = 0;
+  const auto         sequence = fes.get_hierarchy_sequence(contains_fe_index);
+
+  // set up line grid
+  // - refine central cell and flag them for coarsening
+  // - assign highest p-level to leftmost cell
+  //
+  // +-------+---+---+-------+
+  // |       | c | c |       |
+  // |       +---+---+       |
+  // |       | c | c |       |
+  // +-------+---+---+-------+
+  //
+  // after prepare_coarsening_and_refinement(), the p-levels should propagate
+  // through the central cells as if they were already coarsened
+
+  Triangulation<dim> tria;
+  TestGrids::hyper_line(tria, 3);
+
+  // refine the central cell
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->center()[0] > 1. && cell->center()[0] < 2.)
+      cell->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+
+  // now flag these cells for coarsening
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->center()[0] > 1. && cell->center()[0] < 2.)
+      cell->set_coarsen_flag();
+
+  DoFHandler<dim> dofh(tria);
+  for (const auto &cell : dofh.active_cell_iterators())
+    if (cell->center()[0] < 1.)
+      cell->set_active_fe_index(sequence.back());
+  dofh.distribute_dofs(fes);
+
+  const bool fe_indices_changed =
+    hp::Refinement::limit_p_level_difference(dofh,
+                                             max_difference,
+                                             contains_fe_index);
+  (void)fe_indices_changed;
+  Assert(fe_indices_changed, ExcInternalError());
+
+  deallog << "future FE indices before adaptation:" << std::endl;
+  for (const auto &cell : dofh.active_cell_iterators())
+    deallog << " " << cell->id().to_string() << " " << cell->future_fe_index()
+            << std::endl;
+
+  tria.execute_coarsening_and_refinement();
+
+  deallog << "active FE indices after adaptation:" << std::endl;
+  for (const auto &cell : dofh.active_cell_iterators())
+    deallog << " " << cell->id().to_string() << " " << cell->active_fe_index()
+            << std::endl;
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main()
+{
+  initlog();
+
+  test<2>(1);
+  test<2>(2);
+  test<2>(3);
+}
diff --git a/tests/hp/limit_p_level_difference_04.output b/tests/hp/limit_p_level_difference_04.output
new file mode 100644 (file)
index 0000000..6ed8f09
--- /dev/null
@@ -0,0 +1,37 @@
+
+DEAL::future FE indices before adaptation:
+DEAL:: 0_0: 3
+DEAL:: 2_0: 1
+DEAL:: 1_1:0 2
+DEAL:: 1_1:1 2
+DEAL:: 1_1:2 2
+DEAL:: 1_1:3 2
+DEAL::active FE indices after adaptation:
+DEAL:: 0_0: 3
+DEAL:: 1_0: 2
+DEAL:: 2_0: 1
+DEAL::OK
+DEAL::future FE indices before adaptation:
+DEAL:: 0_0: 5
+DEAL:: 2_0: 1
+DEAL:: 1_1:0 3
+DEAL:: 1_1:1 3
+DEAL:: 1_1:2 3
+DEAL:: 1_1:3 3
+DEAL::active FE indices after adaptation:
+DEAL:: 0_0: 5
+DEAL:: 1_0: 3
+DEAL:: 2_0: 1
+DEAL::OK
+DEAL::future FE indices before adaptation:
+DEAL:: 0_0: 7
+DEAL:: 2_0: 1
+DEAL:: 1_1:0 4
+DEAL:: 1_1:1 4
+DEAL:: 1_1:2 4
+DEAL:: 1_1:3 4
+DEAL::active FE indices after adaptation:
+DEAL:: 0_0: 7
+DEAL:: 1_0: 4
+DEAL:: 2_0: 1
+DEAL::OK
diff --git a/tests/mpi/limit_p_level_difference_04.cc b/tests/mpi/limit_p_level_difference_04.cc
new file mode 100644 (file)
index 0000000..28b3753
--- /dev/null
@@ -0,0 +1,159 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// verify restrictions on level differences imposed by
+// DoFHandler::prepare_coarsening_and_refinement()
+// on h-coarsened grids
+
+
+#include <deal.II/distributed/shared_tria.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/distributed/tria_base.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/hp/refinement.h>
+
+#include "../tests.h"
+
+#include "../test_grids.h"
+
+
+template <int dim>
+void
+test(parallel::TriangulationBase<dim> &tria, const unsigned int max_difference)
+{
+  Assert(tria.n_levels() == 0, ExcInternalError());
+  Assert(max_difference > 0, ExcInternalError());
+
+  // setup FE collection
+  const unsigned int fes_size = 2 * max_difference + 2;
+
+  hp::FECollection<dim> fes;
+  while (fes.size() < fes_size)
+    fes.push_back(FE_Q<dim>(1));
+
+  const unsigned int contains_fe_index = 0;
+  const auto         sequence = fes.get_hierarchy_sequence(contains_fe_index);
+
+  // set up line grid
+  // - refine central cell and flag them for coarsening
+  // - assign highest p-level to leftmost cell
+  //
+  // +-------+---+---+-------+
+  // |       | c | c |       |
+  // |       +---+---+       |
+  // |       | c | c |       |
+  // +-------+---+---+-------+
+  //
+  // after prepare_coarsening_and_refinement(), the p-levels should propagate
+  // through the central cells as if they were already coarsened
+
+  TestGrids::hyper_line(tria, 3);
+
+  // refine the central cell
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->center()[0] > 1. && cell->center()[0] < 2.)
+      cell->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+
+  // now flag these cells for coarsening
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->center()[0] > 1. && cell->center()[0] < 2.)
+      cell->set_coarsen_flag();
+
+  DoFHandler<dim> dofh(tria);
+  for (const auto &cell : dofh.active_cell_iterators())
+    if (cell->is_locally_owned() && cell->center()[0] < 1.)
+      cell->set_active_fe_index(sequence.back());
+  dofh.distribute_dofs(fes);
+
+  const bool fe_indices_changed =
+    hp::Refinement::limit_p_level_difference(dofh,
+                                             max_difference,
+                                             contains_fe_index);
+  (void)fe_indices_changed;
+  Assert(fe_indices_changed, ExcInternalError());
+
+  deallog << "future FE indices before adaptation:" << std::endl;
+  for (const auto &cell : dofh.active_cell_iterators())
+    if (cell->is_locally_owned())
+      deallog << " " << cell->id().to_string() << " " << cell->future_fe_index()
+              << std::endl;
+
+  tria.execute_coarsening_and_refinement();
+
+  deallog << "active FE indices after adaptation:" << std::endl;
+  for (const auto &cell : dofh.active_cell_iterators())
+    if (cell->is_locally_owned())
+      deallog << " " << cell->id().to_string() << " " << cell->active_fe_index()
+              << std::endl;
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    log;
+
+  constexpr const unsigned int dim = 2;
+
+  // Known bug:
+  //   Collecting FE indices for parent cells on children fails for p::s::T
+  //   at the moment whenever siblings belong to the different processors.
+  //
+  //  deallog << "parallel::shared::Triangulation" << std::endl;
+  //  {
+  //    parallel::shared::Triangulation<dim> tria(MPI_COMM_WORLD);
+  //
+  //    test<dim>(tria, 1);
+  //    tria.clear();
+  //    test<dim>(tria, 2);
+  //    tria.clear();
+  //    test<dim>(tria, 3);
+  //  }
+  //
+  //  deallog << "parallel::shared::Triangulation with artificial cells"
+  //          << std::endl;
+  //  {
+  //    parallel::shared::Triangulation<dim> tria(MPI_COMM_WORLD,
+  //                                              Triangulation<dim>::none,
+  //                                              /*allow_artificial_cells=*/true);
+  //
+  //    test<dim>(tria, 1);
+  //    tria.clear();
+  //    test<dim>(tria, 2);
+  //    tria.clear();
+  //    test<dim>(tria, 3);
+  //  }
+
+  deallog << "parallel::distributed::Triangulation" << std::endl;
+  {
+    parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+
+    test<dim>(tria, 1);
+    tria.clear();
+    test<dim>(tria, 2);
+    tria.clear();
+    test<dim>(tria, 3);
+  }
+}
diff --git a/tests/mpi/limit_p_level_difference_04.with_p4est=true.mpirun=1.output b/tests/mpi/limit_p_level_difference_04.with_p4est=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..ce0582b
--- /dev/null
@@ -0,0 +1,38 @@
+
+DEAL:0::parallel::distributed::Triangulation
+DEAL:0::future FE indices before adaptation:
+DEAL:0:: 0_0: 3
+DEAL:0:: 2_0: 1
+DEAL:0:: 1_1:0 2
+DEAL:0:: 1_1:1 2
+DEAL:0:: 1_1:2 2
+DEAL:0:: 1_1:3 2
+DEAL:0::active FE indices after adaptation:
+DEAL:0:: 0_0: 3
+DEAL:0:: 1_0: 2
+DEAL:0:: 2_0: 1
+DEAL:0::OK
+DEAL:0::future FE indices before adaptation:
+DEAL:0:: 0_0: 5
+DEAL:0:: 2_0: 1
+DEAL:0:: 1_1:0 3
+DEAL:0:: 1_1:1 3
+DEAL:0:: 1_1:2 3
+DEAL:0:: 1_1:3 3
+DEAL:0::active FE indices after adaptation:
+DEAL:0:: 0_0: 5
+DEAL:0:: 1_0: 3
+DEAL:0:: 2_0: 1
+DEAL:0::OK
+DEAL:0::future FE indices before adaptation:
+DEAL:0:: 0_0: 7
+DEAL:0:: 2_0: 1
+DEAL:0:: 1_1:0 4
+DEAL:0:: 1_1:1 4
+DEAL:0:: 1_1:2 4
+DEAL:0:: 1_1:3 4
+DEAL:0::active FE indices after adaptation:
+DEAL:0:: 0_0: 7
+DEAL:0:: 1_0: 4
+DEAL:0:: 2_0: 1
+DEAL:0::OK
diff --git a/tests/mpi/limit_p_level_difference_04.with_p4est=true.mpirun=4.output b/tests/mpi/limit_p_level_difference_04.with_p4est=true.mpirun=4.output
new file mode 100644 (file)
index 0000000..12eec57
--- /dev/null
@@ -0,0 +1,74 @@
+
+DEAL:0::parallel::distributed::Triangulation
+DEAL:0::future FE indices before adaptation:
+DEAL:0:: 0_0: 3
+DEAL:0::active FE indices after adaptation:
+DEAL:0::OK
+DEAL:0::future FE indices before adaptation:
+DEAL:0:: 0_0: 5
+DEAL:0::active FE indices after adaptation:
+DEAL:0::OK
+DEAL:0::future FE indices before adaptation:
+DEAL:0:: 0_0: 7
+DEAL:0::active FE indices after adaptation:
+DEAL:0::OK
+
+DEAL:1::parallel::distributed::Triangulation
+DEAL:1::future FE indices before adaptation:
+DEAL:1:: 1_1:0 2
+DEAL:1:: 1_1:1 2
+DEAL:1:: 1_1:2 2
+DEAL:1:: 1_1:3 2
+DEAL:1::active FE indices after adaptation:
+DEAL:1:: 0_0: 3
+DEAL:1::OK
+DEAL:1::future FE indices before adaptation:
+DEAL:1:: 1_1:0 3
+DEAL:1:: 1_1:1 3
+DEAL:1:: 1_1:2 3
+DEAL:1:: 1_1:3 3
+DEAL:1::active FE indices after adaptation:
+DEAL:1:: 0_0: 5
+DEAL:1::OK
+DEAL:1::future FE indices before adaptation:
+DEAL:1:: 1_1:0 4
+DEAL:1:: 1_1:1 4
+DEAL:1:: 1_1:2 4
+DEAL:1:: 1_1:3 4
+DEAL:1::active FE indices after adaptation:
+DEAL:1:: 0_0: 7
+DEAL:1::OK
+
+
+DEAL:2::parallel::distributed::Triangulation
+DEAL:2::future FE indices before adaptation:
+DEAL:2::active FE indices after adaptation:
+DEAL:2:: 1_0: 2
+DEAL:2::OK
+DEAL:2::future FE indices before adaptation:
+DEAL:2::active FE indices after adaptation:
+DEAL:2:: 1_0: 3
+DEAL:2::OK
+DEAL:2::future FE indices before adaptation:
+DEAL:2::active FE indices after adaptation:
+DEAL:2:: 1_0: 4
+DEAL:2::OK
+
+
+DEAL:3::parallel::distributed::Triangulation
+DEAL:3::future FE indices before adaptation:
+DEAL:3:: 2_0: 1
+DEAL:3::active FE indices after adaptation:
+DEAL:3:: 2_0: 1
+DEAL:3::OK
+DEAL:3::future FE indices before adaptation:
+DEAL:3:: 2_0: 1
+DEAL:3::active FE indices after adaptation:
+DEAL:3:: 2_0: 1
+DEAL:3::OK
+DEAL:3::future FE indices before adaptation:
+DEAL:3:: 2_0: 1
+DEAL:3::active FE indices after adaptation:
+DEAL:3:: 2_0: 1
+DEAL:3::OK
+

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.