]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Reworked error prediction algorithm.
authorMarc Fehling <marc.fehling@gmx.net>
Fri, 24 Apr 2020 21:24:10 +0000 (23:24 +0200)
committerMatthias Maier <tamiko@43-1.org>
Mon, 11 May 2020 20:23:41 +0000 (15:23 -0500)
doc/doxygen/references.bib
include/deal.II/hp/refinement.h
source/hp/refinement.cc
tests/hp/error_prediction.output [deleted file]
tests/hp/error_prediction_01.cc [moved from tests/hp/error_prediction.cc with 53% similarity]
tests/hp/error_prediction_01.output [new file with mode: 0644]

index b0a4497779f72abd1af43236f497970a40282267..3df58194b8aaf290ef8c368af8239365896680de 100644 (file)
@@ -724,3 +724,27 @@ year = {2008},
   Year                     = {2019},
   Url                      = {https://arxiv.org/abs/1904.03317}
 }
+
+@article{ainsworth1998hp,
+  author    = {Ainsworth, Mark and Senior, Bill},
+  title     = {An adaptive refinement strategy for hp-finite element computations},
+  journal   = {{Applied Numerical Mathematics}},
+  volume    = {26},
+  number    = {1--2},
+  pages     = {165--178},
+  publisher = {Elsevier},
+  year      = {1998},
+  doi       = {10.1016/S0168-9274(97)00083-4}
+}
+
+@article{melenk2001hp,
+  author    = {Melenk, Jens Markus and Wohlmuth, Barbara I.},
+  title     = {{On residual-based a posteriori error estimation in hp-FEM}},
+  journal   = {{Advances in Computational Mathematics}},
+  volume    = {15},
+  number    = {1},
+  pages     = {311--331},
+  publisher = {Springer US},
+  year      = {2001},
+  doi       = {10.1023/A:1014268310921}
+}
index 504cd7b8055888fc61ca989ca505a81ed075da87..10756ba5d0f43766e39830d5c8011af128456afc 100644 (file)
@@ -151,7 +151,7 @@ namespace hp
      * Each cell flagged for h-refinement will also be flagged for p-refinement.
      * The same applies to coarsening.
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
@@ -168,7 +168,7 @@ namespace hp
      * Each entry of the parameter @p p_flags needs to correspond to an active
      * cell.
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
@@ -197,7 +197,7 @@ namespace hp
      * Each entry of the parameter @p criteria needs to correspond to an active
      * cell.
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
@@ -239,7 +239,7 @@ namespace hp
      * cell. Parameters @p p_refine_fraction and @p p_coarsen_fraction need to be
      * in the interval $[0,1]$.
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
@@ -282,7 +282,7 @@ namespace hp
      * cell. Parameters @p p_refine_fraction and @p p_coarsen_fraction need to be
      * in the interval $[0,1]$.
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
@@ -317,23 +317,9 @@ namespace hp
      * Each entry of the parameter @p sobolev_indices needs to correspond
      * to an active cell.
      *
-     * For more theoretical details see
-     * @code{.bib}
-     * @article{Ainsworth1998,
-     *  author    = {Ainsworth, Mark and Senior, Bill},
-     *  title     = {An adaptive refinement strategy for hp-finite element
-     *               computations},
-     *  journal   = {{Applied Numerical Mathematics}},
-     *  volume    = {26},
-     *  number    = {1--2},
-     *  pages     = {165--178},
-     *  publisher = {Elsevier},
-     *  year      = {1998},
-     *  doi       = {10.1016/S0168-9274(97)00083-4}
-     * }
-     * @endcode
+     * For more theoretical details see @cite ainsworth1998hp .
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
@@ -358,7 +344,7 @@ namespace hp
      * Each entry of the parameters @p criteria and @p references needs to
      * correspond to an active cell.
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
@@ -388,18 +374,26 @@ namespace hp
      * Triangulation, thus each container has to be of size
      * Triangulation::n_active_cells(). The errors are interpreted to be
      * measured in the energy norm; this assumption enters the rate of
-     * convergence that is used in the prediction. The `predicted_errors` output
-     * argument has one entry per <i>current</i> cell, with the $2^d$ values for
-     * each cell that will be coarsened away equal, and with the value stored on
-     * a cell to be refined interpreted as applying to each of the future
-     * children.
+     * convergence that is used in the prediction. The $l_2$-norm of the output
+     * argument @p predicted_errors corresponds to the predicted global error
+     * after adaptation will be applied.
+     *
+     * For p-adaptation, the local error is expected to converge exponentially
+     * with the polynomial degree of the assigned finite element. Each increase
+     * or decrease of the degree will thus change its value by a user-defined
+     * control parameter @p gamma_p.
      *
      * For h-adaptation, we expect the local error $\eta_K$ on cell $K$ to be
      * proportional to $(h_K)^{p_K}$ in the energy norm, where $h_K$ denotes the
      * cell diameter and $p_K$ the polynomial degree of the currently assigned
-     * finite element on cell $K$. Here, we assume that the finite element will
-     * not change in the adaptation process so that $p_K = \text{const}$.
-     * However during coarsening, the finite elements on siblings may be
+     * finite element on cell $K$.
+     *
+     * If both h- and p-adaptation are applied simultaneously, we need to
+     * determine the order at which which type of adaptation happensis
+     * performed. We perform p-adaptation first and perform h-adaptation with
+     * the degree of the future finite element $p_{K,\text{future}}$.
+     *
+     * During h-coarsening, the finite elements on siblings may be
      * different, and their parent cell will be assigned to their least
      * dominating finite element that belongs to its most general child. Thus,
      * we will always interpolate on an enclosing finite element space.
@@ -432,25 +426,79 @@ namespace hp
      *
      * The prediction algorithm is formulated as follows with control parameters
      * @p gamma_p, @p gamma_h and @p gamma_n that may be used to influence
-     * prediction for each adaptation type individually.
+     * prediction for each adaptation type individually. The results for each
+     * individual cell are stored in the @p predicted_errors output argument.
      * <table>
      *   <tr><th>Adaptation type <th colspan="2">Prediction formula
      *   <tr><td>no adaptation
-     *       <td>$\eta_{K,\text{pred}} = \eta_{K} \, \gamma_\text{n}$
+     *       <td>$\eta_{K,\text{pred}} = \eta_{K} \,
+     *            \gamma_\text{n}$
      *       <td>$\gamma_\text{n} \in (0,\infty)$
      *   <tr><td>p-adaptation
      *       <td>$\eta_{K,\text{pred}} = \eta_{K} \,
      *            \gamma_\text{p}^{(p_{K,\text{future}} - p_K)}$
      *       <td>$\gamma_\text{p} \in (0,1)$
-     *   <tr><td>h-refinement
-     *       <td>$\eta_{K_c,\text{pred}} = \eta_{K} \,
-     *            \gamma_\text{h} \, 0.5^{p_K} \, 0.5^{\text{dim}}
-     *            \quad \forall K_c \text{ children of } K$
+     *   <tr><td>hp-refinement
+     *       <td>$\eta_{K,\text{pred}} = \eta_{K} \,
+     *            \gamma_\text{h} \, 0.5^{p_{K,\text{future}}} \,
+     *            \gamma_\text{p}^{(p_{K,\text{future}} - p_{K})}$
      *       <td rowspan="2">$\gamma_\text{h} \in (0,\infty)$
-     *   <tr><td>h-coarsening
-     *       <td>$\eta_{K,\text{pred}} = \sum\limits_{K_c} \eta_{K_c} /
-     *            (\gamma_\text{h} \, 0.5^{p_{K_c}})
-     *            \quad \forall K_c \text{ children of } K$
+     *   <tr><td>hp-coarsening
+     *       <td>$\eta_{K,\text{pred}} = \eta_{K} \,
+     *            (\gamma_\text{h} \, 0.5^{p_{K,\text{future}}})^{-1} \,
+     *            \gamma_\text{p}^{(p_{K,\text{future}} - p_{K})}$
+     * </table>
+     *
+     * On basis of the refinement history, we use the predicted error estimates
+     * to decide how cells will be adapted in the next adaptation step.
+     * Comparing the predicted error from the previous adaptation step to the
+     * error estimates of the current step allows us to justify whether our
+     * previous choice of adaptation was justified, and lets us decide how to
+     * adapt in the next one.
+     *
+     * We thus have to transfer the predicted error from the old to the adapted
+     * mesh. When transferring the predicted error to the adapted mesh, make
+     * sure to configure your CellDataTransfer object with
+     * AdaptationStrategies::Refinement::l2_norm() as a refinement strategy and
+     * AdaptationStrategies::Coarsening::l2_norm() as a coarsening strategy.
+     * This ensures that the $l_2$-norm of the predict errors is preserved on
+     * both meshes.
+     *
+     * In the context, we assume that the local error on a cell that will be
+     * h-refined, will lead to errors on the $2^\text{dim}$ children that are
+     * all equal, whereas local errors on siblings will be summed up on the
+     * parent cell in case of h-coarsening. This assumption is often not
+     * satisfied in practice: For example, if a cell is at a corner singularity,
+     * then the one child cell that ends up closest to the singularity will
+     * inherit the majority of the remaining error -- but this function can not
+     * know where the singularity will be, and consequently assumes equal
+     * distribution.
+     *
+     * Incorporating the transfer from the old to the adapted mesh, the complete
+     * error prediction algorithm reads as follows:
+     * <table>
+     *   <tr><th>Adaptation type <th colspan="2">Prediction formula
+     *   <tr><td>no adaptation
+     *       <td>$\eta_{K,\text{pred}} = \eta_{K} \,
+     *            \gamma_\text{n}$
+     *       <td>$\gamma_\text{n} \in (0,\infty)$
+     *   <tr><td>p-adaptation
+     *       <td>$\eta_{K,\text{pred}} = \eta_{K} \,
+     *            \gamma_\text{p}^{(p_{K,\text{future}} - p_K)}$
+     *       <td>$\gamma_\text{p} \in (0,1)$
+     *   <tr><td>hp-refinement
+     *       <td>$\left( \eta_{K_c,\text{pred}} \right)^2 = 0.5^{\text{dim}}
+     *            \left( \eta_{K_p} \,
+     *            \gamma_\text{h} \, 0.5^{p_{K_c,\text{future}}} \,
+     *            \gamma_\text{p}^{(p_{K_c,\text{future}} - p_{K_p})} \right)^2
+     *            \quad \forall K_c \text{ children of } K_p$
+     *       <td rowspan="2">$\gamma_\text{h} \in (0,\infty)$
+     *   <tr><td>hp-coarsening
+     *       <td>$\left( \eta_{K_p,\text{pred}} \right)^2 = \sum\limits_{K_c}
+     *            \left( \eta_{K_c} \,
+     *            (\gamma_\text{h} \, 0.5^{p_{K_p,\text{future}}})^{-1} \,
+     *            \gamma_\text{p}^{(p_{K_p,\text{future}} - p_{K_c})} \right)^2
+     *            \quad \forall K_c \text{ children of } K_p$
      * </table>
      *
      * With these predicted error estimates, we are capable of adapting the
@@ -502,9 +550,11 @@ namespace hp
      *   predicted_error_per_cell);
      *
      * // perform adaptation
-     * CellDataTransfer<dim, spacedim, Vector<float>>
-     *   cell_data_transfer(triangulation);
-     * cell_data_transfer.prepare_for_coarsening_and_refinement();
+     * CellDataTransfer<dim, spacedim, Vector<float>> cell_data_transfer(
+     *   triangulation,
+     *   &AdaptationStrategies::Refinement::l2_norm<dim, spacedim, float>,
+     *   &AdaptationStrategies::Coarsening::l2_norm<dim, spacedim, float>);
+     * cell_data_transfer.prepare_coarsening_and_refinement();
      *
      * triangulation.execute_coarsening_and_refinement();
      *
@@ -513,35 +563,24 @@ namespace hp
      * predicted_error_per_cell = std::move(transferred_errors);
      * @endcode
      *
-     * For more theoretical details see
-     * @code{.bib}
-     * @article{Melenk2001,
-     *  author    = {Melenk, Jens Markus and Wohlmuth, Barbara I.},
-     *  title     = {{On residual-based a posteriori error estimation
-     *                in hp-FEM}},
-     *  journal   = {{Advances in Computational Mathematics}},
-     *  volume    = {15},
-     *  number    = {1},
-     *  pages     = {311--331},
-     *  publisher = {Springer US},
-     *  year      = {2001},
-     *  doi       = {10.1023/A:1014268310921}
-     * }
-     * @endcode
+     * For more theoretical details see @cite melenk2001hp , where the default
+     * parameters for this function come from as well, i.e.
+     * $\gamma_\text{p}^2 = 0.4$, $\gamma_\text{h}^2 = 4$,
+     * $\gamma_\text{n}^2 = 1$.
      *
      * @note This feature is currently only implemented for isotropic refinement.
      *
      * @note We want to predict the error by how adaptation will actually happen.
      *   Thus, this function needs to be called after
-     *   Triangulation::prepare_for_coarsening_and_refinement().
+     *   Triangulation::prepare_coarsening_and_refinement().
      */
     template <int dim, typename Number, int spacedim>
     void
     predict_error(const hp::DoFHandler<dim, spacedim> &dof_handler,
                   const Vector<Number> &               error_indicators,
                   Vector<Number> &                     predicted_errors,
-                  const double                         gamma_p = std::sqrt(0.1),
-                  const double                         gamma_h = 1.,
+                  const double                         gamma_p = std::sqrt(0.4),
+                  const double                         gamma_h = 2.,
                   const double                         gamma_n = 1.);
 
     /**
@@ -559,7 +598,7 @@ namespace hp
      * Removes all refine and coarsen flags on cells that have a
      * @p future_fe_index assigned.
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
@@ -606,7 +645,7 @@ namespace hp
      *   the decision that Triangulation::prepare_coarsening_and_refinement()
      *   would have made later on.
      *
-     * @note Triangulation::prepare_for_coarsening_and_refinement() may change
+     * @note Triangulation::prepare_coarsening_and_refinement() may change
      *   refine and coarsen flags. Avoid calling it before this particular
      *   function.
      */
index d48e0f8a8cf68e4c9f17a5e06b9ca04fa1e6aad8..f47d71f4d3c493dbb1165a60aa7d0f414eddd783 100644 (file)
@@ -511,30 +511,93 @@ namespace hp
       Assert(0 < gamma_h, dealii::GridRefinement::ExcInvalidParameterValue());
       Assert(0 < gamma_n, dealii::GridRefinement::ExcInvalidParameterValue());
 
+      // auxiliary variables
+      unsigned int future_fe_degree       = numbers::invalid_unsigned_int;
+      unsigned int parent_future_fe_index = numbers::invalid_unsigned_int;
+      // store all determined future finite element indices on parent cells for
+      // coarsening
+      std::map<typename hp::DoFHandler<dim, spacedim>::cell_iterator,
+               unsigned int>
+        future_fe_indices_on_coarsened_cells;
+
+      // deep copy error indicators
+      predicted_errors = error_indicators;
+
       for (const auto &cell : dof_handler.active_cell_iterators())
         if (cell->is_locally_owned())
           {
-            if (cell->future_fe_index_set()) // p adaptation
+            // current cell will not be adapted
+            if (!(cell->future_fe_index_set()) && !(cell->refine_flag_set()) &&
+                !(cell->coarsen_flag_set()))
+              {
+                predicted_errors[cell->active_cell_index()] *= gamma_n;
+                continue;
+              }
+
+            // current cell will be adapted
+            // determine degree of its future finite element
+            if (cell->coarsen_flag_set())
+              {
+                // cell will be coarsened, thus determine future finite element
+                // on parent cell
+                const auto &parent = cell->parent();
+                if (future_fe_indices_on_coarsened_cells.find(parent) ==
+                    future_fe_indices_on_coarsened_cells.end())
+                  {
+                    std::set<unsigned int> fe_indices_children;
+                    for (unsigned int child_index = 0;
+                         child_index < parent->n_children();
+                         ++child_index)
+                      {
+                        const auto &child = parent->child(child_index);
+                        Assert(child->is_active() && child->coarsen_flag_set(),
+                               typename dealii::Triangulation<
+                                 dim>::ExcInconsistentCoarseningFlags());
+
+                        fe_indices_children.insert(child->future_fe_index());
+                      }
+                    Assert(!fe_indices_children.empty(), ExcInternalError());
+
+                    parent_future_fe_index =
+                      dof_handler.get_fe_collection()
+                        .find_dominated_fe_extended(fe_indices_children,
+                                                    /*codim=*/0);
+
+                    Assert(
+                      parent_future_fe_index != numbers::invalid_unsigned_int,
+                      typename dealii::hp::FECollection<
+                        dim>::ExcNoDominatedFiniteElementAmongstChildren());
+
+                    future_fe_indices_on_coarsened_cells.insert(
+                      {parent, parent_future_fe_index});
+                  }
+                else
+                  {
+                    parent_future_fe_index =
+                      future_fe_indices_on_coarsened_cells[parent];
+                  }
+
+                future_fe_degree =
+                  dof_handler.get_fe_collection()[parent_future_fe_index]
+                    .degree;
+              }
+            else
               {
-                Assert(cell->future_fe_index_set() && !cell->refine_flag_set(),
-                       ExcMessage(
-                         "For error prediction, a cell marked for p-adaptation "
-                         "should not also be flagged for h-refinement!"));
-                Assert(cell->future_fe_index_set() && !cell->coarsen_flag_set(),
-                       ExcMessage(
-                         "For error prediction, a cell marked for p-adaptation "
-                         "should not also be flagged for h-coarsening!"));
-
-                const int degree_difference =
+                // future finite element on current cell is already set
+                future_fe_degree =
                   dof_handler.get_fe_collection()[cell->future_fe_index()]
-                    .degree -
-                  cell->get_fe().degree;
+                    .degree;
+              }
 
-                predicted_errors[cell->active_cell_index()] =
-                  error_indicators[cell->active_cell_index()] *
-                  std::pow(gamma_p, degree_difference);
+            // step 1: exponential decay with p-adaptation
+            if (cell->future_fe_index_set())
+              {
+                predicted_errors[cell->active_cell_index()] *=
+                  std::pow(gamma_p, future_fe_degree - cell->get_fe().degree);
               }
-            else if (cell->refine_flag_set()) // h refinement
+
+            // step 2: algebraic decay with h-adaptation
+            if (cell->refine_flag_set())
               {
                 Assert(
                   cell->refine_flag_set() ==
@@ -542,20 +605,19 @@ namespace hp
                   ExcMessage(
                     "Error prediction is only valid for isotropic refinement!"));
 
-                predicted_errors[cell->active_cell_index()] =
-                  error_indicators[cell->active_cell_index()] *
-                  (gamma_h * std::pow(.5, dim + cell->get_fe().degree));
-              }
-            else if (cell->coarsen_flag_set()) // h coarsening
-              {
-                predicted_errors[cell->active_cell_index()] =
-                  error_indicators[cell->active_cell_index()] /
-                  (gamma_h * std::pow(.5, cell->get_fe().degree));
+                predicted_errors[cell->active_cell_index()] *=
+                  (gamma_h * std::pow(.5, future_fe_degree));
+
+                // predicted error will be split on children cells
+                // after adaptation via CellDataTransfer
               }
-            else // no changes
+            else if (cell->coarsen_flag_set())
               {
-                predicted_errors[cell->active_cell_index()] =
-                  error_indicators[cell->active_cell_index()] * gamma_n;
+                predicted_errors[cell->active_cell_index()] /=
+                  (gamma_h * std::pow(.5, future_fe_degree));
+
+                // predicted error will be summed up on parent cell
+                // after adaptation via CellDataTransfer
               }
           }
     }
diff --git a/tests/hp/error_prediction.output b/tests/hp/error_prediction.output
deleted file mode 100644 (file)
index 1d11f6b..0000000
+++ /dev/null
@@ -1,19 +0,0 @@
-
-DEAL:1d::ncells:4
-DEAL:1d:: cell:0_0: fe_deg:1 error:10.0 predict:2.50 refining
-DEAL:1d:: cell:1_0: fe_deg:1 error:10.0 predict:20.0 coarsening
-DEAL:1d:: cell:2_0: fe_deg:1 error:10.0 predict:2.50 future_fe_deg:3
-DEAL:1d:: cell:3_0: fe_deg:1 error:10.0 predict:10.0
-DEAL:1d::OK
-DEAL:2d::ncells:4
-DEAL:2d:: cell:0_0: fe_deg:1 error:10.0 predict:1.25 refining
-DEAL:2d:: cell:1_0: fe_deg:1 error:10.0 predict:20.0 coarsening
-DEAL:2d:: cell:2_0: fe_deg:1 error:10.0 predict:2.50 future_fe_deg:3
-DEAL:2d:: cell:3_0: fe_deg:1 error:10.0 predict:10.0
-DEAL:2d::OK
-DEAL:3d::ncells:4
-DEAL:3d:: cell:0_0: fe_deg:1 error:10.0 predict:0.625 refining
-DEAL:3d:: cell:1_0: fe_deg:1 error:10.0 predict:20.0 coarsening
-DEAL:3d:: cell:2_0: fe_deg:1 error:10.0 predict:2.50 future_fe_deg:3
-DEAL:3d:: cell:3_0: fe_deg:1 error:10.0 predict:10.0
-DEAL:3d::OK
similarity index 53%
rename from tests/hp/error_prediction.cc
rename to tests/hp/error_prediction_01.cc
index 8d76702f5683825b61730ae0535f5a65ac2b6a8c..c457fbb45161845c42c36c017bbaa0451d0751b0 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2019 by the deal.II authors
+// Copyright (C) 2019 - 2020 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -15,7 +15,8 @@
 
 
 
-// validate error prediction algorithm for hp adaptive methods
+// validate combination of error prediction and cell data transfer algorithms
+// for hp adaptive methods
 
 
 #include <deal.II/fe/fe_q.h>
@@ -29,6 +30,9 @@
 
 #include <deal.II/lac/vector.h>
 
+#include <deal.II/numerics/adaptation_strategies.h>
+#include <deal.II/numerics/cell_data_transfer.h>
+
 #include "../tests.h"
 
 
@@ -51,44 +55,59 @@ test()
     }
   GridGenerator::subdivided_hyper_rectangle(tria, rep, p1, p2);
 
+  tria.begin_active()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+
   hp::FECollection<dim> fes;
   for (unsigned int d = 1; d <= 3; ++d)
     fes.push_back(FE_Q<dim>(d));
 
   hp::DoFHandler<dim> dh(tria);
   dh.set_fe(fes);
-  for (const auto &cell : dh.active_cell_iterators())
+  for (auto cell = dh.begin(0); cell != dh.end(0); ++cell)
     {
       if (cell->id().to_string() == "0_0:")
-        cell->set_refine_flag();
+        {
+          // h-coarsening
+          for (unsigned int i = 0; i < cell->n_children(); ++i)
+            cell->child(i)->set_coarsen_flag();
+        }
       else if (cell->id().to_string() == "1_0:")
-        cell->set_coarsen_flag();
+        {
+          // h-refinement
+          cell->set_refine_flag();
+        }
       else if (cell->id().to_string() == "2_0:")
-        cell->set_future_fe_index(2);
+        {
+          // p-refinement
+          cell->set_future_fe_index(2);
+        }
     }
 
   // ----- predict -----
-  Vector<float> error_indicators, predicted_errors;
+  Vector<float> error_indicators, predicted_error_indicators;
   error_indicators.reinit(tria.n_active_cells());
-  predicted_errors.reinit(tria.n_active_cells());
+  predicted_error_indicators.reinit(tria.n_active_cells());
   for (unsigned int i = 0; i < tria.n_active_cells(); ++i)
-    error_indicators[i] = 10;
+    error_indicators[i] = 10.;
 
   hp::Refinement::predict_error(dh,
                                 error_indicators,
-                                predicted_errors,
+                                predicted_error_indicators,
                                 /*gamma_p=*/0.5,
                                 /*gamma_h=*/1.,
                                 /*gamma_n=*/1.);
 
   // ----- verify ------
-  deallog << "ncells:" << tria.n_active_cells() << std::endl;
+  deallog << "pre_adaptation" << std::endl
+          << " ncells:" << tria.n_active_cells() << std::endl;
   for (const auto &cell : dh.active_cell_iterators())
     {
       deallog << " cell:" << cell->id().to_string()
               << " fe_deg:" << cell->get_fe().degree
               << " error:" << error_indicators[cell->active_cell_index()]
-              << " predict:" << predicted_errors[cell->active_cell_index()];
+              << " predicted:"
+              << predicted_error_indicators[cell->active_cell_index()];
 
       if (cell->refine_flag_set())
         deallog << " refining";
@@ -101,13 +120,36 @@ test()
       deallog << std::endl;
     }
 
-  // ----- check feature -----
-  hp::Refinement::p_adaptivity_from_reference(
-    dh,
-    error_indicators,
-    predicted_errors,
-    /*compare_refine=*/std::less<float>(),
-    /*compare_coarsen=*/std::less<float>());
+  // ----- execute adaptation -----
+  CellDataTransfer<dim, dim, Vector<float>> cell_data_transfer(
+    tria,
+    &AdaptationStrategies::Refinement::l2_norm<dim, dim, float>,
+    &AdaptationStrategies::Coarsening::l2_norm<dim, dim, float>);
+  cell_data_transfer.prepare_for_coarsening_and_refinement();
+
+  tria.execute_coarsening_and_refinement();
+
+  Vector<float> transferred_indicators(tria.n_active_cells());
+  cell_data_transfer.unpack(predicted_error_indicators, transferred_indicators);
+
+  // ----- verify -----
+  deallog << "post_adaptation" << std::endl
+          << " ncells:" << tria.n_active_cells() << std::endl;
+  for (const auto &cell : dh.active_cell_iterators())
+    deallog << " cell:" << cell->id().to_string()
+            << " fe_deg:" << cell->get_fe().degree << " transferred:"
+            << transferred_indicators[cell->active_cell_index()] << std::endl;
+
+  // ----- verify norms -----
+  const double predicted_error_pre  = predicted_error_indicators.l2_norm(),
+               predicted_error_post = transferred_indicators.l2_norm();
+
+  deallog << "predicted_error_norms" << std::endl
+          << " pre_adaptation:" << predicted_error_pre << std::endl
+          << " post_adaptation:" << predicted_error_post << std::endl;
+
+  Assert(predicted_error_pre == predicted_error_post,
+         ExcMessage("Transfer failed - Results not similar."));
 
   deallog << "OK" << std::endl;
 }
@@ -118,7 +160,6 @@ int
 main(int argc, char *argv[])
 {
   initlog();
-  deallog << std::setprecision(3);
 
   deallog.push("1d");
   test<1>();
diff --git a/tests/hp/error_prediction_01.output b/tests/hp/error_prediction_01.output
new file mode 100644 (file)
index 0000000..540f22e
--- /dev/null
@@ -0,0 +1,71 @@
+
+DEAL:1d::pre_adaptation
+DEAL:1d:: ncells:5
+DEAL:1d:: cell:1_0: fe_deg:1 error:10.0000 predicted:5.00000 refining
+DEAL:1d:: cell:2_0: fe_deg:1 error:10.0000 predicted:2.50000 future_fe_deg:3
+DEAL:1d:: cell:3_0: fe_deg:1 error:10.0000 predicted:10.0000
+DEAL:1d:: cell:0_1:0 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:1d:: cell:0_1:1 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:1d::post_adaptation
+DEAL:1d:: ncells:5
+DEAL:1d:: cell:0_0: fe_deg:1 transferred:28.2843
+DEAL:1d:: cell:2_0: fe_deg:3 transferred:2.50000
+DEAL:1d:: cell:3_0: fe_deg:1 transferred:10.0000
+DEAL:1d:: cell:1_1:0 fe_deg:1 transferred:3.53553
+DEAL:1d:: cell:1_1:1 fe_deg:1 transferred:3.53553
+DEAL:1d::predicted_error_norms
+DEAL:1d:: pre_adaptation:30.5164
+DEAL:1d:: post_adaptation:30.5164
+DEAL:1d::OK
+DEAL:2d::pre_adaptation
+DEAL:2d:: ncells:7
+DEAL:2d:: cell:1_0: fe_deg:1 error:10.0000 predicted:5.00000 refining
+DEAL:2d:: cell:2_0: fe_deg:1 error:10.0000 predicted:2.50000 future_fe_deg:3
+DEAL:2d:: cell:3_0: fe_deg:1 error:10.0000 predicted:10.0000
+DEAL:2d:: cell:0_1:0 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:2d:: cell:0_1:1 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:2d:: cell:0_1:2 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:2d:: cell:0_1:3 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:2d::post_adaptation
+DEAL:2d:: ncells:7
+DEAL:2d:: cell:0_0: fe_deg:1 transferred:40.0000
+DEAL:2d:: cell:2_0: fe_deg:3 transferred:2.50000
+DEAL:2d:: cell:3_0: fe_deg:1 transferred:10.0000
+DEAL:2d:: cell:1_1:0 fe_deg:1 transferred:2.50000
+DEAL:2d:: cell:1_1:1 fe_deg:1 transferred:2.50000
+DEAL:2d:: cell:1_1:2 fe_deg:1 transferred:2.50000
+DEAL:2d:: cell:1_1:3 fe_deg:1 transferred:2.50000
+DEAL:2d::predicted_error_norms
+DEAL:2d:: pre_adaptation:41.6083
+DEAL:2d:: post_adaptation:41.6083
+DEAL:2d::OK
+DEAL:3d::pre_adaptation
+DEAL:3d:: ncells:11
+DEAL:3d:: cell:1_0: fe_deg:1 error:10.0000 predicted:5.00000 refining
+DEAL:3d:: cell:2_0: fe_deg:1 error:10.0000 predicted:2.50000 future_fe_deg:3
+DEAL:3d:: cell:3_0: fe_deg:1 error:10.0000 predicted:10.0000
+DEAL:3d:: cell:0_1:0 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:3d:: cell:0_1:1 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:3d:: cell:0_1:2 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:3d:: cell:0_1:3 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:3d:: cell:0_1:4 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:3d:: cell:0_1:5 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:3d:: cell:0_1:6 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:3d:: cell:0_1:7 fe_deg:1 error:10.0000 predicted:20.0000 coarsening
+DEAL:3d::post_adaptation
+DEAL:3d:: ncells:11
+DEAL:3d:: cell:0_0: fe_deg:1 transferred:56.5685
+DEAL:3d:: cell:2_0: fe_deg:3 transferred:2.50000
+DEAL:3d:: cell:3_0: fe_deg:1 transferred:10.0000
+DEAL:3d:: cell:1_1:0 fe_deg:1 transferred:1.76777
+DEAL:3d:: cell:1_1:1 fe_deg:1 transferred:1.76777
+DEAL:3d:: cell:1_1:2 fe_deg:1 transferred:1.76777
+DEAL:3d:: cell:1_1:3 fe_deg:1 transferred:1.76777
+DEAL:3d:: cell:1_1:4 fe_deg:1 transferred:1.76777
+DEAL:3d:: cell:1_1:5 fe_deg:1 transferred:1.76777
+DEAL:3d:: cell:1_1:6 fe_deg:1 transferred:1.76777
+DEAL:3d:: cell:1_1:7 fe_deg:1 transferred:1.76777
+DEAL:3d::predicted_error_norms
+DEAL:3d:: pre_adaptation:57.7170
+DEAL:3d:: post_adaptation:57.7170
+DEAL:3d::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.