]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use `hp::Refinement::predict_error()` in favor of `parallel::distributed::ErrorPredic...
authorMarc Fehling <mafehling.git@gmail.com>
Thu, 22 Apr 2021 03:18:44 +0000 (21:18 -0600)
committerMarc Fehling <mafehling.git@gmail.com>
Thu, 22 Apr 2021 21:58:46 +0000 (15:58 -0600)
doc/news/changes/incompatibilities/20210421Fehling [new file with mode: 0644]
include/deal.II/hp/refinement.h
tests/mpi/error_prediction_01.cc [moved from tests/mpi/error_predictor_01.cc with 76% similarity]
tests/mpi/error_prediction_01.with_p4est=true.mpirun=1.output [moved from tests/mpi/error_predictor_01.with_p4est=true.mpirun=1.output with 100% similarity]
tests/mpi/error_prediction_01.with_p4est=true.mpirun=4.output [moved from tests/mpi/error_predictor_01.with_p4est=true.mpirun=4.output with 100% similarity]
tests/mpi/error_prediction_02.cc [moved from tests/mpi/error_predictor_02.cc with 80% similarity]
tests/mpi/error_prediction_02.with_p4est=true.mpirun=1.output [moved from tests/mpi/error_predictor_02.with_p4est=true.mpirun=1.output with 100% similarity]
tests/mpi/error_prediction_02.with_p4est=true.mpirun=4.output [moved from tests/mpi/error_predictor_02.with_p4est=true.mpirun=4.output with 100% similarity]

diff --git a/doc/news/changes/incompatibilities/20210421Fehling b/doc/news/changes/incompatibilities/20210421Fehling
new file mode 100644 (file)
index 0000000..e339c8a
--- /dev/null
@@ -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.
+<br>
+(Marc Fehling, 2021/04/21)
index 0b7be5b5c06c24778ba18cd2e272ce2d70376d2a..05726c7b842015007bbc4fb6855bbea708c5baf8 100644 (file)
@@ -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<float> predicted_errors;
+     * triangulation.signals.post_p4est_refinement.connect([&]() {
+     *   const internal::parallel::distributed::TemporarilyMatchRefineFlags<dim>
+     *     refine_modifier(triangulation);
+     *   predicted_errors.reinit(triangulation.n_active_cells());
+     *   hp::Refinement::predict_error(dof_handler,
+     *                                 error_indicators,
+     *                                 predicted_errors);
+     * });
+     * @endcode
+     * The container <code>predicted_errors</code> 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
similarity index 76%
rename from tests/mpi/error_predictor_01.cc
rename to tests/mpi/error_prediction_01.cc
index f531365914551a6f4628f02c53f1aa06b94c9658..4b692b25ccd3463a9c45417d4b41ee584e696f0a 100644 (file)
 // This tests is based on hp/error_prediction.cc
 
 
-#include <deal.II/distributed/error_predictor.h>
+#include <deal.II/distributed/cell_data_transfer.h>
 #include <deal.II/distributed/tria.h>
 
 #include <deal.II/dofs/dof_handler.h>
-#include <deal.II/dofs/dof_tools.h>
 
 #include <deal.II/fe/fe_q.h>
 
+#include <deal.II/hp/refinement.h>
+
 #include <deal.II/lac/vector.h>
 
 #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<float> predicted_errors;
+  tria.signals.post_p4est_refinement.connect([&]() {
+    const internal::parallel::distributed::TemporarilyMatchRefineFlags<dim>
+      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<dim> predictor(dh);
+  parallel::distributed::CellDataTransfer<dim, dim, Vector<float>>
+  data_transfer(tria,
+                /*transfer_variable_size_data=*/false,
+                &AdaptationStrategies::Refinement::l2_norm<dim, dim, float>,
+                &AdaptationStrategies::Coarsening::l2_norm<dim, dim, float>);
 
-  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<float> 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;
similarity index 80%
rename from tests/mpi/error_predictor_02.cc
rename to tests/mpi/error_prediction_02.cc
index 3de8ff0bbd3d776a282fce286ffc0d5b58b49018..e999db84ce923eb3cfd886a30bb0d4f74e8f1183 100644 (file)
 // This tests is based on hp/error_prediction.cc
 
 
-#include <deal.II/distributed/error_predictor.h>
+#include <deal.II/distributed/cell_data_transfer.h>
 #include <deal.II/distributed/tria.h>
 
 #include <deal.II/dofs/dof_handler.h>
-#include <deal.II/dofs/dof_tools.h>
 
 #include <deal.II/fe/fe_q.h>
 
+#include <deal.II/hp/refinement.h>
+
 #include <deal.II/lac/vector.h>
 
 #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<float> predicted_errors;
+  tria.signals.post_p4est_refinement.connect([&]() {
+    const internal::parallel::distributed::TemporarilyMatchRefineFlags<dim>
+      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<dim> predictor(dh);
+  parallel::distributed::CellDataTransfer<dim, dim, Vector<float>>
+  data_transfer(tria,
+                /*transfer_variable_size_data=*/false,
+                &AdaptationStrategies::Refinement::l2_norm<dim, dim, float>,
+                &AdaptationStrategies::Coarsening::l2_norm<dim, dim, float>);
 
-  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<float> 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;

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.