From: Marc Fehling Date: Thu, 22 Apr 2021 03:18:44 +0000 (-0600) Subject: Use `hp::Refinement::predict_error()` in favor of `parallel::distributed::ErrorPredic... X-Git-Tag: v9.3.0-rc1~184^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8d8d37aa02dbebc7842755f61ac698dfdf51acfa;p=dealii.git Use `hp::Refinement::predict_error()` in favor of `parallel::distributed::ErrorPredictor`. --- diff --git a/doc/news/changes/incompatibilities/20210421Fehling b/doc/news/changes/incompatibilities/20210421Fehling new file mode 100644 index 0000000000..e339c8a6c7 --- /dev/null +++ b/doc/news/changes/incompatibilities/20210421Fehling @@ -0,0 +1,6 @@ +Removed: The class parallel::distributed::ErrorPredictor has been +removed. Use the function hp::Refinement::predict_error() in combination +with parallel::distributed::CellDataTransfer instead. Please consult the +documentation of hp::Refinement::predict_error() for instructions. +
+(Marc Fehling, 2021/04/21) diff --git a/include/deal.II/hp/refinement.h b/include/deal.II/hp/refinement.h index 0b7be5b5c0..05726c7b84 100644 --- a/include/deal.II/hp/refinement.h +++ b/include/deal.II/hp/refinement.h @@ -559,6 +559,32 @@ namespace hp * $\gamma_\text{p}^2 = 0.4$, $\gamma_\text{h}^2 = 4$, * $\gamma_\text{n}^2 = 1$. * + * If you are working with parallel::distributed::Triangulation objects, you + * need to pay special attention. Here, p4est determines the details of grid + * refinement, and consequently, it yields more reliable and trustworthy + * results when we determine the predicted errors during the adaptation + * process. We can do exactly this by attaching this function to the signal + * Triangulation::Signals::post_p4est_refinement, which is triggered after + * p4est got refined, but before data is prepared for transfer. Refinement + * and coarsening flags of the Triangulation object need to be matched with + * the already refined p4est oracle using + * internal::parallel::distributed::TemporarilyMatchRefineFlags. + * Thus, a construct like the following is necessary to correctly predict + * errors in parallel distributed applications. + * @code + * Vector predicted_errors; + * triangulation.signals.post_p4est_refinement.connect([&]() { + * const internal::parallel::distributed::TemporarilyMatchRefineFlags + * refine_modifier(triangulation); + * predicted_errors.reinit(triangulation.n_active_cells()); + * hp::Refinement::predict_error(dof_handler, + * error_indicators, + * predicted_errors); + * }); + * @endcode + * The container predicted_errors then needs to follow the + * usual parallel::distributed::CellDataTransfer workflow. + * * @note We want to predict the error by how adaptation will actually happen. * Thus, this function needs to be called after * Triangulation::prepare_coarsening_and_refinement() and diff --git a/tests/mpi/error_predictor_01.cc b/tests/mpi/error_prediction_01.cc similarity index 76% rename from tests/mpi/error_predictor_01.cc rename to tests/mpi/error_prediction_01.cc index f531365914..4b692b25cc 100644 --- a/tests/mpi/error_predictor_01.cc +++ b/tests/mpi/error_prediction_01.cc @@ -19,14 +19,15 @@ // This tests is based on hp/error_prediction.cc -#include +#include #include #include -#include #include +#include + #include #include "../tests.h" @@ -83,6 +84,20 @@ test() for (unsigned int i = 0; i < error_indicators.size(); ++i) error_indicators(i) = 10.; + // ----- connect error predictor ----- + Vector predicted_errors; + tria.signals.post_p4est_refinement.connect([&]() { + const internal::parallel::distributed::TemporarilyMatchRefineFlags + refine_modifier(tria); + predicted_errors.reinit(tria.n_active_cells()); + hp::Refinement::predict_error(dh, + error_indicators, + predicted_errors, + /*gamma_p=*/0.5, + /*gamma_h=*/1., + /*gamma_n=*/1.); + }); + // ----- verify ------ deallog << "pre_adaptation" << std::endl; for (const auto &cell : dh.active_cell_iterators()) @@ -104,16 +119,17 @@ test() } // ----- execute adaptation ----- - parallel::distributed::ErrorPredictor predictor(dh); + parallel::distributed::CellDataTransfer> + data_transfer(tria, + /*transfer_variable_size_data=*/false, + &AdaptationStrategies::Refinement::l2_norm, + &AdaptationStrategies::Coarsening::l2_norm); - predictor.prepare_for_coarsening_and_refinement(error_indicators, - /*gamma_p=*/0.5, - /*gamma_h=*/1., - /*gamma_n=*/1.); + data_transfer.prepare_for_coarsening_and_refinement(predicted_errors); tria.execute_coarsening_and_refinement(); - Vector predicted_errors(tria.n_active_cells()); - predictor.unpack(predicted_errors); + predicted_errors.reinit(tria.n_active_cells()); + data_transfer.unpack(predicted_errors); // ------ verify ------ deallog << "post_adaptation" << std::endl; diff --git a/tests/mpi/error_predictor_01.with_p4est=true.mpirun=1.output b/tests/mpi/error_prediction_01.with_p4est=true.mpirun=1.output similarity index 100% rename from tests/mpi/error_predictor_01.with_p4est=true.mpirun=1.output rename to tests/mpi/error_prediction_01.with_p4est=true.mpirun=1.output diff --git a/tests/mpi/error_predictor_01.with_p4est=true.mpirun=4.output b/tests/mpi/error_prediction_01.with_p4est=true.mpirun=4.output similarity index 100% rename from tests/mpi/error_predictor_01.with_p4est=true.mpirun=4.output rename to tests/mpi/error_prediction_01.with_p4est=true.mpirun=4.output diff --git a/tests/mpi/error_predictor_02.cc b/tests/mpi/error_prediction_02.cc similarity index 80% rename from tests/mpi/error_predictor_02.cc rename to tests/mpi/error_prediction_02.cc index 3de8ff0bbd..e999db84ce 100644 --- a/tests/mpi/error_predictor_02.cc +++ b/tests/mpi/error_prediction_02.cc @@ -19,14 +19,15 @@ // This tests is based on hp/error_prediction.cc -#include +#include #include #include -#include #include +#include + #include #include "../tests.h" @@ -109,6 +110,20 @@ test() for (unsigned int i = 0; i < error_indicators.size(); ++i) error_indicators(i) = 10.; + // ----- connect error predictor ----- + Vector predicted_errors; + tria.signals.post_p4est_refinement.connect([&]() { + const internal::parallel::distributed::TemporarilyMatchRefineFlags + refine_modifier(tria); + predicted_errors.reinit(tria.n_active_cells()); + hp::Refinement::predict_error(dh, + error_indicators, + predicted_errors, + /*gamma_p=*/0.5, + /*gamma_h=*/1., + /*gamma_n=*/1.); + }); + // ----- verify ------ deallog << "pre_adaptation" << std::endl; for (const auto &cell : dh.active_cell_iterators()) @@ -130,16 +145,17 @@ test() } // ----- execute adaptation ----- - parallel::distributed::ErrorPredictor predictor(dh); + parallel::distributed::CellDataTransfer> + data_transfer(tria, + /*transfer_variable_size_data=*/false, + &AdaptationStrategies::Refinement::l2_norm, + &AdaptationStrategies::Coarsening::l2_norm); - predictor.prepare_for_coarsening_and_refinement(error_indicators, - /*gamma_p=*/0.5, - /*gamma_h=*/1., - /*gamma_n=*/1.); + data_transfer.prepare_for_coarsening_and_refinement(predicted_errors); tria.execute_coarsening_and_refinement(); - Vector predicted_errors(tria.n_active_cells()); - predictor.unpack(predicted_errors); + predicted_errors.reinit(tria.n_active_cells()); + data_transfer.unpack(predicted_errors); // ------ verify ------ deallog << "post_adaptation" << std::endl; diff --git a/tests/mpi/error_predictor_02.with_p4est=true.mpirun=1.output b/tests/mpi/error_prediction_02.with_p4est=true.mpirun=1.output similarity index 100% rename from tests/mpi/error_predictor_02.with_p4est=true.mpirun=1.output rename to tests/mpi/error_prediction_02.with_p4est=true.mpirun=1.output diff --git a/tests/mpi/error_predictor_02.with_p4est=true.mpirun=4.output b/tests/mpi/error_prediction_02.with_p4est=true.mpirun=4.output similarity index 100% rename from tests/mpi/error_predictor_02.with_p4est=true.mpirun=4.output rename to tests/mpi/error_prediction_02.with_p4est=true.mpirun=4.output