]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Make access to data members of GridTools::Cache thread-safe.
authorWolfgang Bangerth <bangerth@colostate.edu>
Wed, 19 Jun 2024 23:34:21 +0000 (17:34 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 19 Jun 2024 23:34:39 +0000 (17:34 -0600)
include/deal.II/grid/grid_tools_cache.h
source/grid/grid_tools_cache.cc

index ee158a51f12548b42111334e91106cfa1f5afa79..d5843358028183e38f64e260b39452f4dc922e2a 100644 (file)
@@ -220,7 +220,13 @@ namespace GridTools
 
   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;
 
@@ -241,14 +247,16 @@ namespace GridTools
      */
     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.
@@ -260,18 +268,21 @@ namespace GridTools
      */
     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
@@ -280,7 +291,8 @@ namespace GridTools
     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
@@ -289,21 +301,23 @@ namespace GridTools
     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.
index 04b57a1a89dd5629b8678e05fd5d100b4e68c215..a320b82f39f1e7a6ec1464d8f175375486860514 100644 (file)
@@ -37,7 +37,7 @@ namespace GridTools
   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();
@@ -59,6 +59,13 @@ namespace GridTools
     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);
@@ -73,6 +80,13 @@ namespace GridTools
   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(
@@ -88,6 +102,13 @@ namespace GridTools
   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);
@@ -102,6 +123,13 @@ namespace GridTools
   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();
@@ -124,6 +152,13 @@ namespace GridTools
               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<
@@ -148,6 +183,14 @@ namespace GridTools
               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<
@@ -177,6 +220,13 @@ namespace GridTools
   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())
       {
@@ -202,10 +252,19 @@ namespace GridTools
     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();
@@ -222,10 +281,19 @@ namespace GridTools
     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 =

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.