]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add signals for non nested mg
authorMarco Feder <marco.feder@sissa.it>
Tue, 19 Sep 2023 11:10:56 +0000 (11:10 +0000)
committerMarco Feder <marco.feder@sissa.it>
Tue, 19 Sep 2023 11:10:56 +0000 (11:10 +0000)
include/deal.II/multigrid/mg_transfer_global_coarsening.h
include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h

index 222cb12923560eb87016300b6ec16e3e36124f4c..3b569ce5511da06c3d99675c1fb51b4bc90f816a 100644 (file)
@@ -52,6 +52,49 @@ namespace RepartitioningPolicyTools
 #endif
 
 
+namespace mg
+{
+  /**
+   * A structure with boost::signal objects for optional processing in a
+   * non-nested multigrid solver.
+   *
+   * Similarly to mg::Signals, each signal is called twice: once before and once
+   * after the action is performed. The two calls only differ in the boolean argument @p before, which is true the first time and false the second.
+   *
+   */
+  struct SignalsNonNested
+  {
+    /**
+     * This signal is triggered before and after the call to actual evaluation
+     * function inside RemotePointEvaluation::evaluate_and_process() during
+     * prolongation.
+     */
+    boost::signals2::signal<void(const bool before)> evaluation_prolongation;
+
+    /**
+     * This signal is triggered before and after the call to actual evaluation
+     * function inside RemotePointEvaluation::process_and_evaluate() during
+     * restriction.
+     */
+    boost::signals2::signal<void(const bool before)> evaluation_restriction;
+
+    /**
+     * This signal is triggered before and after the call to
+     * RemotePointEvaluation::evaluate_and_process() used in
+     * MGTwoLevelTransferNonNested::prolongate_and_add(). The difference with the @p evaluation
+     * signal is that also the communication phase is included.
+     */
+    boost::signals2::signal<void(const bool before)> evaluate_and_process;
+
+    /**
+     * This signal is triggered before and after the call to
+     * RemotePointEvaluation::process_and_evaluate() used in
+     * MGTwoLevelTransferNonNested::restrict_and_add(). Similarly to the @p evaluate_and_process signal,
+     * also the communication phase is included.
+     */
+    boost::signals2::signal<void(const bool before)> process_and_evaluate;
+  };
+} // namespace mg
 
 /**
  * Global coarsening utility functions.
@@ -705,6 +748,8 @@ class MGTwoLevelTransferNonNested<dim,
 private:
   using VectorizedArrayType = VectorizedArray<Number, 1>;
 
+  mg::SignalsNonNested signals_non_nested;
+
 public:
   /**
    * AdditionalData structure that can be used to tweak parameters
@@ -798,6 +843,34 @@ public:
   std::size_t
   memory_consumption() const override;
 
+  /**
+   * Connect a function to mg::SignalsNonNested::evaluation_restriction.
+   *
+   */
+  boost::signals2::connection
+  connect_evaluation_prolongation(const std::function<void(const bool)> &slot);
+
+  /**
+   * Connect a function to mg::SignalsNonNested::evaluation_restriction.
+   *
+   */
+  boost::signals2::connection
+  connect_evaluation_restriction(const std::function<void(const bool)> &slot);
+
+  /**
+   * Connect a function to mg::SignalsNonNested::evaluate_and_process.
+   *
+   */
+  boost::signals2::connection
+  connect_evaluate_and_process(const std::function<void(const bool)> &slot);
+
+  /**
+   * Connect a function to mg::SignalsNonNested::process_and_evaluate.
+   *
+   */
+  boost::signals2::connection
+  connect_process_and_evaluate(const std::function<void(const bool)> &slot);
+
 protected:
   AdditionalData additional_data;
   /**
index 55d81d247c30732a438cd1aa606311876cb2b455..b00bcfeea8189794a6af7503c7ba38c112c5fba7 100644 (file)
@@ -4942,6 +4942,7 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
   std::vector<value_type> buffer;
 
   const auto evaluation_function = [&](auto &values, const auto &cell_data) {
+    this->signals_non_nested.evaluation_prolongation(true);
     std::vector<Number> solution_values;
 
     FEPointEvaluation<n_components, dim, dim, Number> evaluator(*mapping_info,
@@ -4971,11 +4972,14 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
           values[q + cell_data.reference_point_ptrs[cell]] =
             evaluator.get_value(q);
       }
+    this->signals_non_nested.evaluation_prolongation(false);
   };
 
+  this->signals_non_nested.evaluate_and_process(true);
   rpe.template evaluate_and_process<value_type>(evaluation_point_results,
                                                 buffer,
                                                 evaluation_function);
+  this->signals_non_nested.evaluate_and_process(false);
 
   // Weight operator in case some points are owned by multiple cells.
   if (rpe.is_map_unique() == false)
@@ -5117,6 +5121,7 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
 
   const auto evaluation_function = [&](const auto &values,
                                        const auto &cell_data) {
+    this->signals_non_nested.evaluation_restriction(true);
     std::vector<Number>                               solution_values;
     FEPointEvaluation<n_components, dim, dim, Number> evaluator(*mapping_info,
                                                                 *fe_coarse);
@@ -5146,11 +5151,14 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
           solution_values.size(),
           true);
       }
+    this->signals_non_nested.evaluation_restriction(false);
   };
 
+  this->signals_non_nested.process_and_evaluate(true);
   rpe.template process_and_evaluate<value_type>(evaluation_point_results,
                                                 buffer,
                                                 evaluation_function);
+  this->signals_non_nested.process_and_evaluate(false);
 }
 
 
@@ -5219,6 +5227,46 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
 }
 
 
+
+template <int dim, typename Number>
+boost::signals2::connection
+MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
+  connect_evaluation_prolongation(const std::function<void(const bool)> &slot)
+{
+  return signals_non_nested.evaluation_prolongation.connect(slot);
+}
+
+
+
+template <int dim, typename Number>
+boost::signals2::connection
+MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
+  connect_evaluation_restriction(const std::function<void(const bool)> &slot)
+{
+  return signals_non_nested.evaluation_restriction.connect(slot);
+}
+
+
+
+template <int dim, typename Number>
+boost::signals2::connection
+MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
+  connect_evaluate_and_process(const std::function<void(const bool)> &slot)
+{
+  return signals_non_nested.evaluate_and_process.connect(slot);
+}
+
+
+
+template <int dim, typename Number>
+boost::signals2::connection
+MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
+  connect_process_and_evaluate(const std::function<void(const bool)> &slot)
+{
+  return signals_non_nested.process_and_evaluate.connect(slot);
+}
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif

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.