]> https://gitweb.dealii.org/ - dealii.git/commitdiff
improve documentation of AdditionaData 15899/head
authorMarco Feder <marco.feder@sissa.it>
Tue, 22 Aug 2023 08:53:36 +0000 (08:53 +0000)
committerMarco Feder <marco.feder@sissa.it>
Wed, 23 Aug 2023 22:29:46 +0000 (22:29 +0000)
include/deal.II/multigrid/mg_transfer_global_coarsening.h
include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h

index 756e75084b6ca9b5c9f3bde12c32058db042f42f..5b66ffc4dc9425434fff8e54c488bfdd63ba97b1 100644 (file)
@@ -707,30 +707,51 @@ private:
 
 public:
   /**
-   * AdditionalData structure with the arguments needed by
-   * RemotePointEvaluation. Default values are the same as the ones described in
-   * the documentation of RemotePointEvaluation. The last boolean parameter, @p enf_all_points_found is true by defaults and
-   * checks if RemotePointEvaluation::all_points_found() evaluates to true, i.e.
-   * all submitted points have been found inside the domain.
+   * AdditionalData structure that can be used to tweak parameters
+   * related to the search procedure (used internally by RemotePointEvaluation)
+   * or, in the future, transfer operators needed by the non-nested multigrid
+   * algorithm.
    */
   struct AdditionalData
   {
-    AdditionalData(const double       tol                = 1e-6,
-                   const bool         enf_unique_mapping = false,
-                   const unsigned int rtree_l            = 0,
-                   const std::function<std::vector<bool>()> &marked_verts = {},
-                   const bool enf_all_points_found = true)
-      : tolerance(tol)
-      , enforce_unique_mapping(enf_unique_mapping)
-      , rtree_level(rtree_l)
-      , marked_vertices(marked_verts)
-      , enforce_all_points_found(enf_all_points_found)
+    /**
+     * Constructor. By default, the @p tolerance and @p rtree_level parameters
+     * are set to the default values used in the constructor of
+     * RemotePointEvaluation, i.e. 1e-6 and 0, respectively. The last Boolean
+     * parameter @p enforce_all_points_found is true by default and checks
+     * that all points submitted internally to RemotePointEvaluation::reinit()
+     * have been found.
+     *
+     */
+    AdditionalData(const double       tolerance                = 1e-6,
+                   const unsigned int rtree_level              = 0,
+                   const bool         enforce_all_points_found = true)
+      : tolerance(tolerance)
+      , rtree_level(rtree_level)
+      , enforce_all_points_found(enforce_all_points_found)
     {}
-    double                             tolerance;
-    bool                               enforce_unique_mapping;
-    unsigned int                       rtree_level;
-    std::function<std::vector<bool>()> marked_vertices;
-    bool                               enforce_all_points_found;
+
+    /**
+     * Tolerance parameter. See the constructor of RemotePointEvaluation for
+     * more details.
+     */
+    double tolerance;
+
+    /**
+     * RTree level parameter. See the constructor of RemotePointEvaluation for
+     * more details.
+     *
+     */
+    unsigned int rtree_level;
+
+    /**
+     * If set to true, it checks if RemotePointEvaluation::all_points_found()
+     * evaluates to true internally during the each call to reinit() from one
+     * level to the next one, ensuring that all submitted points have been found
+     * inside the domain.
+     *
+     */
+    bool enforce_all_points_found;
   };
 
   MGTwoLevelTransferNonNested(const AdditionalData &data = AdditionalData());
@@ -818,7 +839,7 @@ private:
    * Object to evaluate shape functions on one mesh on visited support points of
    * the other mesh.
    */
-  std::shared_ptr<Utilities::MPI::RemotePointEvaluation<dim>> rpe;
+  Utilities::MPI::RemotePointEvaluation<dim> rpe;
 
   /**
    * MappingInfo object needed as Mapping argument by FEPointEvaluation.
index 1de0ac783081ca7f9d02152981d927e02d0ce186..9e23d175a44f5da441f592d2b60e64d7aaf76197 100644 (file)
@@ -4105,13 +4105,8 @@ template <int dim, typename Number>
 MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
   MGTwoLevelTransferNonNested(const AdditionalData &data)
   : additional_data(data)
-{
-  rpe = std::make_shared<Utilities::MPI::RemotePointEvaluation<dim>>(
-    data.tolerance,
-    data.enforce_unique_mapping,
-    data.rtree_level,
-    data.marked_vertices);
-}
+  , rpe(data.tolerance, false, data.rtree_level, {})
+{}
 
 template <int dim, typename Number>
 void
@@ -4183,19 +4178,19 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
       this->partitioner_fine->global_to_local(global_dof_indices[i]);
 
   // hand points over to RPE
-  rpe->reinit(points, dof_handler_coarse.get_triangulation(), mapping_coarse);
+  rpe.reinit(points, dof_handler_coarse.get_triangulation(), mapping_coarse);
 
-  if (additional_data.enforce_all_points_found)
-    AssertThrow(
-      rpe->all_points_found(),
-      ExcMessage(
-        "You requested that all points should be found, but this didn't happen. You can change this option through the AdditionaData struct in the constructor."));
+  AssertThrow(
+    !additional_data.enforce_all_points_found || rpe.all_points_found(),
+    ExcMessage(
+      "You requested that all points should be found, but this didn'thappen."
+      " You can change this option through the AdditionaData struct in the constructor."));
 
   // set up MappingInfo for easier data access
-  mapping_info = internal::fill_mapping_info<dim, Number>(*rpe);
+  mapping_info = internal::fill_mapping_info<dim, Number>(rpe);
 
   // set up constraints
-  const auto &cell_data = rpe->get_cell_data();
+  const auto &cell_data = rpe.get_cell_data();
 
   constraint_info.reinit(dof_handler_coarse,
                          cell_data.cells.size(),
@@ -4204,7 +4199,7 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
   for (unsigned int i = 0; i < cell_data.cells.size(); ++i)
     {
       typename DoFHandler<dim>::active_cell_iterator cell(
-        &rpe->get_triangulation(),
+        &rpe.get_triangulation(),
         cell_data.cells[i].first,
         cell_data.cells[i].second,
         &dof_handler_coarse);
@@ -4314,18 +4309,18 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
       }
   };
 
-  rpe->template evaluate_and_process<value_type>(evaluation_point_results,
-                                                 buffer,
-                                                 evaluation_function);
+  rpe.template evaluate_and_process<value_type>(evaluation_point_results,
+                                                buffer,
+                                                evaluation_function);
 
   // Weight operator in case some points are owned by multiple cells.
-  if (rpe->is_map_unique() == false)
+  if (rpe.is_map_unique() == false)
     {
       const auto evaluation_point_results_temp = evaluation_point_results;
       evaluation_point_results.clear();
-      evaluation_point_results.reserve(rpe->get_point_ptrs().size() - 1);
+      evaluation_point_results.reserve(rpe.get_point_ptrs().size() - 1);
 
-      const auto &ptr = rpe->get_point_ptrs();
+      const auto &ptr = rpe.get_point_ptrs();
 
       for (unsigned int i = 0; i < ptr.size() - 1; ++i)
         {
@@ -4407,7 +4402,7 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
   std::vector<value_type> evaluation_point_results;
   std::vector<value_type> buffer;
 
-  evaluation_point_results.resize(rpe->get_point_ptrs().size() - 1);
+  evaluation_point_results.resize(rpe.get_point_ptrs().size() - 1);
 
   for (unsigned int j = 0; j < evaluation_point_results.size(); ++j)
     {
@@ -4442,9 +4437,9 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
     }
 
   // Weight operator in case some points are owned by multiple cells.
-  if (rpe->is_map_unique() == false)
+  if (rpe.is_map_unique() == false)
     {
-      const auto &ptr = rpe->get_point_ptrs();
+      const auto &ptr = rpe.get_point_ptrs();
 
       for (unsigned int i = 0; i < ptr.size() - 1; ++i)
         {
@@ -4489,9 +4484,9 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
       }
   };
 
-  rpe->template process_and_evaluate<value_type>(evaluation_point_results,
-                                                 buffer,
-                                                 evaluation_function);
+  rpe.template process_and_evaluate<value_type>(evaluation_point_results,
+                                                buffer,
+                                                evaluation_function);
 }
 
 

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.