From 47f9e495f780ec1354a9861222cf15325dccfa25 Mon Sep 17 00:00:00 2001 From: Marco Feder Date: Tue, 19 Sep 2023 11:12:05 +0000 Subject: [PATCH] Testing callbacks for non nested mg --- .../multigrid/mg_transfer_global_coarsening.h | 32 ++-- .../events_non_nested.cc | 155 ++++++++++++++++++ .../events_non_nested.output | 20 +++ 3 files changed, 191 insertions(+), 16 deletions(-) create mode 100644 tests/multigrid-global-coarsening/events_non_nested.cc create mode 100644 tests/multigrid-global-coarsening/events_non_nested.output diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.h index 3b569ce551..99fc5835e2 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.h @@ -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 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 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 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 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 &slot); /** * Connect a function to mg::SignalsNonNested::evaluation_restriction. - * */ boost::signals2::connection connect_evaluation_restriction(const std::function &slot); /** * Connect a function to mg::SignalsNonNested::evaluate_and_process. - * */ boost::signals2::connection connect_evaluate_and_process(const std::function &slot); /** * Connect a function to mg::SignalsNonNested::process_and_evaluate. - * */ boost::signals2::connection connect_process_and_evaluate(const std::function &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 index 0000000000..723a5a593d --- /dev/null +++ b/tests/multigrid-global-coarsening/events_non_nested.cc @@ -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 + +#include "multigrid_util.h" + + +const auto print_evaluate_and_process = [](const bool start) { + deallog << "evaluate_and_process " << (start ? "started" : "finished") + << std::endl; +}; + +template +void +test(const unsigned int n_refinements, + const unsigned int fe_degree_fine, + const bool do_simplex_mesh) +{ + using VectorType = LinearAlgebra::distributed::Vector; + + const unsigned int min_level = 0; + const unsigned int max_level = n_refinements; + + MGLevelObject> triangulations(min_level, max_level); + MGLevelObject> dof_handlers(min_level, max_level); + MGLevelObject> constraints(min_level, max_level); + MGLevelObject>> mappings(min_level, max_level); + MGLevelObject>> + transfers(min_level, max_level); + MGLevelObject> 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> fe; + std::unique_ptr> quad; + std::unique_ptr> mapping; + + if (do_simplex_mesh) + { + fe = std::make_unique>(fe_degree_fine); + quad = std::make_unique>(fe_degree_fine + 1); + mapping = std::make_unique>(FE_SimplexP(1)); + } + else + { + fe = std::make_unique>(fe_degree_fine); + quad = std::make_unique>(fe_degree_fine + 1); + mapping = std::make_unique>(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(), 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>(); + 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 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 index 0000000000..1a294c1e0a --- /dev/null +++ b/tests/multigrid-global-coarsening/events_non_nested.output @@ -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 -- 2.39.5