From: Martin Kronbichler Date: Wed, 6 Jan 2016 11:00:09 +0000 (+0100) Subject: Implement optimized copy for generic vectors X-Git-Tag: v8.4.0-rc2~109^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=6587d41e88245142924a8e24ebe926af8909d8ed;p=dealii.git Implement optimized copy for generic vectors --- diff --git a/include/deal.II/multigrid/mg_transfer.h b/include/deal.II/multigrid/mg_transfer.h index 2dc701135f..b5316b6e68 100644 --- a/include/deal.II/multigrid/mg_transfer.h +++ b/include/deal.II/multigrid/mg_transfer.h @@ -233,6 +233,14 @@ protected: std::vector > > copy_indices_level_mine; + /** + * Stores whether the copy operation from the global to the level vector is + * actually a plain copy to the finest level. This means that the grid has + * no adaptive refinement and the numbering on the finest multigrid level is + * the same as in the global case. + */ + bool perform_plain_copy; + /** * The vector that stores what has been given to the * set_component_to_block_map() function. diff --git a/include/deal.II/multigrid/mg_transfer.templates.h b/include/deal.II/multigrid/mg_transfer.templates.h index a3bd668bdf..59bb8f2413 100644 --- a/include/deal.II/multigrid/mg_transfer.templates.h +++ b/include/deal.II/multigrid/mg_transfer.templates.h @@ -161,6 +161,34 @@ namespace /* ------------------ MGLevelGlobalTransfer ----------------- */ +namespace internal +{ + // generic copy function of two different vectors -> need to access each + // individual entry + template + void copy_vector (const std::vector > ©_indices, + const T &src, + V &dst) + { + // we should have i->second == i->first, therefore we can use the same + // function for both copying to mg as well as copying from mg + for (std::vector >:: + const_iterator i = copy_indices.begin(); i != copy_indices.end(); ++i) + dst(i->first) = src(i->first); + dst.compress(VectorOperation::insert); + } + + // specialized copy function for the same vector + template + void copy_vector (const std::vector > &, + const T &src, + T &dst) + { + dst = src; + } +} + + template template void @@ -175,6 +203,24 @@ MGLevelGlobalTransfer::copy_to_mg std::cout << "copy_to_mg src " << src.l2_norm() << std::endl; MPI_Barrier(MPI_COMM_WORLD); #endif + + if (perform_plain_copy) + { + // if the finest multigrid level covers the whole domain (i.e., no + // adaptive refinement) and the numbering of the finest level DoFs and + // the global DoFs are the same, we can do a plain copy + AssertDimension(dst[dst.max_level()].size(), src.size()); + internal::copy_vector(copy_indices[dst.max_level()], src, dst[dst.max_level()]); + + // do the initial restriction + for (unsigned int level=mg_dof_handler.get_triangulation().n_global_levels()-1; level != 0; ) + { + --level; + this->restrict_and_add (level+1, dst[level], dst[level+1]); + } + return; + } + for (unsigned int level=mg_dof_handler.get_triangulation().n_global_levels(); level != 0;) { --level; @@ -225,6 +271,13 @@ MGLevelGlobalTransfer::copy_from_mg OutVector &dst, const MGLevelObject &src) const { + if (perform_plain_copy) + { + AssertDimension(dst.size(), src[src.max_level()].size()); + internal::copy_vector(copy_indices[src.max_level()], src[src.max_level()], dst); + return; + } + // For non-DG: degrees of freedom in the refinement face may need special // attention, since they belong to the coarse level, but have fine level // basis functions @@ -331,6 +384,8 @@ MGLevelGlobalTransfer >::copy_to_mg VectorView dst_view (src.local_size(), dst[dst.max_level()].begin()); VectorView src_view (src.local_size(), src.begin()); static_cast &>(dst_view) = static_cast &>(src_view); + + // do the initial restriction for (unsigned int level=mg_dof_handler.get_triangulation().n_global_levels()-1; level != 0; ) { --level; @@ -387,10 +442,6 @@ MGLevelGlobalTransfer >::copy_from_mg parallel::distributed::Vector &dst, const MGLevelObject > &src) const { - // For non-DG: degrees of freedom in the refinement face may need special - // attention, since they belong to the coarse level, but have fine level - // basis functions - if (perform_plain_copy) { // In this case, we can simply copy the local range (in parallel by @@ -404,6 +455,10 @@ MGLevelGlobalTransfer >::copy_from_mg return; } + // For non-DG: degrees of freedom in the refinement face may need special + // attention, since they belong to the coarse level, but have fine level + // basis functions + dst = 0; for (unsigned int level=0; level::fill_and_communicate_copy_indices { fill_copy_indices(mg_dof, mg_constrained_dofs, copy_indices, copy_indices_global_mine, copy_indices_level_mine); + + // check if we can run a plain copy operation between the global DoFs and + // the finest level. + perform_plain_copy = + (copy_indices.back().size() == mg_dof.locally_owned_dofs().n_elements()) + && + (mg_dof.locally_owned_dofs().n_elements() == + mg_dof.locally_owned_mg_dofs(mg_dof.get_triangulation().n_global_levels()-1).n_elements()); + if (perform_plain_copy) + { + AssertDimension(copy_indices_global_mine.back().size(), 0); + AssertDimension(copy_indices_level_mine.back().size(), 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 *ptria = + dynamic_cast *> + (&mg_dof.get_tria()); + const MPI_Comm mpi_communicator = ptria != 0 ? ptria->get_communicator() : + MPI_COMM_SELF; + perform_plain_copy = + Utilities::MPI::min(static_cast(perform_plain_copy), + mpi_communicator); + }