From b96849a7044b3ca1ddf0e333a825635133c0b1e2 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Mon, 5 Feb 2018 12:20:16 -0500 Subject: [PATCH] Implement signals for Multigrid --- .../minor/20180331DanielArndtTimoHeister | 11 + include/deal.II/multigrid/multigrid.h | 179 +++++- .../deal.II/multigrid/multigrid.templates.h | 63 ++ .../parallel_multigrid_adaptive_01.cc | 4 +- .../parallel_multigrid_adaptive_02.cc | 4 +- .../parallel_multigrid_adaptive_04.cc | 4 +- .../parallel_multigrid_adaptive_06.cc | 4 +- .../parallel_multigrid_adaptive_06ref.cc | 4 +- tests/multigrid/events_01.cc | 590 ++++++++++++++++++ ...st=true.with_trilinos=true.mpirun=1.output | 82 +++ tests/multigrid/transfer_matrix_free_09.cc | 6 +- tests/multigrid/transfer_matrix_free_10.cc | 6 +- 12 files changed, 933 insertions(+), 24 deletions(-) create mode 100644 doc/news/changes/minor/20180331DanielArndtTimoHeister create mode 100644 tests/multigrid/events_01.cc create mode 100644 tests/multigrid/events_01.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=1.output diff --git a/doc/news/changes/minor/20180331DanielArndtTimoHeister b/doc/news/changes/minor/20180331DanielArndtTimoHeister new file mode 100644 index 0000000000..871060ce1e --- /dev/null +++ b/doc/news/changes/minor/20180331DanielArndtTimoHeister @@ -0,0 +1,11 @@ +New: Multigrid obtained the new signals +Multigrid::Signals::transfer_to_mg, +Multigrid::Signals::transfer_to_global, +Multigrid::Signals::coarse_solve, +Multigrid::Signals::transfer_down_to, +Multigrid::Signals::transfer_up_to, +Multigrid::Signals::pre_smoother_step, +Multigrid::Signals::post_smoother_step +that functions can be connected to. +
+(Daniel Arndt, Timo Heister, 2018/03/31) diff --git a/include/deal.II/multigrid/multigrid.h b/include/deal.II/multigrid/multigrid.h index b9fba95249..dfc811ae2d 100644 --- a/include/deal.II/multigrid/multigrid.h +++ b/include/deal.II/multigrid/multigrid.h @@ -74,6 +74,64 @@ public: f_cycle }; + /** + * A structure that has boost::signal objects for a number of actions that + * happen for the typical usage of this class. + * + * For documentation on signals, see + * http://www.boost.org/doc/libs/release/libs/signals2 . + */ + struct Signals + { + /** + * This signal is triggered before (@p start is true) and after (@p start is + * false) the call to MGTransfer::copy_to_mg which transfers the vector + * given to it to a multi-level vector. + */ + boost::signals2::signal transfer_to_mg; + + /** + * This signal is triggered before (@p start is true) and after (@p start is + * false) the call to MGTransfer::copy_from_mg which transfers the + * multi-level vector given to it to a normal vector. + */ + boost::signals2::signal transfer_to_global; + + /** + * This signal is triggered before (@p start is true) and after (@p start is + * false) the call to the coarse solver on @p level. + */ + boost::signals2::signal coarse_solve; + + /** + * This signal is triggered before (@p start is true) and after (@p start is + * false) the call to MGTransfer::restrict_and_add() which restricts a + * vector from @p level to the next coarser one (@p level - 1). + */ + boost::signals2::signal transfer_down_to; + + /** + * This signal is triggered before (@p start is true) and after (@p start is + * false) the call to MGTransfer::prolongate() which prolongs a vector to + * @p level from the next coarser one (@p level - 1). + */ + boost::signals2::signal transfer_up_to; + + /** + * This signal is triggered before (@p start is true) and after (@p start is + * false) the call to a pre-smoothing step via MGPreSmoother::apply() on + * @p level. + */ + boost::signals2::signal pre_smoother_step; + + /** + * This signal is triggered before (@p start is true) and after (@p start is + * false) the call to a post-smoothing step via MGPostSmoother::apply() on + * @p level. + */ + boost::signals2::signal post_smoother_step; + }; + typedef VectorType vector_type; typedef const VectorType const_vector_type; @@ -224,8 +282,43 @@ public: */ void set_debug (const unsigned int); + /** + * Connect a function to Signals::coarse_solve. + */ + boost::signals2::connection + connect_coarse_solve(const std::function &slot); + + /** + * Connect a function to Signals::transfer_down_to. + */ + boost::signals2::connection + connect_transfer_down_to(const std::function &slot); + + /** + * Connect a function to Signals::transfer_up_to. + */ + boost::signals2::connection + connect_transfer_up_to(const std::function &slot); + + /** + * Connect a function to Signals::pre_smoother_step. + */ + boost::signals2::connection + connect_pre_smoother_step(const std::function &slot); + + /** + * Connect a function to Signals::post_smoother_step. + */ + boost::signals2::connection + connect_post_smoother_step(const std::function &slot); + private: + /** + * Signals for the various actions that the Multigrid algorithm uses. + */ + Signals signals; + /** * The V-cycle multigrid method. level is the level the function * starts on. It will usually be called for the highest level from outside, @@ -405,7 +498,7 @@ public: const OtherVectorType &src) const; /** - * Tranposed preconditioning operator. + * Transposed preconditioning operator. * * Not implemented, but the definition may be needed. */ @@ -414,7 +507,7 @@ public: const OtherVectorType &src) const; /** - * Tranposed preconditioning operator. + * Transposed preconditioning operator. * * Not implemented, but the definition may be needed. */ @@ -443,6 +536,18 @@ public: */ MPI_Comm get_mpi_communicator() const; + /** + * Connect a function to Signals::transfer_to_mg. + */ + boost::signals2::connection + connect_transfer_to_mg(const std::function &slot); + + /** + * Connect a function to Signals::transfer_to_global. + */ + boost::signals2::connection + connect_transfer_to_global(const std::function &slot); + private: /** * Associated @p DoFHandler. @@ -470,6 +575,11 @@ private: * or with one for each block. */ const bool uses_dof_handler_vector; + + /** + * Signals used by this object + */ + typename Multigrid::Signals signals; }; /*@}*/ @@ -577,8 +687,10 @@ namespace internal const TRANSFER &transfer, OtherVectorType &dst, const OtherVectorType &src, - const bool uses_dof_handler_vector, int) + const bool uses_dof_handler_vector, + const typename dealii::Multigrid::Signals &signals, int) { + signals.transfer_to_mg(true); if (uses_dof_handler_vector) transfer.copy_to_mg(dof_handler_vector, multigrid.defect, @@ -587,8 +699,11 @@ namespace internal transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src); + signals.transfer_to_mg(false); multigrid.cycle(); + + signals.transfer_to_global(true); if (uses_dof_handler_vector) transfer.copy_from_mg(dof_handler_vector, dst, @@ -597,6 +712,7 @@ namespace internal transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution); + signals.transfer_to_global(false); } template @@ -606,17 +722,25 @@ namespace internal const TRANSFER &transfer, OtherVectorType &dst, const OtherVectorType &src, - const bool uses_dof_handler_vector,...) + const bool uses_dof_handler_vector, + const typename dealii::Multigrid::Signals &signals, ...) { (void) uses_dof_handler_vector; Assert (!uses_dof_handler_vector, ExcInternalError()); + + signals.transfer_to_mg(true); transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src); + signals.transfer_to_mg(false); + multigrid.cycle(); + + signals.transfer_to_global(true); transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution); + signals.transfer_to_global(false); } template @@ -626,8 +750,10 @@ namespace internal const TRANSFER &transfer, OtherVectorType &dst, const OtherVectorType &src, - const bool uses_dof_handler_vector, int) + const bool uses_dof_handler_vector, + const typename dealii::Multigrid::Signals &signals, int) { + signals.transfer_to_mg(0, true); if (uses_dof_handler_vector) transfer.copy_to_mg(dof_handler_vector, multigrid.defect, @@ -636,8 +762,11 @@ namespace internal transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src); + signals.transfer_to_mg(0, false); multigrid.cycle(); + + signals.transfer_to_global(0, true); if (uses_dof_handler_vector) transfer.copy_from_mg_add(dof_handler_vector, dst, @@ -646,6 +775,7 @@ namespace internal transfer.copy_from_mg_add(*dof_handler_vector[0], dst, multigrid.solution); + signals.transfer_to_global(0, false); } template @@ -655,17 +785,25 @@ namespace internal const TRANSFER &transfer, OtherVectorType &dst, const OtherVectorType &src, - const bool uses_dof_handler_vector,...) + const bool uses_dof_handler_vector, + const typename dealii::Multigrid::Signals &signals, ...) { (void) uses_dof_handler_vector; Assert (!uses_dof_handler_vector, ExcInternalError()); + + signals.transfer_to_mg(true); transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src); + signals.transfer_to_mg(false); + multigrid.cycle(); + + signals.transfer_to_global(true); transfer.copy_from_mg_add(*dof_handler_vector[0], dst, multigrid.solution); + signals.transfer_to_global(false); } } } @@ -717,7 +855,8 @@ PreconditionMG::vmult const OtherVectorType &src) const { internal::PreconditionMGImplementation::vmult(dof_handler_vector_raw,*multigrid,*transfer, - dst,src,uses_dof_handler_vector,0); + dst,src,uses_dof_handler_vector, + this->signals, 0); } @@ -739,6 +878,7 @@ PreconditionMG::locally_owned_domain_indices(const un } + template MPI_Comm PreconditionMG::get_mpi_communicator() const @@ -752,6 +892,28 @@ PreconditionMG::get_mpi_communicator() const } + +template +boost::signals2::connection +PreconditionMG:: +connect_transfer_to_mg(const std::function &slot) +{ + return this->signals.transfer_to_mg.connect(slot); +} + + + +template +boost::signals2::connection +PreconditionMG:: +connect_transfer_to_global( + const std::function &slot) +{ + return this->signals.transfer_to_global.connect(slot); +} + + + template template void @@ -760,7 +922,8 @@ PreconditionMG::vmult_add const OtherVectorType &src) const { internal::PreconditionMGImplementation::vmult_add(dof_handler_vector_raw,*multigrid,*transfer, - dst,src,uses_dof_handler_vector,0); + dst,src,uses_dof_handler_vector, + this->signals, 0); } diff --git a/include/deal.II/multigrid/multigrid.templates.h b/include/deal.II/multigrid/multigrid.templates.h index 77986c0587..3d384e5fd4 100644 --- a/include/deal.II/multigrid/multigrid.templates.h +++ b/include/deal.II/multigrid/multigrid.templates.h @@ -116,16 +116,21 @@ Multigrid::level_v_step (const unsigned int level) if (level == minlevel) { + this->signals.coarse_solve(true, level); if (debug>0) deallog << "Coarse level " << level << std::endl; (*coarse)(level, solution[level], defect[level]); + this->signals.coarse_solve(false, level); return; } // smoothing of the residual if (debug>1) deallog << "Smoothing on level " << level << std::endl; + + this->signals.pre_smoother_step(true, level); pre_smooth->apply(level, solution[level], defect[level]); + this->signals.pre_smoother_step(false, level); if (debug>2) deallog << "Solution norm " << solution[level].l2_norm() @@ -154,13 +159,19 @@ Multigrid::level_v_step (const unsigned int level) edge_down->vmult(level, t[level-1], solution[level]); defect[level-1] -= t[level-1]; } + + this->signals.transfer_down_to(true, level); transfer->restrict_and_add(level, defect[level-1], t[level]); + this->signals.transfer_down_to(false, level); // do recursion level_v_step(level-1); // do coarse grid correction + this->signals.transfer_up_to(true, level); transfer->prolongate(level, t[level], solution[level-1]); + this->signals.transfer_up_to(false, level); + if (debug>2) deallog << "Prolongate norm " << t[level].l2_norm() << std::endl; solution[level] += t[level]; @@ -184,7 +195,9 @@ Multigrid::level_v_step (const unsigned int level) // post-smoothing if (debug>1) deallog << "Smoothing on level " << level << std::endl; + this->signals.post_smoother_step(true, level); post_smooth->smooth(level, solution[level], defect[level]); + this->signals.post_smoother_step(false, level); if (debug>2) deallog << "Solution norm " << solution[level].l2_norm() @@ -367,6 +380,56 @@ Multigrid::vcycle() } + +template +boost::signals2::connection +Multigrid:: +connect_coarse_solve(const std::function &slot) +{ + return this->signals.coarse_solve.connect(slot); +} + + + +template +boost::signals2::connection +Multigrid:: +connect_transfer_down_to(const std::function &slot) +{ + return this->signals.transfer_down_to.connect(slot); +} + + + +template +boost::signals2::connection +Multigrid:: +connect_transfer_up_to(const std::function &slot) +{ + return this->signals.transfer_up_to.connect(slot); +} + + + +template +boost::signals2::connection +Multigrid:: +connect_pre_smoother_step(const std::function &slot) +{ + return this->signals.pre_smoother_step.connect(slot); +} + + + +template +boost::signals2::connection +Multigrid:: +connect_post_smoother_step(const std::function &slot) +{ + return this->signals.post_smoother_step.connect(slot); +} + + DEAL_II_NAMESPACE_CLOSE #endif diff --git a/tests/matrix_free/parallel_multigrid_adaptive_01.cc b/tests/matrix_free/parallel_multigrid_adaptive_01.cc index 06c17e8639..a2dc912116 100644 --- a/tests/matrix_free/parallel_multigrid_adaptive_01.cc +++ b/tests/matrix_free/parallel_multigrid_adaptive_01.cc @@ -556,8 +556,8 @@ void test () for (typename Triangulation::active_cell_iterator cell=tria.begin_active(); cell != tria.end(); ++cell) if (cell->is_locally_owned() && - (cell->center().norm() < 0.5 && (cell->level() < 5 || - cell->center().norm() > 0.45) + ((cell->center().norm() < 0.5 && (cell->level() < 5 || + cell->center().norm() > 0.45)) || (dim == 2 && cell->center().norm() > 1.2))) cell->set_refine_flag(); diff --git a/tests/matrix_free/parallel_multigrid_adaptive_02.cc b/tests/matrix_free/parallel_multigrid_adaptive_02.cc index 006d586259..9787b79451 100644 --- a/tests/matrix_free/parallel_multigrid_adaptive_02.cc +++ b/tests/matrix_free/parallel_multigrid_adaptive_02.cc @@ -283,8 +283,8 @@ void test () for (typename Triangulation::active_cell_iterator cell=tria.begin_active(); cell != tria.end(); ++cell) if (cell->is_locally_owned() && - (cell->center().norm() < 0.5 && (cell->level() < 5 || - cell->center().norm() > 0.45) + ((cell->center().norm() < 0.5 && (cell->level() < 5 || + cell->center().norm() > 0.45)) || (dim == 2 && cell->center().norm() > 1.2))) cell->set_refine_flag(); diff --git a/tests/matrix_free/parallel_multigrid_adaptive_04.cc b/tests/matrix_free/parallel_multigrid_adaptive_04.cc index 43c2b9b0bf..ec176b4251 100644 --- a/tests/matrix_free/parallel_multigrid_adaptive_04.cc +++ b/tests/matrix_free/parallel_multigrid_adaptive_04.cc @@ -283,8 +283,8 @@ void test () for (typename Triangulation::active_cell_iterator cell=tria.begin_active(); cell != tria.end(); ++cell) if (cell->is_locally_owned() && - (cell->center().norm() < 0.5 && (cell->level() < 5 || - cell->center().norm() > 0.45) + ((cell->center().norm() < 0.5 && (cell->level() < 5 || + cell->center().norm() > 0.45)) || (dim == 2 && cell->center().norm() > 1.2))) cell->set_refine_flag(); diff --git a/tests/matrix_free/parallel_multigrid_adaptive_06.cc b/tests/matrix_free/parallel_multigrid_adaptive_06.cc index 30ad01e7d0..a1b7f56ea1 100644 --- a/tests/matrix_free/parallel_multigrid_adaptive_06.cc +++ b/tests/matrix_free/parallel_multigrid_adaptive_06.cc @@ -346,8 +346,8 @@ void test (const unsigned int nbands = 1) for (typename Triangulation::active_cell_iterator cell=tria.begin_active(); cell != tria.end(); ++cell) if (cell->is_locally_owned() && - (cell->center().norm() < 0.5 && (cell->level() < 5 || - cell->center().norm() > 0.45) + ((cell->center().norm() < 0.5 && (cell->level() < 5 || + cell->center().norm() > 0.45)) || (dim == 2 && cell->center().norm() > 1.2))) cell->set_refine_flag(); diff --git a/tests/matrix_free/parallel_multigrid_adaptive_06ref.cc b/tests/matrix_free/parallel_multigrid_adaptive_06ref.cc index 2823d82a01..1374e7d6c1 100644 --- a/tests/matrix_free/parallel_multigrid_adaptive_06ref.cc +++ b/tests/matrix_free/parallel_multigrid_adaptive_06ref.cc @@ -275,8 +275,8 @@ void test () for (typename Triangulation::active_cell_iterator cell=tria.begin_active(); cell != tria.end(); ++cell) if (cell->is_locally_owned() && - (cell->center().norm() < 0.5 && (cell->level() < 5 || - cell->center().norm() > 0.45) + ((cell->center().norm() < 0.5 && (cell->level() < 5 || + cell->center().norm() > 0.45)) || (dim == 2 && cell->center().norm() > 1.2))) cell->set_refine_flag(); diff --git a/tests/multigrid/events_01.cc b/tests/multigrid/events_01.cc new file mode 100644 index 0000000000..3f00a66802 --- /dev/null +++ b/tests/multigrid/events_01.cc @@ -0,0 +1,590 @@ +/* --------------------------------------------------------------------- + * + * Copyright (C) 2016 - 2017 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 at + * the top level of the deal.II distribution. + * + * --------------------------------------------------------------------- + + * + * Test different coarse grid solvers in parallel (based on step-50) + */ + +#include "../tests.h" + +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +#include +#include + +#include +#include +#include + +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include + + +#include + +namespace LA +{ + using namespace dealii::LinearAlgebraTrilinos; +} + +#include +#include + +namespace Step50 +{ + using namespace dealii; + + template + class LaplaceProblem + { + public: + LaplaceProblem (const unsigned int deg); + void run (); + + private: + typedef LA::MPI::SparseMatrix matrix_t; + typedef LA::MPI::Vector vector_t; + + void setup_system (); + void assemble_system (); + void assemble_multigrid (); + void vcycle (const MGCoarseGridBase &coarse_grid_solver); + + void solve (); + + void refine_grid (); + void output_results (const unsigned int cycle) const; + + parallel::distributed::Triangulation triangulation; + FE_Q fe; + DoFHandler mg_dof_handler; + + + matrix_t system_matrix; + + IndexSet locally_relevant_set; + + ConstraintMatrix constraints; + + vector_t solution; + vector_t system_rhs; + + const unsigned int degree; + + MGLevelObject mg_matrices; + MGLevelObject mg_interface_matrices; + + MGConstrainedDoFs mg_constrained_dofs; + }; + + + template + LaplaceProblem::LaplaceProblem (const unsigned int degree) + : + triangulation (MPI_COMM_WORLD,Triangulation:: + limit_level_difference_at_vertices, + parallel::distributed::Triangulation::construct_multigrid_hierarchy), + fe (degree), + mg_dof_handler (triangulation), + degree(degree) + {} + + + template + void LaplaceProblem::setup_system () + { + mg_dof_handler.distribute_dofs (fe); + mg_dof_handler.distribute_mg_dofs (fe); + + DoFTools::extract_locally_relevant_dofs (mg_dof_handler, + locally_relevant_set); + + solution.reinit(mg_dof_handler.locally_owned_dofs(), MPI_COMM_WORLD); + system_rhs.reinit(mg_dof_handler.locally_owned_dofs(), MPI_COMM_WORLD); + constraints.reinit (locally_relevant_set); + DoFTools::make_hanging_node_constraints (mg_dof_handler, constraints); + + typename FunctionMap::type dirichlet_boundary; + Functions::ConstantFunction homogeneous_dirichlet_bc (0.0); + dirichlet_boundary[0] = &homogeneous_dirichlet_bc; + VectorTools::interpolate_boundary_values (mg_dof_handler, + dirichlet_boundary, + constraints); + constraints.close (); + + DynamicSparsityPattern dsp(mg_dof_handler.n_dofs(), mg_dof_handler.n_dofs()); + DoFTools::make_sparsity_pattern (mg_dof_handler, dsp, constraints); + system_matrix.reinit (mg_dof_handler.locally_owned_dofs(), dsp, MPI_COMM_WORLD, true); + + + mg_constrained_dofs.clear(); + mg_constrained_dofs.initialize(mg_dof_handler, dirichlet_boundary); + + const unsigned int n_levels = triangulation.n_global_levels(); + + mg_interface_matrices.resize(0, n_levels-1); + mg_interface_matrices.clear_elements (); + mg_matrices.resize(0, n_levels-1); + mg_matrices.clear_elements (); + + for (unsigned int level=0; level + void LaplaceProblem::assemble_system () + { + const QGauss quadrature_formula(degree+1); + + FEValues fe_values (fe, quadrature_formula, + update_values | update_gradients | + update_quadrature_points | update_JxW_values); + + const unsigned int dofs_per_cell = fe.dofs_per_cell; + const unsigned int n_q_points = quadrature_formula.size(); + + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); + Vector cell_rhs (dofs_per_cell); + + std::vector local_dof_indices (dofs_per_cell); + + typename DoFHandler::active_cell_iterator + cell = mg_dof_handler.begin_active(), + endc = mg_dof_handler.end(); + for (; cell!=endc; ++cell) + if (cell->is_locally_owned()) + { + cell_matrix = 0; + cell_rhs = 0; + + fe_values.reinit (cell); + + for (unsigned int q_point=0; q_pointget_dof_indices (local_dof_indices); + constraints.distribute_local_to_global (cell_matrix, cell_rhs, + local_dof_indices, + system_matrix, system_rhs); + } + + system_matrix.compress(VectorOperation::add); + system_rhs.compress(VectorOperation::add); + } + + + template + void LaplaceProblem::assemble_multigrid () + { + QGauss quadrature_formula(1+degree); + + FEValues fe_values (fe, quadrature_formula, + update_values | update_gradients | + update_quadrature_points | update_JxW_values); + + const unsigned int dofs_per_cell = fe.dofs_per_cell; + const unsigned int n_q_points = quadrature_formula.size(); + + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); + + std::vector local_dof_indices (dofs_per_cell); + + std::vector boundary_constraints (triangulation.n_global_levels()); + ConstraintMatrix empty_constraints; + for (unsigned int level=0; level::cell_iterator cell = mg_dof_handler.begin(), + endc = mg_dof_handler.end(); + + for (; cell!=endc; ++cell) + if (cell->level_subdomain_id()==triangulation.locally_owned_subdomain()) + { + cell_matrix = 0; + fe_values.reinit (cell); + + for (unsigned int q_point=0; q_pointget_mg_dof_indices (local_dof_indices); + + boundary_constraints[cell->level()] + .distribute_local_to_global (cell_matrix, + local_dof_indices, + mg_matrices[cell->level()]); + const IndexSet &interface_dofs_on_level + = mg_constrained_dofs.get_refinement_edge_indices(cell->level()); + const unsigned int lvl = cell->level(); + + + for (unsigned int i=0; ilevel()]); + } + + for (unsigned int i=0; i + class MGCoarseAMG : public MGCoarseGridBase + { + public: + MGCoarseAMG(const LA::MPI::SparseMatrix &coarse_matrix, + const LA::MPI::PreconditionAMG::AdditionalData additional_data) + : + count(0) + { + precondition_amg.initialize(coarse_matrix, additional_data); + } + + virtual void operator() (const unsigned int, + VECTOR &dst, + const VECTOR &src) const + { + ++count; + precondition_amg.vmult(dst, src); + } + + mutable int count; + + private: + LA::MPI::PreconditionAMG precondition_amg; + }; + + template + void LaplaceProblem::vcycle (const MGCoarseGridBase &coarse_grid_solver) + { + MGTransferPrebuilt mg_transfer(mg_constrained_dofs); + mg_transfer.build_matrices(mg_dof_handler); + + typedef LA::MPI::PreconditionJacobi Smoother; + MGSmootherPrecondition mg_smoother; + mg_smoother.initialize(mg_matrices, Smoother::AdditionalData(0.5)); + mg_smoother.set_steps(2); + mg::Matrix mg_matrix(mg_matrices); + mg::Matrix mg_interface_up(mg_interface_matrices); + mg::Matrix mg_interface_down(mg_interface_matrices); + + Multigrid mg(mg_dof_handler, + mg_matrix, + coarse_grid_solver, + mg_transfer, + mg_smoother, + mg_smoother); + + auto print_transfer_to_mg = [](const bool start) + { + deallog << "transfer_to_mg " << (start?"started":"finished") << std::endl; + }; + + auto print_transfer_to_global = [](const bool start) + { + deallog << "transfer_to_global " << (start?"started":"finished") << std::endl; + }; + + auto print_coarse_solve = [](const bool start, const unsigned int level) + { + deallog << "coarse_solve " << (start?"started":"finished") + << " on level " << level << std::endl; + + }; + + auto print_transfer_down_to = [](const bool start, const unsigned int level) + { + deallog << "transfer_down_to " << (start?"started":"finished") + << " on level " << level << std::endl; + + }; + + auto print_transfer_up_to = [](const bool start, const unsigned int level) + { + deallog << "transfer_up_to " << (start?"started":"finished") + << " on level " << level << std::endl; + }; + + auto print_pre_smoother_step = [](const bool start, const unsigned int level) + { + deallog << "pre_smoother_step " << (start?"started":"finished") + << " on level " << level << std::endl; + }; + + auto print_post_smoother_step = [](const bool start, const unsigned int level) + { + deallog << "post_smoother_step " << (start?"started":"finished") + << " on level " << level << std::endl; + }; + + mg.set_edge_matrices(mg_interface_down, mg_interface_up); + + PreconditionMG > + preconditioner(mg_dof_handler, mg, mg_transfer); + + preconditioner.connect_transfer_to_mg(print_transfer_to_mg); + preconditioner.connect_transfer_to_global(print_transfer_to_global); + + mg.connect_coarse_solve(print_coarse_solve); + mg.connect_transfer_down_to(print_transfer_down_to); + mg.connect_transfer_up_to(print_transfer_up_to); + mg.connect_pre_smoother_step(print_pre_smoother_step); + mg.connect_post_smoother_step(print_post_smoother_step); + + LA::MPI::Vector check, tmp; + check.reinit(mg_dof_handler.locally_owned_dofs(),MPI_COMM_WORLD); + tmp=check; + preconditioner.vmult(check, system_rhs); + + double residual = system_matrix.residual(tmp, check, system_rhs); + deallog << "check residual: " << residual << std::endl; + + } + + + + + template + void LaplaceProblem::solve () + { + matrix_t &coarse_matrix = mg_matrices[0]; + SolverControl coarse_solver_control (1000, 1e-8, false, false); + SolverCG coarse_solver(coarse_solver_control); + + PreconditionIdentity id; + + MGCoarseGridIterativeSolver, + matrix_t, + PreconditionIdentity> coarse_grid_solver(coarse_solver, + coarse_matrix, + id); + + vcycle(coarse_grid_solver); + deallog << "coarse iterations: " << coarse_solver_control.last_step() << std::endl; + + } + + template + void LaplaceProblem::refine_grid () + { + Vector estimated_error_per_cell (triangulation.n_active_cells()); + + LA::MPI::Vector temp_solution; + temp_solution.reinit(locally_relevant_set, MPI_COMM_WORLD); + temp_solution = solution; + + KellyErrorEstimator::estimate (static_cast&>(mg_dof_handler), + QGauss(degree+1), + typename FunctionMap::type(), + temp_solution, + estimated_error_per_cell); + + parallel::distributed::GridRefinement:: + refine_and_coarsen_fixed_fraction (triangulation, + estimated_error_per_cell, + 0.3, 0.0); + + triangulation.execute_coarsening_and_refinement (); + } + + + + template + void LaplaceProblem::output_results (const unsigned int cycle) const + { + } + + + template + void LaplaceProblem::run () + { + for (unsigned int cycle=0; cycle<3; ++cycle) + { + deallog << "Cycle " << cycle << ':' << std::endl; + + { + unsigned int n_subdiv = 6+5*cycle; + triangulation.clear(); + GridGenerator::subdivided_hyper_cube (triangulation, n_subdiv, 0, 1); + triangulation.refine_global(1); + if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0) + triangulation.begin_active()->set_refine_flag(); + triangulation.execute_coarsening_and_refinement(); + } + + deallog << " Number of active cells: " + << triangulation.n_global_active_cells() + << std::endl; + + setup_system (); + + deallog << " Number of degrees of freedom: " + << mg_dof_handler.n_dofs() + << " (by level: "; + for (unsigned int level=0; level laplace_problem(3/*degree*/); + laplace_problem.run (); + } + catch (std::exception &exc) + { + std::cerr << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + throw; + } + catch (...) + { + std::cerr << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + throw; + } + + return 0; +} diff --git a/tests/multigrid/events_01.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=1.output b/tests/multigrid/events_01.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=1.output new file mode 100644 index 0000000000..bf4004d436 --- /dev/null +++ b/tests/multigrid/events_01.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=1.output @@ -0,0 +1,82 @@ + +DEAL::Cycle 0: +DEAL:: Number of active cells: 147 +DEAL:: Number of degrees of freedom: 1406 (by level: 361, 1369, 49) +DEAL::transfer_to_mg started +DEAL::transfer_to_mg finished +DEAL::pre_smoother_step started on level 2 +DEAL::pre_smoother_step finished on level 2 +DEAL::transfer_down_to started on level 2 +DEAL::transfer_down_to finished on level 2 +DEAL::pre_smoother_step started on level 1 +DEAL::pre_smoother_step finished on level 1 +DEAL::transfer_down_to started on level 1 +DEAL::transfer_down_to finished on level 1 +DEAL::coarse_solve started on level 0 +DEAL::coarse_solve finished on level 0 +DEAL::transfer_up_to started on level 1 +DEAL::transfer_up_to finished on level 1 +DEAL::post_smoother_step started on level 1 +DEAL::post_smoother_step finished on level 1 +DEAL::transfer_up_to started on level 2 +DEAL::transfer_up_to finished on level 2 +DEAL::post_smoother_step started on level 2 +DEAL::post_smoother_step finished on level 2 +DEAL::transfer_to_global started +DEAL::transfer_to_global finished +DEAL::check residual: 0.00154307 +DEAL::coarse iterations: 31 +DEAL::Cycle 1: +DEAL:: Number of active cells: 487 +DEAL:: Number of degrees of freedom: 4526 (by level: 1156, 4489, 49) +DEAL::transfer_to_mg started +DEAL::transfer_to_mg finished +DEAL::pre_smoother_step started on level 2 +DEAL::pre_smoother_step finished on level 2 +DEAL::transfer_down_to started on level 2 +DEAL::transfer_down_to finished on level 2 +DEAL::pre_smoother_step started on level 1 +DEAL::pre_smoother_step finished on level 1 +DEAL::transfer_down_to started on level 1 +DEAL::transfer_down_to finished on level 1 +DEAL::coarse_solve started on level 0 +DEAL::coarse_solve finished on level 0 +DEAL::transfer_up_to started on level 1 +DEAL::transfer_up_to finished on level 1 +DEAL::post_smoother_step started on level 1 +DEAL::post_smoother_step finished on level 1 +DEAL::transfer_up_to started on level 2 +DEAL::transfer_up_to finished on level 2 +DEAL::post_smoother_step started on level 2 +DEAL::post_smoother_step finished on level 2 +DEAL::transfer_to_global started +DEAL::transfer_to_global finished +DEAL::check residual: 0.000788123 +DEAL::coarse iterations: 56 +DEAL::Cycle 2: +DEAL:: Number of active cells: 1027 +DEAL:: Number of degrees of freedom: 9446 (by level: 2401, 9409, 49) +DEAL::transfer_to_mg started +DEAL::transfer_to_mg finished +DEAL::pre_smoother_step started on level 2 +DEAL::pre_smoother_step finished on level 2 +DEAL::transfer_down_to started on level 2 +DEAL::transfer_down_to finished on level 2 +DEAL::pre_smoother_step started on level 1 +DEAL::pre_smoother_step finished on level 1 +DEAL::transfer_down_to started on level 1 +DEAL::transfer_down_to finished on level 1 +DEAL::coarse_solve started on level 0 +DEAL::coarse_solve finished on level 0 +DEAL::transfer_up_to started on level 1 +DEAL::transfer_up_to finished on level 1 +DEAL::post_smoother_step started on level 1 +DEAL::post_smoother_step finished on level 1 +DEAL::transfer_up_to started on level 2 +DEAL::transfer_up_to finished on level 2 +DEAL::post_smoother_step started on level 2 +DEAL::post_smoother_step finished on level 2 +DEAL::transfer_to_global started +DEAL::transfer_to_global finished +DEAL::check residual: 0.000522658 +DEAL::coarse iterations: 81 diff --git a/tests/multigrid/transfer_matrix_free_09.cc b/tests/multigrid/transfer_matrix_free_09.cc index 34ca08e0e7..2c2cc16cd7 100644 --- a/tests/multigrid/transfer_matrix_free_09.cc +++ b/tests/multigrid/transfer_matrix_free_09.cc @@ -45,9 +45,9 @@ void check(const FiniteElement &fe) for (typename Triangulation::active_cell_iterator cell=tr.begin_active(); cell != tr.end(); ++cell) if ((cell->center().norm() < 0.5 && (cell->level() < 5 || - cell->center().norm() > 0.45) - || - (dim == 2 && cell->center().norm() > 1.2))) + cell->center().norm() > 0.45)) + || + (dim == 2 && cell->center().norm() > 1.2)) cell->set_refine_flag(); tr.execute_coarsening_and_refinement (); diff --git a/tests/multigrid/transfer_matrix_free_10.cc b/tests/multigrid/transfer_matrix_free_10.cc index d7daa71503..cf1b72fff8 100644 --- a/tests/multigrid/transfer_matrix_free_10.cc +++ b/tests/multigrid/transfer_matrix_free_10.cc @@ -45,9 +45,9 @@ void check(const FiniteElement &fe) for (typename Triangulation::active_cell_iterator cell=tr.begin_active(); cell != tr.end(); ++cell) if ((cell->center().norm() < 0.5 && (cell->level() < 5 || - cell->center().norm() > 0.45) - || - (dim == 2 && cell->center().norm() > 1.2))) + cell->center().norm() > 0.45)) + || + (dim == 2 && cell->center().norm() > 1.2)) cell->set_refine_flag(); tr.execute_coarsening_and_refinement (); -- 2.39.5