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}
+}
* 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.
*/
* 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.
*/
* 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.
*/
* 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.
*/
* 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.
*/
* 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.
*/
* 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.
*/
* 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.
*
* 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
* 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();
*
* 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.);
/**
* 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.
*/
* 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.
*/
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() ==
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
}
}
}
+++ /dev/null
-
-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
// ---------------------------------------------------------------------
//
-// 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.
//
-// 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>
#include <deal.II/lac/vector.h>
+#include <deal.II/numerics/adaptation_strategies.h>
+#include <deal.II/numerics/cell_data_transfer.h>
+
#include "../tests.h"
}
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";
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;
}
main(int argc, char *argv[])
{
initlog();
- deallog << std::setprecision(3);
deallog.push("1d");
test<1>();
--- /dev/null
+
+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