From 0ca09c50447cdb7c8a0da7e50cd78480f3bc8516 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Thu, 23 Apr 2020 08:45:25 +0200 Subject: [PATCH] Create GridTools::compute_vertices_with_ghost_neighbors() --- include/deal.II/distributed/p4est_wrappers.h | 12 ---- include/deal.II/distributed/tria_base.h | 7 ++- include/deal.II/grid/grid_tools.h | 13 ++++- source/distributed/p4est_wrappers.cc | 52 ----------------- source/distributed/p4est_wrappers.inst.in | 27 --------- source/distributed/tria.cc | 2 +- source/distributed/tria_base.cc | 47 +-------------- source/dofs/dof_handler_policy.cc | 5 +- source/grid/grid_tools.cc | 60 ++++++++++++++++++++ source/grid/grid_tools.inst.in | 4 ++ 10 files changed, 87 insertions(+), 142 deletions(-) diff --git a/include/deal.II/distributed/p4est_wrappers.h b/include/deal.II/distributed/p4est_wrappers.h index ff796047d8..77a763ba8b 100644 --- a/include/deal.II/distributed/p4est_wrappers.h +++ b/include/deal.II/distributed/p4est_wrappers.h @@ -542,18 +542,6 @@ namespace internal tree_exists_locally(const typename types::forest *parallel_forest, const typename types::topidx coarse_grid_cell); - - - /** - * Compute the ghost neighbors surrounding each vertex by querying p4est - */ - template - std::map> - compute_vertices_with_ghost_neighbors( - const dealii::parallel::distributed::Triangulation &tria, - typename dealii::internal::p4est::types::forest *parallel_forest, - typename dealii::internal::p4est::types::ghost * parallel_ghost); - } // namespace p4est } // namespace internal diff --git a/include/deal.II/distributed/tria_base.h b/include/deal.II/distributed/tria_base.h index 6b6efd4b38..a4276bfbb6 100644 --- a/include/deal.II/distributed/tria_base.h +++ b/include/deal.II/distributed/tria_base.h @@ -213,8 +213,13 @@ namespace parallel /** * Return a map that, for each vertex, lists all the processors whose * subdomains are adjacent to that vertex. + * + * @deprecated Use GridTools::compute_vertices_with_ghost_neighbors() + * instead of + * parallel::TriangulationBase::compute_vertices_with_ghost_neighbors(). */ - virtual std::map> + DEAL_II_DEPRECATED virtual std::map> compute_vertices_with_ghost_neighbors() const; /** diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index b8cb641a6e..03aa9b8fa2 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -2986,6 +2986,17 @@ namespace GridTools std::map> &coinciding_vertex_groups, std::map &vertex_to_coinciding_vertex_group); + /** + * Return a map that, for each vertex, lists all the processes whose + * subdomains are adjacent to that vertex. + * + * @param[in] tria Triangulation. + */ + template + std::map> + compute_vertices_with_ghost_neighbors( + const Triangulation &tria); + /** * A structure that allows the transfer of cell data of type @p T from one processor * to another. It corresponds to a packed buffer that stores a vector of @@ -3927,7 +3938,7 @@ namespace GridTools std::map> vertices_with_ghost_neighbors = - tria->compute_vertices_with_ghost_neighbors(); + GridTools::compute_vertices_with_ghost_neighbors(*tria); for (const auto &cell : tria->active_cell_iterators()) if (cell->is_locally_owned()) diff --git a/source/distributed/p4est_wrappers.cc b/source/distributed/p4est_wrappers.cc index d71f6fc9b6..bbaf79f4cc 100644 --- a/source/distributed/p4est_wrappers.cc +++ b/source/distributed/p4est_wrappers.cc @@ -473,58 +473,6 @@ namespace internal std::size_t (&functions<2>::connectivity_memory_used)( types<2>::connectivity *p4est) = p4est_connectivity_memory_used; - template - std::map> - compute_vertices_with_ghost_neighbors( - const typename dealii::parallel::distributed::Triangulation - & tria, - typename dealii::internal::p4est::types::forest *parallel_forest, - typename dealii::internal::p4est::types::ghost * parallel_ghost) - { - std::map> - vertices_with_ghost_neighbors; - - dealii::internal::p4est::FindGhosts fg; - fg.subids = sc_array_new(sizeof(dealii::types::subdomain_id)); - fg.triangulation = &tria; - fg.vertices_with_ghost_neighbors = &vertices_with_ghost_neighbors; - - switch (dim) - { - case 2: - p4est_iterate( - reinterpret_cast::forest *>( - parallel_forest), - reinterpret_cast::ghost *>( - parallel_ghost), - static_cast(&fg), - nullptr, - find_ghosts_face<2, spacedim>, - find_ghosts_corner<2, spacedim>); - break; - - case 3: - p8est_iterate( - reinterpret_cast::forest *>( - parallel_forest), - reinterpret_cast::ghost *>( - parallel_ghost), - static_cast(&fg), - nullptr, - find_ghosts_face<3, 3>, - find_ghosts_edge<3, 3>, - find_ghosts_corner<3, 3>); - break; - - default: - Assert(false, ExcNotImplemented()); - } - - sc_array_destroy(fg.subids); - - return vertices_with_ghost_neighbors; - } - constexpr unsigned int functions<2>::max_level; void (&functions<2>::transfer_fixed)(const types<2>::gloidx *dest_gfq, diff --git a/source/distributed/p4est_wrappers.inst.in b/source/distributed/p4est_wrappers.inst.in index 81ae074778..47b1445833 100644 --- a/source/distributed/p4est_wrappers.inst.in +++ b/source/distributed/p4est_wrappers.inst.in @@ -54,30 +54,3 @@ for (deal_II_dimension : DIMENSIONS) \} #endif // DEAL_II_WITH_P4EST } - - -for (deal_II_dimension, deal_II_space_dimension : DIMENSIONS) - { -#ifdef DEAL_II_WITH_P4EST - - namespace internal - \{ - namespace p4est - \{ -# if deal_II_dimension > 1 -# if deal_II_space_dimension >= deal_II_dimension - template std::map> - compute_vertices_with_ghost_neighbors( - const dealii::parallel::distributed:: - Triangulation &tria, - typename dealii::internal::p4est::types::forest - *parallel_forest, - typename dealii::internal::p4est::types::ghost - *parallel_ghost); - -# endif -# endif - \} - \} -#endif // DEAL_II_WITH_P4EST - } diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index 8aa5670873..4739542ca4 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -4198,7 +4198,7 @@ namespace parallel // at that boundary. const std::map> vertices_with_ghost_neighbors = - this->compute_vertices_with_ghost_neighbors(); + GridTools::compute_vertices_with_ghost_neighbors(*this); // now collect cells and their vertices // for the interested neighbors diff --git a/source/distributed/tria_base.cc b/source/distributed/tria_base.cc index beecfd9486..e12ddedbf9 100644 --- a/source/distributed/tria_base.cc +++ b/source/distributed/tria_base.cc @@ -353,52 +353,7 @@ namespace parallel TriangulationBase::compute_vertices_with_ghost_neighbors() const { - // 1) collect for each vertex on periodic faces all vertices it coincides - // with - std::map> coinciding_vertex_groups; - std::map vertex_to_coinciding_vertex_group; - - GridTools::collect_coinciding_vertices(*this, - coinciding_vertex_groups, - vertex_to_coinciding_vertex_group); - - // 2) collect vertices belonging to local cells - std::vector vertex_of_own_cell(this->n_vertices(), false); - for (const auto &cell : this->active_cell_iterators()) - if (cell->is_locally_owned()) - for (const unsigned int v : GeometryInfo::vertex_indices()) - vertex_of_own_cell[cell->vertex_index(v)] = true; - - // 3) for each vertex belonging to a locally owned cell all ghost - // neighbors (including the periodic own) - std::map> result; - - // loop over all active ghost cells - for (const auto &cell : this->active_cell_iterators()) - if (cell->is_ghost()) - { - const types::subdomain_id owner = cell->subdomain_id(); - - // loop over all its vertices - for (const unsigned int v : GeometryInfo::vertex_indices()) - { - // set owner if vertex belongs to a local cell - if (vertex_of_own_cell[cell->vertex_index(v)]) - result[cell->vertex_index(v)].insert(owner); - - // mark also nodes coinciding due to periodicity - auto coinciding_vertex_group = - vertex_to_coinciding_vertex_group.find(cell->vertex_index(v)); - if (coinciding_vertex_group != - vertex_to_coinciding_vertex_group.end()) - for (auto coinciding_vertex : - coinciding_vertex_groups[coinciding_vertex_group->second]) - if (vertex_of_own_cell[coinciding_vertex]) - result[coinciding_vertex].insert(owner); - } - } - - return result; + return GridTools::compute_vertices_with_ghost_neighbors(*this); } diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc index a21950148c..e8078d8611 100644 --- a/source/dofs/dof_handler_policy.cc +++ b/source/dofs/dof_handler_policy.cc @@ -4550,7 +4550,7 @@ namespace internal // to exchange DoF indices const std::map> vertices_with_ghost_neighbors = - triangulation->compute_vertices_with_ghost_neighbors(); + GridTools::compute_vertices_with_ghost_neighbors(*triangulation); // mark all cells that either have to send data (locally // owned cells that are adjacent to ghost neighbors in some @@ -5102,7 +5102,8 @@ namespace internal const std::map> vertices_with_ghost_neighbors = - triangulation->compute_vertices_with_ghost_neighbors(); + GridTools::compute_vertices_with_ghost_neighbors( + *triangulation); // Send and receive cells. After this, only the local cells diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 377d3d7ba6..0a576d5e5f 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -5501,6 +5501,66 @@ namespace GridTools coinciding_vertex_groups[p.second].push_back(p.first); } } + + + + template + std::map> + compute_vertices_with_ghost_neighbors( + const Triangulation &tria) + { + if (dynamic_cast *>( + &tria) == nullptr) // nothing to do for a serial triangulation + return {}; + + // 1) collect for each vertex on periodic faces all vertices it coincides + // with + std::map> coinciding_vertex_groups; + std::map vertex_to_coinciding_vertex_group; + + GridTools::collect_coinciding_vertices(tria, + coinciding_vertex_groups, + vertex_to_coinciding_vertex_group); + + // 2) collect vertices belonging to local cells + std::vector vertex_of_own_cell(tria.n_vertices(), false); + for (const auto &cell : tria.active_cell_iterators()) + if (cell->is_locally_owned()) + for (const unsigned int v : GeometryInfo::vertex_indices()) + vertex_of_own_cell[cell->vertex_index(v)] = true; + + // 3) for each vertex belonging to a locally owned cell all ghost + // neighbors (including the periodic own) + std::map> result; + + // loop over all active ghost cells + for (const auto &cell : tria.active_cell_iterators()) + if (cell->is_ghost()) + { + const types::subdomain_id owner = cell->subdomain_id(); + + // loop over all its vertices + for (const unsigned int v : GeometryInfo::vertex_indices()) + { + // set owner if vertex belongs to a local cell + if (vertex_of_own_cell[cell->vertex_index(v)]) + result[cell->vertex_index(v)].insert(owner); + + // mark also nodes coinciding due to periodicity + auto coinciding_vertex_group = + vertex_to_coinciding_vertex_group.find(cell->vertex_index(v)); + if (coinciding_vertex_group != + vertex_to_coinciding_vertex_group.end()) + for (auto coinciding_vertex : + coinciding_vertex_groups[coinciding_vertex_group->second]) + if (vertex_of_own_cell[coinciding_vertex]) + result[coinciding_vertex].insert(owner); + } + } + + return result; + } + } /* namespace GridTools */ diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index 824c610b40..ac471d4370 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -410,6 +410,10 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) const Triangulation &, std::map> &, std::map &); + + template std::map> + compute_vertices_with_ghost_neighbors( + const Triangulation &); \} #endif } -- 2.39.5