From 876862db9b80b2ad18d178297e74bac187f7799f Mon Sep 17 00:00:00 2001 From: David Wells Date: Fri, 6 Sep 2024 18:52:00 -0400 Subject: [PATCH] Triangulation::get_boundary_ids(): don't return duplicates. In 1d we support, in spacedim == 2 and spacedim == 3, more than two boundary nodes, so we need to remove duplicates from this vector. --- doc/news/changes/minor/20240906Wells | 4 ++++ source/grid/tria.cc | 32 +++++++--------------------- tests/grid/grid_in_vtk_1d_3d.cc | 6 ++++++ tests/grid/grid_in_vtk_1d_3d.output | 4 ++++ 4 files changed, 22 insertions(+), 24 deletions(-) create mode 100644 doc/news/changes/minor/20240906Wells diff --git a/doc/news/changes/minor/20240906Wells b/doc/news/changes/minor/20240906Wells new file mode 100644 index 0000000000..e80e66bc09 --- /dev/null +++ b/doc/news/changes/minor/20240906Wells @@ -0,0 +1,4 @@ +Fixed: Triangulation::get_boundary_ids() does not return duplicates when +`dim == 1` any more. +
+(David Wells, 2024/09/06) diff --git a/source/grid/tria.cc b/source/grid/tria.cc index 446727d20e..a496a30a66 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -12331,30 +12331,14 @@ DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim)) std::vector Triangulation::get_boundary_ids() const { - // in 1d, we store a map of all used boundary indicators. use it for - // our purposes - if (dim == 1) - { - std::vector boundary_ids; - for (std::map::const_iterator p = - vertex_to_boundary_id_map_1d->begin(); - p != vertex_to_boundary_id_map_1d->end(); - ++p) - boundary_ids.push_back(p->second); - - return boundary_ids; - } - else - { - std::set b_ids; - for (auto cell : active_cell_iterators()) - if (cell->is_locally_owned()) - for (const unsigned int face : cell->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; - } + std::set boundary_ids; + for (const auto &cell : active_cell_iterators()) + if (cell->is_locally_owned()) + for (const auto &face : cell->face_indices()) + if (cell->at_boundary(face)) + boundary_ids.insert(cell->face(face)->boundary_id()); + + return {boundary_ids.begin(), boundary_ids.end()}; } diff --git a/tests/grid/grid_in_vtk_1d_3d.cc b/tests/grid/grid_in_vtk_1d_3d.cc index d8eef4a0e1..cfddaeb55c 100644 --- a/tests/grid/grid_in_vtk_1d_3d.cc +++ b/tests/grid/grid_in_vtk_1d_3d.cc @@ -58,6 +58,9 @@ main() } } + for (const auto &boundary_id : triangulation1.get_boundary_ids()) + deallog << " boundary_id = " << boundary_id << std::endl; + deallog << "Triangulation 2:\n"; for (const auto &cell : triangulation2.active_cell_iterators()) { @@ -69,6 +72,9 @@ main() } } + for (const auto &boundary_id : triangulation1.get_boundary_ids()) + deallog << " boundary_id = " << boundary_id << std::endl; + GridGenerator::merge_triangulations(triangulation2, triangulation1, triangulation2); diff --git a/tests/grid/grid_in_vtk_1d_3d.output b/tests/grid/grid_in_vtk_1d_3d.output index f1d25018eb..8f505ac215 100644 --- a/tests/grid/grid_in_vtk_1d_3d.output +++ b/tests/grid/grid_in_vtk_1d_3d.output @@ -12,6 +12,8 @@ DEAL:: vertex: 0.00000 2.01000 -1.00000 DEAL:: cell=0.3 DEAL:: vertex: 0.00000 2.01000 0.00000 DEAL:: vertex: 0.00000 3.01000 0.00000 +DEAL:: boundary_id = 0 +DEAL:: boundary_id = 1 DEAL::Triangulation 2: cell=0.0 DEAL:: vertex: 0.00000 0.00000 0.00000 @@ -25,6 +27,8 @@ DEAL:: vertex: 0.00000 0.00000 -1.00000 DEAL:: cell=0.3 DEAL:: vertex: 0.00000 0.00000 0.00000 DEAL:: vertex: 0.00000 1.00000 0.00000 +DEAL:: boundary_id = 0 +DEAL:: boundary_id = 1 0.00000 0.00000 0.00000 0 0 0.00000 0.00000 1.00000 0 0 -- 2.39.5