private:
/**
- * Keep track of what needs to be updated next.
+ * Keep track of what needs to be updated every time the triangulation
+ * is changed. Each of the get_*() functions above checks whether a
+ * data structure still needs to be updated before returning a reference
+ * to it; if it needs to be updated, it does so and resets the
+ * corresponding flag. The constructor as well as the signal slot
+ * that is called whenever the triangulation changes the re-sets this
+ * flag to `update_all` again.
*/
mutable CacheUpdateFlags update_flags;
*/
mutable std::vector<
std::set<typename Triangulation<dim, spacedim>::active_cell_iterator>>
- vertex_to_cells;
+ vertex_to_cells;
+ mutable std::mutex vertex_to_cells_mutex;
/**
* Store vertex to cell center directions, as generated by
* GridTools::vertex_to_cell_centers_directions().
*/
mutable std::vector<std::vector<Tensor<1, spacedim>>>
- vertex_to_cell_centers;
+ vertex_to_cell_centers;
+ mutable std::mutex vertex_to_cell_centers_mutex;
/**
* A collection of rtree objects covering the whole mesh.
*/
mutable std::map<unsigned int,
RTree<std::pair<BoundingBox<spacedim>, unsigned int>>>
- covering_rtree;
+ covering_rtree;
+ mutable std::mutex covering_rtree_mutex;
/**
* Store the used vertices of the Triangulation, as generated by
* GridTools::extract_used_vertices().
*/
mutable std::map<unsigned int, Point<spacedim>> used_vertices;
+ mutable std::mutex used_vertices_mutex;
/**
* Store an RTree object, containing the used vertices of the triangulation.
*/
mutable RTree<std::pair<Point<spacedim>, unsigned int>> used_vertices_rtree;
+ mutable std::mutex used_vertices_rtree_mutex;
/**
* Store an RTree object, containing the bounding boxes of the cells of the
mutable RTree<
std::pair<BoundingBox<spacedim>,
typename Triangulation<dim, spacedim>::active_cell_iterator>>
- cell_bounding_boxes_rtree;
+ cell_bounding_boxes_rtree;
+ mutable std::mutex cell_bounding_boxes_rtree_mutex;
/**
* Store an RTree object, containing the bounding boxes of the locally owned
mutable RTree<
std::pair<BoundingBox<spacedim>,
typename Triangulation<dim, spacedim>::active_cell_iterator>>
- locally_owned_cell_bounding_boxes_rtree;
-
+ locally_owned_cell_bounding_boxes_rtree;
+ mutable std::mutex locally_owned_cell_bounding_boxes_rtree_mutex;
/**
* Store an std::vector of std::set of integer containing the id of all
* subdomain to which a vertex is connected to.
*/
mutable std::vector<std::set<unsigned int>> vertex_to_neighbor_subdomain;
+ mutable std::mutex vertex_to_neighbor_subdomain_mutex;
/**
* Store an std::map of unsigned integer containing
* the set of subdomains connected to each vertex.
*/
mutable std::map<unsigned int, std::set<dealii::types::subdomain_id>>
- vertices_with_ghost_neighbors;
+ vertices_with_ghost_neighbors;
+ mutable std::mutex vertices_with_ghost_neighbors_mutex;
/**
* Storage for the status of the triangulation signal.
template <int dim, int spacedim>
Cache<dim, spacedim>::~Cache()
{
- // Make sure that the signals that was attached to the triangulation
+ // Make sure that the signal that was attached to the triangulation
// is removed here.
if (tria_signal.connected())
tria_signal.disconnect();
std::set<typename Triangulation<dim, spacedim>::active_cell_iterator>> &
Cache<dim, spacedim>::get_vertex_to_cell_map() const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(vertex_to_cells_mutex);
+
if (update_flags & update_vertex_to_cell_map)
{
vertex_to_cells = GridTools::vertex_to_cell_map(*tria);
const std::vector<std::vector<Tensor<1, spacedim>>> &
Cache<dim, spacedim>::get_vertex_to_cell_centers_directions() const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(vertex_to_cell_centers_mutex);
+
if (update_flags & update_vertex_to_cell_centers_directions)
{
vertex_to_cell_centers = GridTools::vertex_to_cell_centers_directions(
const std::map<unsigned int, Point<spacedim>> &
Cache<dim, spacedim>::get_used_vertices() const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(used_vertices_mutex);
+
if (update_flags & update_used_vertices)
{
used_vertices = GridTools::extract_used_vertices(*tria, *mapping);
const RTree<std::pair<Point<spacedim>, unsigned int>> &
Cache<dim, spacedim>::get_used_vertices_rtree() const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(used_vertices_rtree_mutex);
+
if (update_flags & update_used_vertices_rtree)
{
const auto &used_vertices = get_used_vertices();
typename Triangulation<dim, spacedim>::active_cell_iterator>> &
Cache<dim, spacedim>::get_cell_bounding_boxes_rtree() const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(cell_bounding_boxes_rtree_mutex);
+
if (update_flags & update_cell_bounding_boxes_rtree)
{
std::vector<std::pair<
typename Triangulation<dim, spacedim>::active_cell_iterator>> &
Cache<dim, spacedim>::get_locally_owned_cell_bounding_boxes_rtree() const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(
+ locally_owned_cell_bounding_boxes_rtree_mutex);
+
if (update_flags & update_locally_owned_cell_bounding_boxes_rtree)
{
std::vector<std::pair<
const RTree<std::pair<BoundingBox<spacedim>, unsigned int>> &
Cache<dim, spacedim>::get_covering_rtree(const unsigned int level) const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(covering_rtree_mutex);
+
if (update_flags & update_covering_rtree ||
covering_rtree.find(level) == covering_rtree.end())
{
return covering_rtree[level];
}
+
+
template <int dim, int spacedim>
const std::vector<std::set<unsigned int>> &
Cache<dim, spacedim>::get_vertex_to_neighbor_subdomain() const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(vertex_to_neighbor_subdomain_mutex);
+
if (update_flags & update_vertex_to_neighbor_subdomain)
{
vertex_to_neighbor_subdomain.clear();
return vertex_to_neighbor_subdomain;
}
+
+
template <int dim, int spacedim>
const std::map<unsigned int, std::set<types::subdomain_id>> &
Cache<dim, spacedim>::get_vertices_with_ghost_neighbors() const
{
+ // In the following, we will first check whether the data structure
+ // in question needs to be updated (in which case we update it, and
+ // reset the flag that indices that this needs to happen to zero), and
+ // then return it. Make this thread-safe by using a mutex to guard
+ // all of this:
+ std::lock_guard<std::mutex> lock(vertices_with_ghost_neighbors_mutex);
+
if (update_flags & update_vertex_with_ghost_neighbors)
{
vertices_with_ghost_neighbors =