From a0137a6e6b3a7ceef164582007ff4ca6f2ad472b Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Wed, 18 Sep 2019 07:47:25 +0200 Subject: [PATCH] Clean up compute_vertices_with_ghost_neighbors --- .../distributed/fully_distributed_tria.h | 7 - include/deal.II/distributed/tria.h | 17 --- source/distributed/fully_distributed_tria.cc | 124 ------------------ source/distributed/tria.cc | 26 +--- source/distributed/tria_base.cc | 105 +++++++++++++-- 5 files changed, 99 insertions(+), 180 deletions(-) diff --git a/include/deal.II/distributed/fully_distributed_tria.h b/include/deal.II/distributed/fully_distributed_tria.h index eef0b55d2d..4c7843278c 100644 --- a/include/deal.II/distributed/fully_distributed_tria.h +++ b/include/deal.II/distributed/fully_distributed_tria.h @@ -327,13 +327,6 @@ namespace parallel virtual std::size_t memory_consumption() const override; - /** - * This function determines the neighboring subdomains that are adjacent - * to each vertex. - */ - virtual std::map> - compute_vertices_with_ghost_neighbors() const override; - virtual bool is_multilevel_hierarchy_constructed() const override; diff --git a/include/deal.II/distributed/tria.h b/include/deal.II/distributed/tria.h index 467a21489c..8376b466f6 100644 --- a/include/deal.II/distributed/tria.h +++ b/include/deal.II/distributed/tria.h @@ -1193,16 +1193,6 @@ namespace parallel std::vector get_cell_weights() const; - /** - * Override the implementation in parallel::TriangulationBase because - * we can ask p4est about ghost neighbors across periodic boundaries. - * - * Specifically, this function determines the neighboring subdomains that - * are adjacent to each vertex. - */ - virtual std::map> - compute_vertices_with_ghost_neighbors() const override; - /** * This method returns a bit vector of length tria.n_vertices() * indicating the locally active vertices on a level, i.e., the vertices @@ -1368,13 +1358,6 @@ namespace parallel */ Settings settings; - /** - * Like above, this method, which is only implemented for dim = 2 or 3, - * needs a stub because it is used in dof_handler_policy.cc - */ - virtual std::map> - compute_vertices_with_ghost_neighbors() const override; - /** * Like above, this method, which is only implemented for dim = 2 or 3, * needs a stub because it is used in dof_handler_policy.cc diff --git a/source/distributed/fully_distributed_tria.cc b/source/distributed/fully_distributed_tria.cc index f4fe3720d4..b6625bc541 100644 --- a/source/distributed/fully_distributed_tria.cc +++ b/source/distributed/fully_distributed_tria.cc @@ -346,130 +346,6 @@ namespace parallel - template - std::map> - Triangulation::compute_vertices_with_ghost_neighbors() const - { - // 1) collect for each vertex on periodic faces all vertices it coincides - // with - std::map> - vertex_to_coinciding_vertices; - { - // 1a) collect nodes coinciding due to periodicity - std::map vertex_to_coinciding_vertex; - for (auto &cell : this->active_cell_iterators()) - if (cell->is_locally_owned() || cell->is_ghost()) - for (auto i = 0u; i < GeometryInfo::faces_per_cell; ++i) - if (cell->has_periodic_neighbor(i) && - cell->periodic_neighbor(i)->active()) - { - auto face_t = cell->face(i); - auto face_n = cell->periodic_neighbor(i)->face( - cell->periodic_neighbor_face_no(i)); - for (auto j = 0u; j < GeometryInfo::vertices_per_face; - ++j) - { - auto v_t = face_t->vertex_index(j); - auto v_n = face_n->vertex_index(j); - unsigned int temp = std::min(v_t, v_n); - { - auto it = vertex_to_coinciding_vertex.find(v_t); - if (it != vertex_to_coinciding_vertex.end()) - temp = std::min(temp, it->second); - } - { - auto it = vertex_to_coinciding_vertex.find(v_n); - if (it != vertex_to_coinciding_vertex.end()) - temp = std::min(temp, it->second); - } - vertex_to_coinciding_vertex[v_t] = temp; - vertex_to_coinciding_vertex[v_n] = temp; - } - } - - // 1b) compress map: let vertices point to the coinciding vertex with - // the smallest index - for (auto &p : vertex_to_coinciding_vertex) - { - if (p.first == p.second) - continue; - unsigned int temp = p.second; - while (temp != vertex_to_coinciding_vertex[temp]) - temp = vertex_to_coinciding_vertex[temp]; - p.second = temp; - } - -#ifdef DEBUG - // check if map is actually compressed - for (auto p : vertex_to_coinciding_vertex) - { - if (p.first == p.second) - continue; - auto pp = vertex_to_coinciding_vertex.find(p.second); - if (pp->first == pp->second) - continue; - AssertThrow(false, ExcMessage("Map has to be compressed!")); - } -#endif - - // 1c) create a map: smallest index of coinciding index -> all - // coinciding indices - std::map> - smallest_coinciding_vertex_to_coinciding_vertices; - for (auto p : vertex_to_coinciding_vertex) - smallest_coinciding_vertex_to_coinciding_vertices[p.second] = - std::set(); - - for (auto p : vertex_to_coinciding_vertex) - smallest_coinciding_vertex_to_coinciding_vertices[p.second].insert( - p.first); - - // 1d) create a map: vertex -> all coinciding indices - for (auto &s : smallest_coinciding_vertex_to_coinciding_vertices) - for (auto &ss : s.second) - vertex_to_coinciding_vertices[ss] = s.second; - } - - // 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 (auto v = 0u; v < GeometryInfo::vertices_per_cell; ++v) - vertex_of_own_cell[cell->vertex_index(v)] = true; - - // 3) for 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 auto owner = cell->subdomain_id(); - - // loop over all its vertices - for (unsigned int v = 0; v < GeometryInfo::vertices_per_cell; - ++v) - { - // 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_vertices = - vertex_to_coinciding_vertices.find(cell->vertex_index(v)); - if (coinciding_vertices != vertex_to_coinciding_vertices.end()) - for (auto coinciding_vertex : coinciding_vertices->second) - if (vertex_of_own_cell[coinciding_vertex]) - result[coinciding_vertex].insert(owner); - } - } - - return result; - } - - - template bool Triangulation::is_multilevel_hierarchy_constructed() const diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index aef7015732..528e8234df 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -4194,7 +4194,8 @@ namespace parallel // Here, it is sufficient to collect all vertices that are located // at that boundary. const std::map> - vertices_with_ghost_neighbors = compute_vertices_with_ghost_neighbors(); + vertices_with_ghost_neighbors = + this->compute_vertices_with_ghost_neighbors(); // now collect cells and their vertices // for the interested neighbors @@ -4470,19 +4471,6 @@ namespace parallel - template - std::map> - Triangulation::compute_vertices_with_ghost_neighbors() const - { - Assert(dim > 1, ExcNotImplemented()); - - return dealii::internal::p4est::compute_vertices_with_ghost_neighbors< - dim, - spacedim>(*this, this->parallel_forest, this->parallel_ghost); - } - - - template std::vector Triangulation::mark_locally_active_vertices_on_level( @@ -5035,16 +5023,6 @@ namespace parallel - template - std::map> - Triangulation<1, spacedim>::compute_vertices_with_ghost_neighbors() const - { - Assert(false, ExcNotImplemented()); - return std::map>(); - } - - - template std::map> Triangulation<1, spacedim>::compute_level_vertices_with_ghost_neighbors( diff --git a/source/distributed/tria_base.cc b/source/distributed/tria_base.cc index 44252fa12f..283029a696 100644 --- a/source/distributed/tria_base.cc +++ b/source/distributed/tria_base.cc @@ -362,36 +362,125 @@ namespace parallel TriangulationBase::compute_vertices_with_ghost_neighbors() const { - // TODO: we are not treating periodic neighbors correctly here. If we do - // we can remove the overriding implementation for p::d::Triangulation - // that is currently using a p4est callback to get correct ghost neighbors - // over periodic faces. - Assert(this->get_periodic_face_map().size() == 0, ExcNotImplemented()); + // 1) collect for each vertex on periodic faces all vertices it coincides + // with + std::map> + vertex_to_coinciding_vertices; + { + // 1a) collect nodes coinciding due to periodicity + std::map vertex_to_coinciding_vertex; + for (auto &cell : this->active_cell_iterators()) + if (cell->is_locally_owned() || cell->is_ghost()) + for (auto i = 0u; i < GeometryInfo::faces_per_cell; ++i) + if (cell->has_periodic_neighbor(i) && + cell->periodic_neighbor(i)->active()) + { + auto face_t = cell->face(i); + auto face_n = cell->periodic_neighbor(i)->face( + cell->periodic_neighbor_face_no(i)); + for (auto j = 0u; j < GeometryInfo::vertices_per_face; ++j) + { + auto v_t = face_t->vertex_index(j); + auto v_n = face_n->vertex_index(j); + unsigned int temp = std::min(v_t, v_n); + { + auto it = vertex_to_coinciding_vertex.find(v_t); + if (it != vertex_to_coinciding_vertex.end()) + temp = std::min(temp, it->second); + } + { + auto it = vertex_to_coinciding_vertex.find(v_n); + if (it != vertex_to_coinciding_vertex.end()) + temp = std::min(temp, it->second); + } + vertex_to_coinciding_vertex[v_t] = temp; + vertex_to_coinciding_vertex[v_n] = temp; + } + } + + // 1b) compress map: let vertices point to the coinciding vertex with + // the smallest index + for (auto &p : vertex_to_coinciding_vertex) + { + if (p.first == p.second) + continue; + unsigned int temp = p.second; + while (temp != vertex_to_coinciding_vertex[temp]) + temp = vertex_to_coinciding_vertex[temp]; + p.second = temp; + } +#ifdef DEBUG + // check if map is actually compressed + for (auto p : vertex_to_coinciding_vertex) + { + if (p.first == p.second) + continue; + auto pp = vertex_to_coinciding_vertex.find(p.second); + if (pp->first == pp->second) + continue; + AssertThrow(false, ExcMessage("Map has to be compressed!")); + } +#endif - std::vector vertex_of_own_cell(this->n_vertices(), false); + // 1c) create a map: smallest index of coinciding index -> all + // coinciding indices + std::map> + smallest_coinciding_vertex_to_coinciding_vertices; + for (auto p : vertex_to_coinciding_vertex) + smallest_coinciding_vertex_to_coinciding_vertices[p.second] = + std::set(); + + for (auto p : vertex_to_coinciding_vertex) + smallest_coinciding_vertex_to_coinciding_vertices[p.second].insert( + p.first); + + // 1d) create a map: vertex -> all coinciding indices + for (auto &s : smallest_coinciding_vertex_to_coinciding_vertices) + for (auto &ss : s.second) + vertex_to_coinciding_vertices[ss] = s.second; + } + // 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 (unsigned int v = 0; v < GeometryInfo::vertices_per_cell; ++v) + for (auto v = 0u; v < GeometryInfo::vertices_per_cell; ++v) vertex_of_own_cell[cell->vertex_index(v)] = true; + // 3) for 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(); + const auto owner = cell->subdomain_id(); + + // loop over all its vertices for (unsigned int v = 0; v < GeometryInfo::vertices_per_cell; ++v) { + // 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_vertices = + vertex_to_coinciding_vertices.find(cell->vertex_index(v)); + if (coinciding_vertices != vertex_to_coinciding_vertices.end()) + for (auto coinciding_vertex : coinciding_vertices->second) + if (vertex_of_own_cell[coinciding_vertex]) + result[coinciding_vertex].insert(owner); } } return result; } + + template DistributedTriangulationBase::DistributedTriangulationBase( MPI_Comm mpi_communicator, -- 2.39.5