From 2236087c5807a0841f1fabc194f49a02b690958a Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Tue, 4 May 2021 11:44:34 +0200 Subject: [PATCH] Implement cache to speed up TriaAccessor::vertex_index() --- include/deal.II/grid/tria.h | 7 +++++ .../deal.II/grid/tria_accessor.templates.h | 22 +++++++++++++--- include/deal.II/grid/tria_levels.h | 14 +++++++--- source/grid/tria.cc | 26 +++++++++++++++++++ 4 files changed, 63 insertions(+), 6 deletions(-) diff --git a/include/deal.II/grid/tria.h b/include/deal.II/grid/tria.h index 9fbe2107da..deea012397 100644 --- a/include/deal.II/grid/tria.h +++ b/include/deal.II/grid/tria.h @@ -3897,6 +3897,12 @@ private: void reset_global_cell_indices(); + /** + * Reset cache for the cells' vertex indices. + */ + void + reset_cell_vertex_indices_cache(); + /** * Refine all cells on all levels which were previously flagged for * refinement. @@ -4255,6 +4261,7 @@ Triangulation::load(Archive &ar, const unsigned int) level->global_active_cell_indices.resize(level->refine_flags.size()); level->global_level_cell_indices.resize(level->refine_flags.size()); } + reset_cell_vertex_indices_cache(); reset_active_cell_indices(); reset_global_cell_indices(); } diff --git a/include/deal.II/grid/tria_accessor.templates.h b/include/deal.II/grid/tria_accessor.templates.h index 1524e05e9f..5e249501a2 100644 --- a/include/deal.II/grid/tria_accessor.templates.h +++ b/include/deal.II/grid/tria_accessor.templates.h @@ -1112,8 +1112,23 @@ TriaAccessor::vertex_index( { AssertIndexRange(corner, this->n_vertices()); - return dealii::internal::TriaAccessorImplementation::Implementation:: - vertex_index(*this, corner); + if (structdim == dim) + { + constexpr unsigned int max_vertices_per_cell = 1 << dim; + const std::size_t my_index = + static_cast(this->present_index) * max_vertices_per_cell; + AssertIndexRange(my_index + corner, + this->tria->levels[this->present_level] + ->cell_vertex_indices_cache.size()); + const unsigned int vertex_index = + this->tria->levels[this->present_level] + ->cell_vertex_indices_cache[my_index + corner]; + Assert(vertex_index != numbers::invalid_unsigned_int, ExcInternalError()); + return vertex_index; + } + else + return dealii::internal::TriaAccessorImplementation::Implementation:: + vertex_index(*this, corner); } @@ -3206,7 +3221,8 @@ namespace internal ((i == 1) && cell.at_boundary(1) ? dealii::TriaAccessor<0, 1, spacedim>::right_vertex : dealii::TriaAccessor<0, 1, spacedim>::interior_vertex)), - cell.vertex_index(i)); + dealii::internal::TriaAccessorImplementation::Implementation:: + vertex_index(cell, i)); return dealii::TriaIterator>(a); } diff --git a/include/deal.II/grid/tria_levels.h b/include/deal.II/grid/tria_levels.h index 39969f954e..64f774477f 100644 --- a/include/deal.II/grid/tria_levels.h +++ b/include/deal.II/grid/tria_levels.h @@ -224,6 +224,14 @@ namespace internal */ std::vector reference_cell; + /** + * A cache for the vertex indices of the cells (`structdim == dim`), in + * order to more quickly retrieve these frequently accessed quantities. + * For simplified addressing, the information is indexed by the maximum + * number of vertices possible for a cell (quadrilateral/hexahedron). + */ + std::vector cell_vertex_indices_cache; + /** * Determine an estimate for the memory consumption (in bytes) of this * object. @@ -250,9 +258,9 @@ namespace internal ar &refine_flags &coarsen_flags; - // do not serialize 'active_cell_indices' here. instead of storing them - // to the stream and re-reading them again later, we just rebuild them - // in Triangulation::load() + // do not serialize `active_cell_indices` and `vertex_indices_cache` + // here. instead of storing them to the stream and re-reading them again + // later, we just rebuild them in Triangulation::load() ar &neighbors; ar &subdomain_ids; diff --git a/source/grid/tria.cc b/source/grid/tria.cc index a5a42d27c1..4f8c7d1b83 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -10539,6 +10539,7 @@ Triangulation::create_triangulation( // update our counts of the various elements of a triangulation, and set // active_cell_indices of all cells + reset_cell_vertex_indices_cache(); internal::TriangulationImplementation::Implementation::compute_number_cache( *this, levels.size(), number_cache); reset_active_cell_indices(); @@ -13390,6 +13391,8 @@ Triangulation::execute_coarsening_and_refinement() const DistortedCellList cells_with_distorted_children = execute_refinement(); + reset_cell_vertex_indices_cache(); + // verify a case with which we have had // some difficulty in the past (see the // deal.II/coarsening_* tests) @@ -13453,6 +13456,29 @@ Triangulation::reset_global_cell_indices() +template +void +Triangulation::reset_cell_vertex_indices_cache() +{ + for (unsigned int l = 0; l < levels.size(); ++l) + { + constexpr unsigned int max_vertices_per_cell = 1 << dim; + std::vector &cache = levels[l]->cell_vertex_indices_cache; + cache.clear(); + cache.resize(levels[l]->refine_flags.size() * max_vertices_per_cell, + numbers::invalid_unsigned_int); + for (const auto &cell : cell_iterators_on_level(l)) + { + const unsigned int my_index = cell->index() * max_vertices_per_cell; + for (const unsigned int i : cell->vertex_indices()) + cache[my_index + i] = internal::TriaAccessorImplementation:: + Implementation::vertex_index(*cell, i); + } + } +} + + + template void Triangulation::update_periodic_face_map() -- 2.39.5