From d8de253be31f13e3fb7e7ea779dff5c5fac9aad0 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Fri, 24 Apr 2020 23:24:10 +0200 Subject: [PATCH] Reworked error prediction algorithm. --- doc/doxygen/references.bib | 24 +++ include/deal.II/hp/refinement.h | 165 +++++++++++------- source/hp/refinement.cc | 120 ++++++++++--- tests/hp/error_prediction.output | 19 -- ...r_prediction.cc => error_prediction_01.cc} | 81 ++++++--- tests/hp/error_prediction_01.output | 71 ++++++++ 6 files changed, 349 insertions(+), 131 deletions(-) delete mode 100644 tests/hp/error_prediction.output rename tests/hp/{error_prediction.cc => error_prediction_01.cc} (53%) create mode 100644 tests/hp/error_prediction_01.output diff --git a/doc/doxygen/references.bib b/doc/doxygen/references.bib index b0a4497779..3df58194b8 100644 --- a/doc/doxygen/references.bib +++ b/doc/doxygen/references.bib @@ -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} +} diff --git a/include/deal.II/hp/refinement.h b/include/deal.II/hp/refinement.h index 504cd7b805..10756ba5d0 100644 --- a/include/deal.II/hp/refinement.h +++ b/include/deal.II/hp/refinement.h @@ -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 current 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. * *
Adaptation type Prediction formula *
no adaptation - * $\eta_{K,\text{pred}} = \eta_{K} \, \gamma_\text{n}$ + * $\eta_{K,\text{pred}} = \eta_{K} \, + * \gamma_\text{n}$ * $\gamma_\text{n} \in (0,\infty)$ *
p-adaptation * $\eta_{K,\text{pred}} = \eta_{K} \, * \gamma_\text{p}^{(p_{K,\text{future}} - p_K)}$ * $\gamma_\text{p} \in (0,1)$ - *
h-refinement - * $\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$ + *
hp-refinement + * $\eta_{K,\text{pred}} = \eta_{K} \, + * \gamma_\text{h} \, 0.5^{p_{K,\text{future}}} \, + * \gamma_\text{p}^{(p_{K,\text{future}} - p_{K})}$ * $\gamma_\text{h} \in (0,\infty)$ - *
h-coarsening - * $\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$ + *
hp-coarsening + * $\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})}$ + *
+ * + * 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: + * + *
Adaptation type Prediction formula + *
no adaptation + * $\eta_{K,\text{pred}} = \eta_{K} \, + * \gamma_\text{n}$ + * $\gamma_\text{n} \in (0,\infty)$ + *
p-adaptation + * $\eta_{K,\text{pred}} = \eta_{K} \, + * \gamma_\text{p}^{(p_{K,\text{future}} - p_K)}$ + * $\gamma_\text{p} \in (0,1)$ + *
hp-refinement + * $\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$ + * $\gamma_\text{h} \in (0,\infty)$ + *
hp-coarsening + * $\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$ *
* * With these predicted error estimates, we are capable of adapting the @@ -502,9 +550,11 @@ namespace hp * predicted_error_per_cell); * * // perform adaptation - * CellDataTransfer> - * cell_data_transfer(triangulation); - * cell_data_transfer.prepare_for_coarsening_and_refinement(); + * CellDataTransfer> cell_data_transfer( + * triangulation, + * &AdaptationStrategies::Refinement::l2_norm, + * &AdaptationStrategies::Coarsening::l2_norm); + * 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 void predict_error(const hp::DoFHandler &dof_handler, const Vector & error_indicators, Vector & 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. */ diff --git a/source/hp/refinement.cc b/source/hp/refinement.cc index d48e0f8a8c..f47d71f4d3 100644 --- a/source/hp/refinement.cc +++ b/source/hp/refinement.cc @@ -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::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 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 index 1d11f6b76b..0000000000 --- a/tests/hp/error_prediction.output +++ /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 diff --git a/tests/hp/error_prediction.cc b/tests/hp/error_prediction_01.cc similarity index 53% rename from tests/hp/error_prediction.cc rename to tests/hp/error_prediction_01.cc index 8d76702f56..c457fbb451 100644 --- a/tests/hp/error_prediction.cc +++ b/tests/hp/error_prediction_01.cc @@ -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 @@ -29,6 +30,9 @@ #include +#include +#include + #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 fes; for (unsigned int d = 1; d <= 3; ++d) fes.push_back(FE_Q(d)); hp::DoFHandler 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 error_indicators, predicted_errors; + Vector 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(), - /*compare_coarsen=*/std::less()); + // ----- execute adaptation ----- + CellDataTransfer> cell_data_transfer( + tria, + &AdaptationStrategies::Refinement::l2_norm, + &AdaptationStrategies::Coarsening::l2_norm); + cell_data_transfer.prepare_for_coarsening_and_refinement(); + + tria.execute_coarsening_and_refinement(); + + Vector 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 index 0000000000..540f22e27e --- /dev/null +++ b/tests/hp/error_prediction_01.output @@ -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 -- 2.39.5