From 412a627db00783cf1a33297c287aa7117c10ced2 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Fri, 20 Sep 2019 09:46:57 +0200 Subject: [PATCH] Change to Table for copy in MGLevelGlobalTransfer --- include/deal.II/base/table.h | 1 - include/deal.II/multigrid/mg_transfer.h | 27 ++--- .../deal.II/multigrid/mg_transfer.templates.h | 86 +++++++------- .../deal.II/multigrid/mg_transfer_internal.h | 3 +- source/multigrid/mg_level_global_transfer.cc | 108 ++++++++---------- source/multigrid/mg_transfer_internal.cc | 15 ++- source/multigrid/mg_transfer_internal.inst.in | 2 +- 7 files changed, 114 insertions(+), 128 deletions(-) diff --git a/include/deal.II/base/table.h b/include/deal.II/base/table.h index af1d41a065..7e2b3d498c 100644 --- a/include/deal.II/base/table.h +++ b/include/deal.II/base/table.h @@ -1254,7 +1254,6 @@ public: typename AlignedVector::const_reference operator()(const size_type i, const size_type j) const; - /** * Direct access to one element of the table by specifying all indices at * the same time. Range checks are performed. diff --git a/include/deal.II/multigrid/mg_transfer.h b/include/deal.II/multigrid/mg_transfer.h index ba66e59b4a..d5537f61df 100644 --- a/include/deal.II/multigrid/mg_transfer.h +++ b/include/deal.II/multigrid/mg_transfer.h @@ -541,49 +541,44 @@ protected: * Mapping for the copy_to_mg() and copy_from_mg() functions. Here only * index pairs locally owned is stored. * - * The data is organized as follows: one vector per level. Each element of - * these vectors contains first the global index, then the level index. + * The data is organized as follows: one table per level. This table has two + * rows. The first row contains the global index, the second one the level + * index. */ - std::vector>> copy_indices; - + std::vector> copy_indices; /** * Same as above, but used to transfer solution vectors. */ - std::vector>> - solution_copy_indices; + std::vector> solution_copy_indices; /** * Additional degrees of freedom for the copy_to_mg() function. These are * the ones where the global degree of freedom is locally owned and the * level degree of freedom is not. * - * Organization of the data is like for @p copy_indices_mine. + * Organization of the data is like for @p copy_indices. */ - std::vector>> - copy_indices_global_mine; + std::vector> copy_indices_global_mine; /** * Same as above, but used to transfer solution vectors. */ - std::vector>> - solution_copy_indices_global_mine; + std::vector> solution_copy_indices_global_mine; /** * Additional degrees of freedom for the copy_from_mg() function. These are * the ones where the level degree of freedom is locally owned and the * global degree of freedom is not. * - * Organization of the data is like for @p copy_indices_mine. + * Organization of the data is like for @p copy_indices. */ - std::vector>> - copy_indices_level_mine; + std::vector> copy_indices_level_mine; /** * Same as above, but used to transfer solution vectors. */ - std::vector>> - solution_copy_indices_level_mine; + std::vector> solution_copy_indices_level_mine; /** * This variable stores whether the copy operation from the global to the diff --git a/include/deal.II/multigrid/mg_transfer.templates.h b/include/deal.II/multigrid/mg_transfer.templates.h index 62ca7118df..533d35f68c 100644 --- a/include/deal.II/multigrid/mg_transfer.templates.h +++ b/include/deal.II/multigrid/mg_transfer.templates.h @@ -417,13 +417,11 @@ MGLevelGlobalTransfer>::copy_to_mg( { LinearAlgebra::distributed::Vector &this_ghosted_global_vector = solution_transfer ? solution_ghosted_global_vector : ghosted_global_vector; - const std::vector>> - &this_copy_indices = - solution_transfer ? solution_copy_indices : copy_indices; - const std::vector>> - &this_copy_indices_level_mine = - solution_transfer ? solution_copy_indices_level_mine : - copy_indices_level_mine; + const std::vector> &this_copy_indices = + solution_transfer ? solution_copy_indices : copy_indices; + const std::vector> &this_copy_indices_level_mine = + solution_transfer ? solution_copy_indices_level_mine : + copy_indices_level_mine; (void)mg_dof_handler; @@ -466,11 +464,17 @@ MGLevelGlobalTransfer>::copy_to_mg( } else if (perform_renumbered_plain_copy) { - Assert(copy_indices_level_mine.back().empty(), ExcInternalError()); + AssertDimension(dst[dst.max_level()].local_size(), src.local_size()); + AssertDimension(this_copy_indices.back().n_cols(), src.local_size()); + Assert(copy_indices_level_mine.back().n_rows() == 0, ExcInternalError()); LinearAlgebra::distributed::Vector &dst_level = dst[dst.max_level()]; - for (std::pair i : this_copy_indices.back()) - dst_level.local_element(i.second) = src.local_element(i.first); + // as opposed to the copy_unknowns lambda below, we here know that all + // src elements will be touched, so we only need to do indirect + // addressing on one index + for (unsigned int i = 0; i < this_copy_indices.back().n_cols(); ++i) + dst_level.local_element(this_copy_indices.back()(1, i)) = + src.local_element(i); return; } @@ -487,20 +491,20 @@ MGLevelGlobalTransfer>::copy_to_mg( { --level; - using dof_pair_iterator = - std::vector>::const_iterator; LinearAlgebra::distributed::Vector &dst_level = dst[level]; + auto copy_unknowns = [&](const Table<2, unsigned int> &indices) { + for (unsigned int i = 0; i < indices.n_cols(); ++i) + dst_level.local_element(indices(1, i)) = + this_ghosted_global_vector.local_element(indices(0, i)); + }; + // first copy local unknowns - for (const auto &indices : this_copy_indices[level]) - dst_level.local_element(indices.second) = - this_ghosted_global_vector.local_element(indices.first); + copy_unknowns(this_copy_indices[level]); // Do the same for the indices where the level index is local, but the // global index is not - for (const auto &indices : this_copy_indices_level_mine[level]) - dst_level.local_element(indices.second) = - this_ghosted_global_vector.local_element(indices.first); + copy_unknowns(this_copy_indices_level_mine[level]); dst_level.compress(VectorOperation::insert); } @@ -532,12 +536,16 @@ MGLevelGlobalTransfer>::copy_from_mg( } else if (perform_renumbered_plain_copy) { + AssertDimension(src[src.max_level()].local_size(), dst.local_size()); + AssertDimension(copy_indices.back().n_cols(), dst.local_size()); + Assert(copy_indices_global_mine.back().n_rows() == 0, ExcInternalError()); Assert(copy_indices_global_mine.back().empty(), ExcInternalError()); const LinearAlgebra::distributed::Vector &src_level = src[src.max_level()]; dst.zero_out_ghosts(); - for (std::pair i : copy_indices.back()) - dst.local_element(i.first) = src_level.local_element(i.second); + for (unsigned int i = 0; i < copy_indices.back().n_cols(); ++i) + dst.local_element(i) = + src_level.local_element(copy_indices.back()(1, i)); return; } @@ -548,9 +556,6 @@ MGLevelGlobalTransfer>::copy_from_mg( dst = 0; for (unsigned int level = src.min_level(); level <= src.max_level(); ++level) { - using dof_pair_iterator = - std::vector>::const_iterator; - // the ghosted vector should already have the correct local size (but // different parallel layout) AssertDimension(ghosted_level_vector[level].local_size(), @@ -563,18 +568,18 @@ MGLevelGlobalTransfer>::copy_from_mg( ghosted_vector = src[level]; ghosted_vector.update_ghost_values(); + auto copy_unknowns = [&](const Table<2, unsigned int> &indices) { + for (unsigned int i = 0; i < indices.n_cols(); ++i) + dst.local_element(indices(0, i)) = + ghosted_vector.local_element(indices(1, i)); + }; + // first copy local unknowns - for (dof_pair_iterator i = copy_indices[level].begin(); - i != copy_indices[level].end(); - ++i) - dst.local_element(i->first) = ghosted_vector.local_element(i->second); + copy_unknowns(copy_indices[level]); // Do the same for the indices where the level index is local, but the // global index is not - for (dof_pair_iterator i = copy_indices_global_mine[level].begin(); - i != copy_indices_global_mine[level].end(); - ++i) - dst.local_element(i->first) = ghosted_vector.local_element(i->second); + copy_unknowns(copy_indices_global_mine[level]); } dst.compress(VectorOperation::insert); } @@ -597,9 +602,6 @@ MGLevelGlobalTransfer>:: dst.zero_out_ghosts(); for (unsigned int level = src.min_level(); level <= src.max_level(); ++level) { - using dof_pair_iterator = - std::vector>::const_iterator; - // the ghosted vector should already have the correct local size (but // different parallel layout) AssertDimension(ghosted_level_vector[level].local_size(), @@ -612,18 +614,18 @@ MGLevelGlobalTransfer>:: ghosted_vector = src[level]; ghosted_vector.update_ghost_values(); + auto copy_unknowns = [&](const Table<2, unsigned int> &indices) { + for (unsigned int i = 0; i < indices.n_cols(); ++i) + dst.local_element(indices(0, i)) += + ghosted_vector.local_element(indices(1, i)); + }; + // first add local unknowns - for (dof_pair_iterator i = copy_indices[level].begin(); - i != copy_indices[level].end(); - ++i) - dst.local_element(i->first) += ghosted_vector.local_element(i->second); + copy_unknowns(copy_indices[level]); // Do the same for the indices where the level index is local, but the // global index is not - for (dof_pair_iterator i = copy_indices_global_mine[level].begin(); - i != copy_indices_global_mine[level].end(); - ++i) - dst.local_element(i->first) += ghosted_vector.local_element(i->second); + copy_unknowns(copy_indices_global_mine[level]); } dst.compress(VectorOperation::add); } diff --git a/include/deal.II/multigrid/mg_transfer_internal.h b/include/deal.II/multigrid/mg_transfer_internal.h index 804348ccf1..e90e7ca3da 100644 --- a/include/deal.II/multigrid/mg_transfer_internal.h +++ b/include/deal.II/multigrid/mg_transfer_internal.h @@ -129,8 +129,7 @@ namespace internal std::vector &n_owned_level_cells, std::vector>> &dirichlet_indices, std::vector> & weights_on_refined, - std::vector>> - ©_indices_global_mine, + std::vector> ©_indices_global_mine, MGLevelObject> &ghosted_level_vector); diff --git a/source/multigrid/mg_level_global_transfer.cc b/source/multigrid/mg_level_global_transfer.cc index 5ac77d336d..633b4ff980 100644 --- a/source/multigrid/mg_level_global_transfer.cc +++ b/source/multigrid/mg_level_global_transfer.cc @@ -164,15 +164,12 @@ namespace SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer>> - mg_constrained_dofs, - const MPI_Comm mpi_communicator, - const bool transfer_solution_vectors, - std::vector>> - ©_indices, - std::vector>> - ©_indices_global_mine, - std::vector>> - & copy_indices_level_mine, + mg_constrained_dofs, + const MPI_Comm mpi_communicator, + const bool transfer_solution_vectors, + std::vector> & copy_indices, + std::vector> & copy_indices_global_mine, + std::vector> & copy_indices_level_mine, LinearAlgebra::distributed::Vector &ghosted_global_vector, MGLevelObject> &ghosted_level_vector) @@ -251,40 +248,32 @@ namespace *ghosted_global_vector.get_partitioner(); const Utilities::MPI::Partitioner &level_partitioner = *ghosted_level_vector[level].get_partitioner(); - // owned-owned case: the locally owned indices are going to control - // the local index - copy_indices[level].resize(my_copy_indices[level].size()); - for (unsigned int i = 0; i < my_copy_indices[level].size(); ++i) - copy_indices[level][i] = std::pair( - global_partitioner.global_to_local(my_copy_indices[level][i].first), - level_partitioner.global_to_local( - my_copy_indices[level][i].second)); - - // remote-owned case: the locally owned indices for the level and the - // ghost dofs for the global indices set the local index - copy_indices_level_mine[level].resize( - my_copy_indices_level_mine[level].size()); - for (unsigned int i = 0; i < my_copy_indices_level_mine[level].size(); - ++i) - copy_indices_level_mine[level][i] = - std::pair( - global_partitioner.global_to_local( - my_copy_indices_level_mine[level][i].first), - level_partitioner.global_to_local( - my_copy_indices_level_mine[level][i].second)); - - // owned-remote case: the locally owned indices for the global dofs - // and the ghost dofs for the level indices set the local index - copy_indices_global_mine[level].resize( - my_copy_indices_global_mine[level].size()); - for (unsigned int i = 0; i < my_copy_indices_global_mine[level].size(); - ++i) - copy_indices_global_mine[level][i] = - std::pair( - global_partitioner.global_to_local( - my_copy_indices_global_mine[level][i].first), - level_partitioner.global_to_local( - my_copy_indices_global_mine[level][i].second)); + + auto translate_indices = + [&](const std::vector< + std::pair> + & global_copy_indices, + Table<2, unsigned int> &local_copy_indices) { + local_copy_indices.reinit(2, global_copy_indices.size()); + for (unsigned int i = 0; i < global_copy_indices.size(); ++i) + { + local_copy_indices(0, i) = global_partitioner.global_to_local( + global_copy_indices[i].first); + local_copy_indices(1, i) = level_partitioner.global_to_local( + global_copy_indices[i].second); + } + }; + + // owned-owned case + translate_indices(my_copy_indices[level], copy_indices[level]); + + // remote-owned case + translate_indices(my_copy_indices_level_mine[level], + copy_indices_level_mine[level]); + + // owned-remote case + translate_indices(my_copy_indices_global_mine[level], + copy_indices_global_mine[level]); } } } // namespace @@ -322,21 +311,20 @@ MGLevelGlobalTransfer>:: solution_ghosted_level_vector); bool my_perform_renumbered_plain_copy = - (this->copy_indices.back().size() == + (this->copy_indices.back().n_cols() == mg_dof.locally_owned_dofs().n_elements()); bool my_perform_plain_copy = false; if (my_perform_renumbered_plain_copy) { my_perform_plain_copy = true; - AssertDimension(this->copy_indices_global_mine.back().size(), 0); - AssertDimension(this->copy_indices_level_mine.back().size(), 0); + AssertDimension(this->copy_indices_global_mine.back().n_rows(), 0); + AssertDimension(this->copy_indices_level_mine.back().n_rows(), 0); // check whether there is a renumbering of degrees of freedom on // either the finest level or the global dofs, which means that we // cannot apply a plain copy - for (unsigned int i = 0; i < this->copy_indices.back().size(); ++i) - if (this->copy_indices.back()[i].first != - this->copy_indices.back()[i].second) + for (unsigned int i = 0; i < this->copy_indices.back().n_cols(); ++i) + if (this->copy_indices.back()(0, i) != this->copy_indices.back()(1, i)) { my_perform_plain_copy = false; break; @@ -355,6 +343,9 @@ MGLevelGlobalTransfer>:: // if we do a plain copy, no need to hold additional ghosted vectors if (perform_renumbered_plain_copy) { + for (unsigned int i = 0; i < this->copy_indices.back().n_cols(); ++i) + AssertDimension(this->copy_indices.back()(0, i), i); + ghosted_global_vector.reinit(0); ghosted_level_vector.resize(0, 0); solution_ghosted_global_vector.reinit(0); @@ -389,24 +380,25 @@ MGLevelGlobalTransfer>:: { for (unsigned int level = 0; level < copy_indices.size(); ++level) { - for (unsigned int i = 0; i < copy_indices[level].size(); ++i) - os << "copy_indices[" << level << "]\t" << copy_indices[level][i].first - << '\t' << copy_indices[level][i].second << std::endl; + for (unsigned int i = 0; i < copy_indices[level].n_cols(); ++i) + os << "copy_indices[" << level << "]\t" << copy_indices[level](0, i) + << '\t' << copy_indices[level](1, i) << std::endl; } for (unsigned int level = 0; level < copy_indices_level_mine.size(); ++level) { - for (unsigned int i = 0; i < copy_indices_level_mine[level].size(); ++i) + for (unsigned int i = 0; i < copy_indices_level_mine[level].n_cols(); ++i) os << "copy_ifrom [" << level << "]\t" - << copy_indices_level_mine[level][i].first << '\t' - << copy_indices_level_mine[level][i].second << std::endl; + << copy_indices_level_mine[level](0, i) << '\t' + << copy_indices_level_mine[level](1, i) << std::endl; } for (unsigned int level = 0; level < copy_indices_global_mine.size(); ++level) { - for (unsigned int i = 0; i < copy_indices_global_mine[level].size(); ++i) + for (unsigned int i = 0; i < copy_indices_global_mine[level].n_cols(); + ++i) os << "copy_ito [" << level << "]\t" - << copy_indices_global_mine[level][i].first << '\t' - << copy_indices_global_mine[level][i].second << std::endl; + << copy_indices_global_mine[level](0, i) << '\t' + << copy_indices_global_mine[level](1, i) << std::endl; } } diff --git a/source/multigrid/mg_transfer_internal.cc b/source/multigrid/mg_transfer_internal.cc index 8366291d47..13d19f3b1a 100644 --- a/source/multigrid/mg_transfer_internal.cc +++ b/source/multigrid/mg_transfer_internal.cc @@ -377,8 +377,7 @@ namespace internal std::vector & ghosted_level_dofs, const MPI_Comm & communicator, LinearAlgebra::distributed::Vector &ghosted_level_vector, - std::vector> - ©_indices_global_mine) + Table<2, unsigned int> & copy_indices_global_mine) { std::sort(ghosted_level_dofs.begin(), ghosted_level_dofs.end()); IndexSet ghosted_dofs(locally_owned.size()); @@ -394,10 +393,11 @@ namespace internal // partitioner that we are going to use for the vector const auto &part = ghosted_level_vector.get_partitioner(); ghosted_dofs.add_indices(part->ghost_indices()); - for (auto &indices : copy_indices_global_mine) - indices.second = locally_owned.n_elements() + - ghosted_dofs.index_within_set( - part->local_to_global(indices.second)); + for (unsigned int i = 0; i < copy_indices_global_mine.n_cols(); ++i) + copy_indices_global_mine(1, i) = + locally_owned.n_elements() + + ghosted_dofs.index_within_set( + part->local_to_global(copy_indices_global_mine(1, i))); } ghosted_level_vector.reinit(locally_owned, ghosted_dofs, communicator); } @@ -592,8 +592,7 @@ namespace internal std::vector &n_owned_level_cells, std::vector>> &dirichlet_indices, std::vector> & weights_on_refined, - std::vector>> - ©_indices_global_mine, + std::vector> ©_indices_global_mine, MGLevelObject> &ghosted_level_vector) { diff --git a/source/multigrid/mg_transfer_internal.inst.in b/source/multigrid/mg_transfer_internal.inst.in index ea67a0e379..204503303a 100644 --- a/source/multigrid/mg_transfer_internal.inst.in +++ b/source/multigrid/mg_transfer_internal.inst.in @@ -85,7 +85,7 @@ for (deal_II_dimension : DIMENSIONS; S : REAL_SCALARS) std::vector &, std::vector>> &, std::vector> &, - std::vector>> &, + std::vector> &, MGLevelObject> &); \} \} -- 2.39.5