Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3075-gc235bd4825 2025-04-15 08:10:00+00:00
\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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
grid_tools_cache.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2017 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_grid_grid_tools_cache_h
16#define dealii_grid_grid_tools_cache_h
17
18
19#include <deal.II/base/config.h>
20
23#include <deal.II/base/point.h>
24
25#include <deal.II/fe/mapping.h>
26
28#include <deal.II/grid/tria.h>
31
33
34#include <boost/signals2.hpp>
35
36#include <atomic>
37#include <cmath>
38#include <set>
39
40
42
43namespace GridTools
44{
66 template <int dim, int spacedim = dim>
68 {
69 public:
81
89
93 ~Cache() override;
94
105 void
107
108
113 const std::vector<
114 std::set<typename Triangulation<dim, spacedim>::active_cell_iterator>> &
116
121 const std::vector<std::vector<Tensor<1, spacedim>>> &
123
128 const std::map<unsigned int, Point<spacedim>> &
129 get_used_vertices() const;
130
135 const RTree<std::pair<Point<spacedim>, unsigned int>> &
137
144 const RTree<
145 std::pair<BoundingBox<spacedim>,
148
159 const RTree<
160 std::pair<BoundingBox<spacedim>,
163
164
171 const std::vector<std::set<unsigned int>> &
173
178 const std::map<unsigned int, std::set<types::subdomain_id>> &
180
186
191 get_mapping() const;
192
193
222 const RTree<std::pair<BoundingBox<spacedim>, unsigned int>> &
223 get_covering_rtree(const unsigned int level = 0) const;
224
225 private:
240 mutable std::atomic<std::underlying_type_t<CacheUpdateFlags>> update_flags;
241
247
252
253
258 mutable std::vector<
259 std::set<typename Triangulation<dim, spacedim>::active_cell_iterator>>
261 mutable std::mutex vertex_to_cells_mutex;
262
267 mutable std::vector<std::vector<Tensor<1, spacedim>>>
270
279 mutable std::map<unsigned int,
282 mutable std::mutex covering_rtree_mutex;
283
288 mutable std::map<unsigned int, Point<spacedim>> used_vertices;
289 mutable std::mutex used_vertices_mutex;
290
295 mutable std::mutex used_vertices_rtree_mutex;
296
301 mutable RTree<
302 std::pair<BoundingBox<spacedim>,
306
311 mutable RTree<
312 std::pair<BoundingBox<spacedim>,
316
321 mutable std::vector<std::set<unsigned int>> vertex_to_neighbor_subdomain;
323
328 mutable std::map<unsigned int, std::set<::types::subdomain_id>>
331
335 boost::signals2::connection tria_change_signal;
336
340 boost::signals2::connection tria_create_signal;
341 };
342
343
344
345 // Inline functions
346 template <int dim, int spacedim>
347 inline const Triangulation<dim, spacedim> &
349 {
350 return *tria;
351 }
352
353
354
355 template <int dim, int spacedim>
356 inline const Mapping<dim, spacedim> &
358 {
359 Assert(mapping, ExcNotInitialized());
360 return *mapping;
361 }
362} // namespace GridTools
363
364
365
367
368#endif
const std::vector< std::set< typename Triangulation< dim, spacedim >::active_cell_iterator > > & get_vertex_to_cell_map() const
std::vector< std::vector< Tensor< 1, spacedim > > > vertex_to_cell_centers
std::mutex cell_bounding_boxes_rtree_mutex
std::mutex covering_rtree_mutex
const RTree< std::pair< BoundingBox< spacedim >, unsigned int > > & get_covering_rtree(const unsigned int level=0) const
ObserverPointer< const Mapping< dim, spacedim >, Cache< dim, spacedim > > mapping
std::vector< std::set< unsigned int > > vertex_to_neighbor_subdomain
const std::vector< std::vector< Tensor< 1, spacedim > > > & get_vertex_to_cell_centers_directions() const
boost::signals2::connection tria_change_signal
boost::signals2::connection tria_create_signal
RTree< std::pair< BoundingBox< spacedim >, typename Triangulation< dim, spacedim >::active_cell_iterator > > cell_bounding_boxes_rtree
std::mutex vertex_to_cells_mutex
RTree< std::pair< BoundingBox< spacedim >, typename Triangulation< dim, spacedim >::active_cell_iterator > > locally_owned_cell_bounding_boxes_rtree
const std::vector< std::set< unsigned int > > & get_vertex_to_neighbor_subdomain() const
std::mutex vertex_to_neighbor_subdomain_mutex
std::map< unsigned int, Point< spacedim > > used_vertices
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
std::mutex vertex_to_cell_centers_mutex
const Mapping< dim, spacedim > & get_mapping() const
RTree< std::pair< Point< spacedim >, unsigned int > > used_vertices_rtree
const Triangulation< dim, spacedim > & get_triangulation() const
std::mutex used_vertices_rtree_mutex
std::atomic< std::underlying_type_t< CacheUpdateFlags > > update_flags
const RTree< std::pair< Point< spacedim >, unsigned int > > & get_used_vertices_rtree() const
std::map< unsigned int, std::set<::types::subdomain_id > > vertices_with_ghost_neighbors
const RTree< std::pair< BoundingBox< spacedim >, typename Triangulation< dim, spacedim >::active_cell_iterator > > & get_cell_bounding_boxes_rtree() const
std::vector< std::set< typename Triangulation< dim, spacedim >::active_cell_iterator > > vertex_to_cells
std::mutex locally_owned_cell_bounding_boxes_rtree_mutex
const std::map< unsigned int, Point< spacedim > > & get_used_vertices() const
std::mutex vertices_with_ghost_neighbors_mutex
std::map< unsigned int, RTree< std::pair< BoundingBox< spacedim >, unsigned int > > > covering_rtree
ObserverPointer< const Triangulation< dim, spacedim >, Cache< dim, spacedim > > tria
void mark_for_update(const CacheUpdateFlags &flags=update_all)
std::mutex used_vertices_mutex
Abstract base class for mapping classes.
Definition mapping.h:320
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
unsigned int level
Definition grid_out.cc:4635
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotInitialized()
boost::geometry::index::rtree< LeafType, IndexType, IndexableGetter > RTree
Definition rtree.h:145