Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
grid_tools_cache.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2017 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
18
22
24
25namespace GridTools
26{
27 template <int dim, int spacedim>
29 const Mapping<dim, spacedim> & mapping)
30 : update_flags(update_all)
31 , tria(&tria)
32 , mapping(&mapping)
33 {
35 tria.signals.any_change.connect([&]() { mark_for_update(update_all); });
36 }
37
38 template <int dim, int spacedim>
40 {
41 // Make sure that the signals that was attached to the triangulation
42 // is removed here.
43 if (tria_signal.connected())
44 tria_signal.disconnect();
45 }
46
47
48
49 template <int dim, int spacedim>
50 void
52 {
53 update_flags |= flags;
54 }
55
56
57
58 template <int dim, int spacedim>
59 const std::vector<
60 std::set<typename Triangulation<dim, spacedim>::active_cell_iterator>> &
62 {
63 if (update_flags & update_vertex_to_cell_map)
64 {
65 vertex_to_cells = GridTools::vertex_to_cell_map(*tria);
66 update_flags = update_flags & ~update_vertex_to_cell_map;
67 }
68 return vertex_to_cells;
69 }
70
71
72
73 template <int dim, int spacedim>
74 const std::vector<std::vector<Tensor<1, spacedim>>> &
76 {
78 {
79 vertex_to_cell_centers = GridTools::vertex_to_cell_centers_directions(
80 *tria, get_vertex_to_cell_map());
81 update_flags = update_flags & ~update_vertex_to_cell_centers_directions;
82 }
83 return vertex_to_cell_centers;
84 }
85
86
87
88 template <int dim, int spacedim>
89 const std::map<unsigned int, Point<spacedim>> &
91 {
92 if (update_flags & update_used_vertices)
93 {
94 used_vertices = GridTools::extract_used_vertices(*tria, *mapping);
95 update_flags = update_flags & ~update_used_vertices;
96 }
97 return used_vertices;
98 }
99
100
101
102 template <int dim, int spacedim>
103 const RTree<std::pair<Point<spacedim>, unsigned int>> &
105 {
106 if (update_flags & update_used_vertices_rtree)
107 {
108 const auto &used_vertices = get_used_vertices();
109 std::vector<std::pair<Point<spacedim>, unsigned int>> vertices(
110 used_vertices.size());
111 unsigned int i = 0;
112 for (const auto &it : used_vertices)
113 vertices[i++] = std::make_pair(it.second, it.first);
114 used_vertices_rtree = pack_rtree(vertices);
115 update_flags = update_flags & ~update_used_vertices_rtree;
116 }
117 return used_vertices_rtree;
118 }
120
121
122 template <int dim, int spacedim>
123 const RTree<
124 std::pair<BoundingBox<spacedim>,
127 {
128 if (update_flags & update_cell_bounding_boxes_rtree)
129 {
130 std::vector<std::pair<
133 boxes;
134 boxes.reserve(tria->n_active_cells());
135 for (const auto &cell : tria->active_cell_iterators())
136 boxes.emplace_back(mapping->get_bounding_box(cell), cell);
137
138 cell_bounding_boxes_rtree = pack_rtree(boxes);
139 update_flags = update_flags & ~update_cell_bounding_boxes_rtree;
140 }
141 return cell_bounding_boxes_rtree;
142 }
143
145
146 template <int dim, int spacedim>
147 const RTree<
148 std::pair<BoundingBox<spacedim>,
151 {
153 {
154 std::vector<std::pair<
157 boxes;
158 if (const parallel::TriangulationBase<dim, spacedim> *parallel_tria =
160 &*tria))
161 boxes.reserve(parallel_tria->n_locally_owned_active_cells());
162 else
163 boxes.reserve(tria->n_active_cells());
164 for (const auto &cell : tria->active_cell_iterators() |
165 IteratorFilters::LocallyOwnedCell())
166 boxes.emplace_back(mapping->get_bounding_box(cell), cell);
167
168 locally_owned_cell_bounding_boxes_rtree = pack_rtree(boxes);
169 update_flags =
170 update_flags & ~update_locally_owned_cell_bounding_boxes_rtree;
171 }
172 return locally_owned_cell_bounding_boxes_rtree;
173 }
174
175
177 template <int dim, int spacedim>
178 const RTree<std::pair<BoundingBox<spacedim>, unsigned int>> &
180 {
181 if (update_flags & update_covering_rtree ||
182 covering_rtree.find(level) == covering_rtree.end())
183 {
184 const auto boxes =
185 extract_rtree_level(get_locally_owned_cell_bounding_boxes_rtree(),
186 level);
187
188 if (const auto tria_mpi =
189 dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
190 &(*tria)))
191 {
193 boxes, tria_mpi->get_communicator());
194 }
195 else
196 {
197 covering_rtree[level] =
198 GridTools::build_global_description_tree(boxes, MPI_COMM_SELF);
199 }
200 update_flags = update_flags & ~update_covering_rtree;
201 }
202
203 return covering_rtree[level];
204 }
205
206 template <int dim, int spacedim>
207 const std::vector<std::set<unsigned int>> &
209 {
210 if (update_flags & update_vertex_to_neighbor_subdomain)
211 {
212 vertex_to_neighbor_subdomain.clear();
213 vertex_to_neighbor_subdomain.resize(tria->n_vertices());
214 for (const auto &cell : tria->active_cell_iterators())
215 {
216 if (cell->is_ghost())
217 for (const unsigned int v : cell->vertex_indices())
218 vertex_to_neighbor_subdomain[cell->vertex_index(v)].insert(
219 cell->subdomain_id());
221 update_flags = update_flags & ~update_vertex_to_neighbor_subdomain;
222 }
223 return vertex_to_neighbor_subdomain;
224 }
225
226 template <int dim, int spacedim>
227 const std::map<unsigned int, std::set<types::subdomain_id>> &
229 {
230 if (update_flags & update_vertex_with_ghost_neighbors)
231 {
234
235 update_flags = update_flags & ~update_vertex_with_ghost_neighbors;
236 }
237
239 }
240
241#include "grid_tools_cache.inst"
242
243} // namespace GridTools
244
SmartPointer< const Triangulation< dim, spacedim >, Cache< dim, spacedim > > tria
const std::vector< std::set< typename Triangulation< dim, spacedim >::active_cell_iterator > > & get_vertex_to_cell_map() const
Cache(const Triangulation< dim, spacedim > &tria, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()))
const RTree< std::pair< BoundingBox< spacedim >, unsigned int > > & get_covering_rtree(const unsigned int level=0) const
const std::vector< std::vector< Tensor< 1, spacedim > > > & get_vertex_to_cell_centers_directions() const
const std::vector< std::set< unsigned int > > & get_vertex_to_neighbor_subdomain() const
const std::map< unsigned int, std::set< types::subdomain_id > > & get_vertices_with_ghost_neighbors() const
const RTree< std::pair< BoundingBox< spacedim >, typename Triangulation< dim, spacedim >::active_cell_iterator > > & get_locally_owned_cell_bounding_boxes_rtree() const
const RTree< std::pair< Point< spacedim >, unsigned int > > & get_used_vertices_rtree() const
const RTree< std::pair< BoundingBox< spacedim >, typename Triangulation< dim, spacedim >::active_cell_iterator > > & get_cell_bounding_boxes_rtree() const
const std::map< unsigned int, Point< spacedim > > & get_used_vertices() const
void mark_for_update(const CacheUpdateFlags &flags=update_all)
boost::signals2::connection tria_signal
Abstract base class for mapping classes.
Definition mapping.h:317
unsigned int n_active_cells() const
unsigned int n_vertices() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 3 > vertices[4]
Point< 2 > second
Definition grid_out.cc:4616
Point< 2 > first
Definition grid_out.cc:4615
unsigned int level
Definition grid_out.cc:4618
IteratorRange< active_cell_iterator > active_cell_iterators() const
RTree< std::pair< BoundingBox< spacedim >, unsigned int > > build_global_description_tree(const std::vector< BoundingBox< spacedim > > &local_description, const MPI_Comm mpi_communicator)
std::map< unsigned int, Point< spacedim > > extract_used_vertices(const Triangulation< dim, spacedim > &container, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()))
std::map< unsigned int, std::set<::types::subdomain_id > > compute_vertices_with_ghost_neighbors(const Triangulation< dim, spacedim > &tria)
std::vector< std::set< typename Triangulation< dim, spacedim >::active_cell_iterator > > vertex_to_cell_map(const Triangulation< dim, spacedim > &triangulation)
std::vector< std::vector< Tensor< 1, spacedim > > > vertex_to_cell_centers_directions(const Triangulation< dim, spacedim > &mesh, const std::vector< std::set< typename Triangulation< dim, spacedim >::active_cell_iterator > > &vertex_to_cells)
STL namespace.
std::map< unsigned int, std::set<::types::subdomain_id > > * vertices_with_ghost_neighbors
std::vector< BoundingBox< boost::geometry::dimension< typename Rtree::indexable_type >::value > > extract_rtree_level(const Rtree &tree, const unsigned int level)
boost::geometry::index::rtree< LeafType, IndexType, IndexableGetter > RTree
Definition rtree.h:144
RTree< typename LeafTypeIterator::value_type, IndexType, IndexableGetter > pack_rtree(const LeafTypeIterator &begin, const LeafTypeIterator &end)
const ::Triangulation< dim, spacedim > & tria