]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Testing callbacks for non nested mg
authorMarco Feder <marco.feder@sissa.it>
Tue, 19 Sep 2023 11:12:05 +0000 (11:12 +0000)
committerMarco Feder <marco.feder@sissa.it>
Wed, 20 Sep 2023 16:38:46 +0000 (16:38 +0000)
include/deal.II/multigrid/mg_transfer_global_coarsening.h
tests/multigrid-global-coarsening/events_non_nested.cc [new file with mode: 0644]
tests/multigrid-global-coarsening/events_non_nested.output [new file with mode: 0644]

index 3b569ce5511da06c3d99675c1fb51b4bc90f816a..99fc5835e2f3eff31a79751febc58641dd33e273 100644 (file)
@@ -58,39 +58,43 @@ 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.
+   * 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
+   * booleanargument @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.
+     * This signal is triggered before and after the call to the actual
+     * evaluationfunction 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.
+     * This signal is triggered before and after the call to the 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.
+     * MGTwoLevelTransferNonNested::prolongate_and_add(). The difference
+     * with the @p evaluation_prolongation 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.
+     * MGTwoLevelTransferNonNested::restrict_and_add(). Similarly to
+     * the @p evaluation_restriction signal, also the communication phase is
+     * included.
      */
     boost::signals2::signal<void(const bool before)> process_and_evaluate;
   };
@@ -845,28 +849,24 @@ public:
 
   /**
    * 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);
diff --git a/tests/multigrid-global-coarsening/events_non_nested.cc b/tests/multigrid-global-coarsening/events_non_nested.cc
new file mode 100644 (file)
index 0000000..723a5a5
--- /dev/null
@@ -0,0 +1,155 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 - 2023 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+/**
+ * Test signals for prolongation and restriction within the non-nested mg
+ * infrastructure.
+ */
+
+#include <deal.II/grid/grid_tools.h>
+
+#include "multigrid_util.h"
+
+
+const auto print_evaluate_and_process = [](const bool start) {
+  deallog << "evaluate_and_process " << (start ? "started" : "finished")
+          << std::endl;
+};
+
+template <int dim, typename Number = double>
+void
+test(const unsigned int n_refinements,
+     const unsigned int fe_degree_fine,
+     const bool         do_simplex_mesh)
+{
+  using VectorType = LinearAlgebra::distributed::Vector<Number>;
+
+  const unsigned int min_level = 0;
+  const unsigned int max_level = n_refinements;
+
+  MGLevelObject<Triangulation<dim>>        triangulations(min_level, max_level);
+  MGLevelObject<DoFHandler<dim>>           dof_handlers(min_level, max_level);
+  MGLevelObject<AffineConstraints<Number>> constraints(min_level, max_level);
+  MGLevelObject<std::unique_ptr<Mapping<dim>>> mappings(min_level, max_level);
+  MGLevelObject<std::shared_ptr<MGTwoLevelTransferNonNested<dim, VectorType>>>
+                                       transfers(min_level, max_level);
+  MGLevelObject<Operator<dim, Number>> operators(min_level, max_level);
+
+
+  // set up levels
+  for (auto l = min_level; l <= max_level; ++l)
+    {
+      auto &tria        = triangulations[l];
+      auto &dof_handler = dof_handlers[l];
+      auto &constraint  = constraints[l];
+      auto &op          = operators[l];
+
+      std::unique_ptr<FiniteElement<dim>> fe;
+      std::unique_ptr<Quadrature<dim>>    quad;
+      std::unique_ptr<Mapping<dim>>       mapping;
+
+      if (do_simplex_mesh)
+        {
+          fe      = std::make_unique<FE_SimplexP<dim>>(fe_degree_fine);
+          quad    = std::make_unique<QGaussSimplex<dim>>(fe_degree_fine + 1);
+          mapping = std::make_unique<MappingFE<dim>>(FE_SimplexP<dim>(1));
+        }
+      else
+        {
+          fe      = std::make_unique<FE_Q<dim>>(fe_degree_fine);
+          quad    = std::make_unique<QGauss<dim>>(fe_degree_fine + 1);
+          mapping = std::make_unique<MappingQ<dim>>(1);
+        }
+
+      mappings[l] = mapping->clone();
+
+      // set up triangulation
+      if (do_simplex_mesh)
+        GridGenerator::subdivided_hyper_cube_with_simplices(tria, 2);
+      else
+        GridGenerator::subdivided_hyper_cube(tria, 2);
+      tria.refine_global(l);
+
+      // set up dofhandler
+      dof_handler.reinit(tria);
+      dof_handler.distribute_dofs(*fe);
+
+      // set up constraints
+      const IndexSet locally_relevant_dofs =
+        DoFTools::extract_locally_relevant_dofs(dof_handler);
+      constraint.reinit(locally_relevant_dofs);
+      VectorTools::interpolate_boundary_values(
+        *mapping, dof_handler, 0, Functions::ZeroFunction<dim>(), constraint);
+      constraint.close();
+
+      // set up operator
+      op.reinit(*mapping, dof_handler, *quad, constraint);
+    }
+
+  // set up transfer operator and test connections
+  for (unsigned int l = min_level; l < max_level; ++l)
+    {
+      transfers[l + 1] =
+        std::make_shared<MGTwoLevelTransferNonNested<dim, VectorType>>();
+      transfers[l + 1]->connect_evaluate_and_process(
+        print_evaluate_and_process);
+      transfers[l + 1]->reinit(dof_handlers[l + 1],
+                               dof_handlers[l],
+                               *mappings[l + 1],
+                               *mappings[l],
+                               constraints[l + 1],
+                               constraints[l]);
+    }
+
+  MGTransferGlobalCoarsening<dim, VectorType> transfer(
+    transfers,
+    [&](const auto l, auto &vec) { operators[l].initialize_dof_vector(vec); });
+
+
+  GMGParameters mg_data; // TODO
+
+  VectorType dst, src;
+  operators[max_level].initialize_dof_vector(dst);
+  operators[max_level].initialize_dof_vector(src);
+
+  operators[max_level].rhs(src);
+
+  ReductionControl solver_control(
+    mg_data.maxiter, mg_data.abstol, mg_data.reltol, false, false);
+
+  mg_solve(solver_control,
+           dst,
+           src,
+           mg_data,
+           dof_handlers[max_level],
+           operators[max_level],
+           operators,
+           transfer);
+
+  deallog << dim << ' ' << fe_degree_fine << ' ' << n_refinements << ' '
+          << (do_simplex_mesh ? "tri " : "quad") << ' '
+          << solver_control.last_step() << std::endl;
+}
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    all;
+
+  deallog.precision(8);
+  test<2>(3, 2, false);
+}
diff --git a/tests/multigrid-global-coarsening/events_non_nested.output b/tests/multigrid-global-coarsening/events_non_nested.output
new file mode 100644 (file)
index 0000000..1a294c1
--- /dev/null
@@ -0,0 +1,20 @@
+
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0:cg::evaluate_and_process started
+DEAL:0:cg::evaluate_and_process finished
+DEAL:0::2 2 3 quad 3

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.