From b35a2736487aa567defedb762b98dabf47db2051 Mon Sep 17 00:00:00 2001 From: David Wells Date: Sun, 5 Aug 2018 20:50:12 -0400 Subject: [PATCH] Add Triangulation::active_face_iterators(). This function is like Triangulation::active_cell_iterators(), but for faces. --- doc/news/changes/minor/20180805DavidWells | 4 + include/deal.II/grid/tria.h | 21 +++++ source/grid/tria.cc | 10 +++ tests/dofs/range_based_for_tria_face.cc | 58 ++++++++++++++ tests/dofs/range_based_for_tria_face.output | 4 + tests/grid/filtered_iterator_06.cc | 88 +++++++++++++++++++++ tests/grid/filtered_iterator_06.output | 2 + 7 files changed, 187 insertions(+) create mode 100644 doc/news/changes/minor/20180805DavidWells create mode 100644 tests/dofs/range_based_for_tria_face.cc create mode 100644 tests/dofs/range_based_for_tria_face.output create mode 100644 tests/grid/filtered_iterator_06.cc create mode 100644 tests/grid/filtered_iterator_06.output diff --git a/doc/news/changes/minor/20180805DavidWells b/doc/news/changes/minor/20180805DavidWells new file mode 100644 index 0000000000..9d133e018b --- /dev/null +++ b/doc/news/changes/minor/20180805DavidWells @@ -0,0 +1,4 @@ +New: A convenience function Triangulation::active_face_iterators() (which is like +Triangulation::active_cell_iterators(), but for faces) has been added. +
+(David Wells, 2018/08/05) diff --git a/include/deal.II/grid/tria.h b/include/deal.II/grid/tria.h index bbba763462..fea118a11a 100644 --- a/include/deal.II/grid/tria.h +++ b/include/deal.II/grid/tria.h @@ -2823,6 +2823,27 @@ public: face_iterator end_face() const; + /** + * Return an iterator range that contains all active faces that make up this + * triangulation. This function is the face version of + * Triangulation::active_cell_iterators(), and allows one to write code + * like, e.g., + * + * @code + * Triangulation triangulation; + * ... + * for (auto &face : triangulation.active_face_iterators()) + * face->set_manifold_id(42); + * @endcode + * + * @return The half open range [this->begin_active_face(), + * this->end_face()) + * + * @ingroup CPP11 + */ + IteratorRange + active_face_iterators() const; + /* * @} */ diff --git a/source/grid/tria.cc b/source/grid/tria.cc index e2325b8b61..947a120e0b 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -12153,6 +12153,16 @@ Triangulation::end_face() const } + +template +IteratorRange::active_face_iterator> +Triangulation::active_face_iterators() const +{ + return IteratorRange< + typename Triangulation::active_face_iterator>( + begin_active_face(), end_face()); +} + /*------------------------ Vertex iterator functions ------------------------*/ diff --git a/tests/dofs/range_based_for_tria_face.cc b/tests/dofs/range_based_for_tria_face.cc new file mode 100644 index 0000000000..6ce8a7e44b --- /dev/null +++ b/tests/dofs/range_based_for_tria_face.cc @@ -0,0 +1,58 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 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. +// +// --------------------------------------------------------------------- + + + +// Check face range-based for loops for triangulations + +#include +#include +#include +#include + +#include + +#include "../tests.h" + + +template +void +check() +{ + Triangulation tr; + GridGenerator::hyper_cube(tr); + tr.refine_global(2); + + for (auto &face : tr.active_face_iterators()) + face->set_manifold_id(42); + + for (const auto &cell : tr.active_cell_iterators()) + for (unsigned int face_n = 0; face_n < GeometryInfo::faces_per_cell; + ++face_n) + Assert(cell->face(face_n)->manifold_id() == 42, ExcInternalError()); + + deallog << "OK" << std::endl; +} + + +int +main() +{ + initlog(); + + check<2, 2>(); + check<2, 3>(); + check<3, 3>(); +} diff --git a/tests/dofs/range_based_for_tria_face.output b/tests/dofs/range_based_for_tria_face.output new file mode 100644 index 0000000000..fb71de2867 --- /dev/null +++ b/tests/dofs/range_based_for_tria_face.output @@ -0,0 +1,4 @@ + +DEAL::OK +DEAL::OK +DEAL::OK diff --git a/tests/grid/filtered_iterator_06.cc b/tests/grid/filtered_iterator_06.cc new file mode 100644 index 0000000000..0ff0657e0c --- /dev/null +++ b/tests/grid/filtered_iterator_06.cc @@ -0,0 +1,88 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 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. +// +// --------------------------------------------------------------------- + +// check that we can create filtered face iterators from a face iterator +// range. Similar to filtered_iterator_05. + +#include +#include +#include +#include +#include + +#include + +#include "../tests.h" + +template +void +test() +{ + Triangulation tria; + GridGenerator::hyper_cube(tria); + // refine the boundary cells a few times + for (unsigned int i = 0; i < 5; ++i) + { + for (auto &cell : tria.active_cell_iterators()) + cell->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + } + + // we now have a number of cells, flag them with some manifold ids + std::size_t i = 0; + for (auto &cell : tria.active_cell_iterators()) + { + cell->set_all_manifold_ids(i % 3); + ++i; + } + + // Count the faces that are on the boundary and have a manifold_id of 0 + const types::manifold_id manifold_id = 0; + std::set::active_face_iterator> + face_set; + for (const auto &cell : tria.active_cell_iterators()) + for (unsigned int face_n = 0; face_n < GeometryInfo::faces_per_cell; + ++face_n) + if (cell->face(face_n)->at_boundary() && + cell->face(face_n)->manifold_id() == manifold_id) + face_set.insert(cell->face(face_n)); + + std::size_t n_filtered_cells = 0; + for (const auto filtered_cell : filter_iterators( + tria.active_face_iterators(), + IteratorFilters::AtBoundary(), + [=](const typename Triangulation::active_face_iterator + &face) { return face->manifold_id() == manifold_id; })) + { + AssertThrow(face_set.count(filtered_cell) == 1, + ExcMessage("Wrong cell filtered.")); + ++n_filtered_cells; + } + AssertThrow(n_filtered_cells == face_set.size(), + ExcMessage("Filtered cells missing.")); +} + +int +main() +{ + initlog(); + deallog << std::setprecision(4); + + test<2, 2>(); + test<2, 3>(); + test<3, 3>(); + + deallog << "OK" << std::endl; +} diff --git a/tests/grid/filtered_iterator_06.output b/tests/grid/filtered_iterator_06.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/grid/filtered_iterator_06.output @@ -0,0 +1,2 @@ + +DEAL::OK -- 2.39.5