From e32c19412545693869d3f618348e041eb13765f7 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sat, 4 Apr 2020 00:30:24 +0200 Subject: [PATCH] Fix get_boundary_ids() and get_manifold_ids() for parallel triangulations --- include/deal.II/distributed/tria_base.h | 18 ++++++ include/deal.II/grid/tria.h | 18 +++--- source/distributed/tria_base.cc | 22 +++++++ source/grid/tria.cc | 28 ++++----- .../distributed_grids/get_boundary_ids_01.cc | 58 +++++++++++++++++++ .../get_boundary_ids_01.mpirun=2.output | 5 ++ 6 files changed, 127 insertions(+), 22 deletions(-) create mode 100644 tests/distributed_grids/get_boundary_ids_01.cc create mode 100644 tests/distributed_grids/get_boundary_ids_01.mpirun=2.output diff --git a/include/deal.II/distributed/tria_base.h b/include/deal.II/distributed/tria_base.h index 458c72699c..0203f6cc77 100644 --- a/include/deal.II/distributed/tria_base.h +++ b/include/deal.II/distributed/tria_base.h @@ -211,6 +211,24 @@ namespace parallel virtual std::map> compute_vertices_with_ghost_neighbors() const; + /** + * @copydoc dealii::Triangulation::get_boundary_ids() + * + * @note This function involves a global communication gathering all current + * IDs from all processes. + */ + virtual std::vector + get_boundary_ids() const override; + + /** + * @copydoc dealii::Triangulation::get_manifold_ids() + * + * @note This function involves a global communication gathering all current + * IDs from all processes. + */ + virtual std::vector + get_manifold_ids() const override; + protected: /** * MPI communicator to be used for the triangulation. We create a unique diff --git a/include/deal.II/grid/tria.h b/include/deal.II/grid/tria.h index 5866fd2188..0c40ea2ef4 100644 --- a/include/deal.II/grid/tria.h +++ b/include/deal.II/grid/tria.h @@ -1708,30 +1708,32 @@ public: /** * Return a vector containing all boundary indicators assigned to boundary - * faces of this Triangulation object. Note, that each boundary indicator is - * reported only once. The size of the return vector will represent the - * number of different indicators (which is greater or equal one). + * faces of active cells of this Triangulation object. Note, that each + * boundary indicator is reported only once. The size of the return vector + * will represent the number of different indicators (which is greater or + * equal one). * * @ingroup boundary * * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - std::vector + virtual std::vector get_boundary_ids() const; /** * Return a vector containing all manifold indicators assigned to the - * objects of this Triangulation. Note, that each manifold indicator is - * reported only once. The size of the return vector will represent the - * number of different indicators (which is greater or equal one). + * objects of the active cells of this Triangulation. Note, that each + * manifold indicator is reported only once. The size of the return vector + * will represent the number of different indicators (which is greater or + * equal one). * * @ingroup manifold * * @see * @ref GlossManifoldIndicator "Glossary entry on manifold indicators" */ - std::vector + virtual std::vector get_manifold_ids() const; /** diff --git a/source/distributed/tria_base.cc b/source/distributed/tria_base.cc index 7117d6f2ac..beecfd9486 100644 --- a/source/distributed/tria_base.cc +++ b/source/distributed/tria_base.cc @@ -403,6 +403,28 @@ namespace parallel + template + std::vector + TriangulationBase::get_boundary_ids() const + { + return Utilities::MPI::compute_set_union( + dealii::Triangulation::get_boundary_ids(), + this->mpi_communicator); + } + + + + template + std::vector + TriangulationBase::get_manifold_ids() const + { + return Utilities::MPI::compute_set_union( + dealii::Triangulation::get_manifold_ids(), + this->mpi_communicator); + } + + + template DistributedTriangulationBase::DistributedTriangulationBase( MPI_Comm mpi_communicator, diff --git a/source/grid/tria.cc b/source/grid/tria.cc index 98f8d49917..3bc8fa488a 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -10397,11 +10397,11 @@ Triangulation::get_boundary_ids() const else { std::set b_ids; - active_cell_iterator cell = begin_active(); - for (; cell != end(); ++cell) - for (const unsigned int face : GeometryInfo::face_indices()) - if (cell->at_boundary(face)) - b_ids.insert(cell->face(face)->boundary_id()); + for (auto cell : active_cell_iterators()) + if (cell->is_locally_owned()) + for (const unsigned int face : GeometryInfo::face_indices()) + if (cell->at_boundary(face)) + b_ids.insert(cell->face(face)->boundary_id()); std::vector boundary_ids(b_ids.begin(), b_ids.end()); return boundary_ids; } @@ -10414,15 +10414,15 @@ std::vector Triangulation::get_manifold_ids() const { std::set m_ids; - active_cell_iterator cell = begin_active(); - for (; cell != end(); ++cell) - { - m_ids.insert(cell->manifold_id()); - if (dim > 1) - for (const unsigned int face : GeometryInfo::face_indices()) - if (cell->at_boundary(face)) - m_ids.insert(cell->face(face)->manifold_id()); - } + for (auto cell : active_cell_iterators()) + if (cell->is_locally_owned()) + { + m_ids.insert(cell->manifold_id()); + if (dim > 1) + for (const unsigned int face : GeometryInfo::face_indices()) + if (cell->at_boundary(face)) + m_ids.insert(cell->face(face)->manifold_id()); + } std::vector manifold_indicators(m_ids.begin(), m_ids.end()); return manifold_indicators; diff --git a/tests/distributed_grids/get_boundary_ids_01.cc b/tests/distributed_grids/get_boundary_ids_01.cc new file mode 100644 index 0000000000..3020c84b20 --- /dev/null +++ b/tests/distributed_grids/get_boundary_ids_01.cc @@ -0,0 +1,58 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// --------------------------------------------------------------------- + + +// create a tria mesh and copy it + +#include + +#include + +#include "../tests.h" + +template +void +test() +{ + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + GridGenerator::subdivided_hyper_cube(tria, 4); + + for (auto cell : tria.active_cell_iterators()) + if (cell->is_locally_owned()) + for (const unsigned int face : GeometryInfo::face_indices()) + { + if (std::fabs(cell->face(face)->center()(0) - 0.0) < 1e-12) + cell->face(face)->set_all_boundary_ids(1); + if (std::fabs(cell->face(face)->center()(0) - 1.0) < 1e-12) + cell->face(face)->set_all_boundary_ids(2); + if (std::fabs(cell->face(face)->center()(1) - 0.0) < 1e-12) + cell->face(face)->set_all_boundary_ids(3); + if (std::fabs(cell->face(face)->center()(1) - 1.0) < 1e-12) + cell->face(face)->set_all_boundary_ids(4); + } + + for (auto i : tria.get_boundary_ids()) + deallog << i << " "; + deallog << std::endl; +} + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + test<2>(); +} diff --git a/tests/distributed_grids/get_boundary_ids_01.mpirun=2.output b/tests/distributed_grids/get_boundary_ids_01.mpirun=2.output new file mode 100644 index 0000000000..8bdee82a30 --- /dev/null +++ b/tests/distributed_grids/get_boundary_ids_01.mpirun=2.output @@ -0,0 +1,5 @@ + +DEAL:0::1 2 3 4 + +DEAL:1::1 2 3 4 + -- 2.39.5