From: tcclevenger Date: Sat, 3 Aug 2019 01:06:36 +0000 (-0600) Subject: add function X-Git-Tag: v9.2.0-rc1~1306^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F8444%2Fhead;p=dealii.git add function --- diff --git a/include/deal.II/multigrid/mg_tools.h b/include/deal.II/multigrid/mg_tools.h index 99048d1ecf..c38bfaf58b 100644 --- a/include/deal.II/multigrid/mg_tools.h +++ b/include/deal.II/multigrid/mg_tools.h @@ -279,6 +279,19 @@ namespace MGTools template unsigned int max_level_for_coarse_mesh(const Triangulation &tria); + + /** + * Return the imbalance of the parallel distribution of the multigrid + * mesh hierarchy. Ideally this value is equal to 1 (every processor owns + * the same number of cells on each level, approximately true for most + * globally refined meshes). Values greater than 1 estimate the slowdown + * one should see in a geometric multigrid v-cycle as compared with the same + * computation on a perfectly distributed mesh hierarchy. + */ + template + double + workload_imbalance(const Triangulation &tria); + } // namespace MGTools /* @} */ diff --git a/source/multigrid/mg_tools.cc b/source/multigrid/mg_tools.cc index d107f40bda..4d3c8c04f2 100644 --- a/source/multigrid/mg_tools.cc +++ b/source/multigrid/mg_tools.cc @@ -1661,6 +1661,58 @@ namespace MGTools return global_min; } + + + + template + double + workload_imbalance(const Triangulation &tria) + { + double workload_imbalance = 1.0; + + // It is only necessary to calculate the imbalance + // on a distributed mesh. The imbalance is always + // 1.0 for the serial case. + if (const parallel::Triangulation *tr = + dynamic_cast *>(&tria)) + { + const unsigned int n_proc = + Utilities::MPI::n_mpi_processes(tr->get_communicator()); + const unsigned int n_global_levels = tr->n_global_levels(); + + // This value will represent the sum over the multigrid + // levels of the maximum number of cells owned by any + // one processesor on that level. + types::global_dof_index work_estimate = 0; + + // Sum of all cells in the multigrid hierarchy + types::global_dof_index total_cells_in_hierarchy = 0; + + for (int lvl = n_global_levels - 1; lvl >= 0; --lvl) + { + // Number of cells this processor owns on this level + types::global_dof_index n_owned_cells_on_lvl = 0; + + for (const auto &cell : tr->cell_iterators_on_level(lvl)) + if (cell->is_locally_owned_on_level()) + ++n_owned_cells_on_lvl; + + work_estimate += + dealii::Utilities::MPI::max(n_owned_cells_on_lvl, + tr->get_communicator()); + + total_cells_in_hierarchy += + dealii::Utilities::MPI::sum(n_owned_cells_on_lvl, + tr->get_communicator()); + } + + const double ideal_work = + total_cells_in_hierarchy / static_cast(n_proc); + workload_imbalance = work_estimate / ideal_work; + } + + return workload_imbalance; + } } // namespace MGTools diff --git a/source/multigrid/mg_tools.inst.in b/source/multigrid/mg_tools.inst.in index 1a78b8e2ef..e1252193fb 100644 --- a/source/multigrid/mg_tools.inst.in +++ b/source/multigrid/mg_tools.inst.in @@ -162,6 +162,10 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) template unsigned int max_level_for_coarse_mesh( const Triangulation &tria); + + template double + workload_imbalance( + const Triangulation &tria); #endif \} } diff --git a/tests/multigrid/workload_imbalance.cc b/tests/multigrid/workload_imbalance.cc new file mode 100644 index 0000000000..20492c9d99 --- /dev/null +++ b/tests/multigrid/workload_imbalance.cc @@ -0,0 +1,81 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2006 - 2018 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. +// +// --------------------------------------------------------------------- + + +// check function MGTools::workload_imbalance() + +#include + +#include +#include + +#include +#include +#include +#include + +#include +#include + +#include + +#include "../tests.h" + + +template +void +test() +{ + parallel::distributed::Triangulation tria( + MPI_COMM_WORLD, + Triangulation::limit_level_difference_at_vertices, + parallel::distributed::Triangulation::construct_multigrid_hierarchy); + GridGenerator::hyper_cube(tria, 0, 1); + tria.refine_global(1); + + unsigned int counter = 0; + for (auto cell : tria.active_cell_iterators()) + { + cell->set_refine_flag(); + ++counter; + + if (dim == 2 && counter == 1) + break; + if (dim == 3 && counter == 4) + break; + } + tria.execute_coarsening_and_refinement(); + + const double workload_imbalance = MGTools::workload_imbalance(tria); + deallog << "Workload imbalance: " << workload_imbalance << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + deallog.push("2d"); + test<2>(); + deallog.pop(); + deallog.push("3d"); + test<3>(); + deallog.pop(); + + deallog << "OK" << std::endl; +} diff --git a/tests/multigrid/workload_imbalance.with_mpi=true.with_p4est=true.mpirun=3.output b/tests/multigrid/workload_imbalance.with_mpi=true.with_p4est=true.mpirun=3.output new file mode 100644 index 0000000000..9d5b5fc113 --- /dev/null +++ b/tests/multigrid/workload_imbalance.with_mpi=true.with_p4est=true.mpirun=3.output @@ -0,0 +1,14 @@ + +DEAL:0:2d::Workload imbalance: 2.66667 +DEAL:0:3d::Workload imbalance: 1.60976 +DEAL:0::OK + +DEAL:1:2d::Workload imbalance: 2.66667 +DEAL:1:3d::Workload imbalance: 1.60976 +DEAL:1::OK + + +DEAL:2:2d::Workload imbalance: 2.66667 +DEAL:2:3d::Workload imbalance: 1.60976 +DEAL:2::OK +