]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Revise some operations in compute_index_owner_internal 14286/head
authorMartin Kronbichler <martin.kronbichler@uni-a.de>
Mon, 19 Sep 2022 10:32:27 +0000 (12:32 +0200)
committerMartin Kronbichler <martin.kronbichler@uni-a.de>
Mon, 19 Sep 2022 20:23:56 +0000 (22:23 +0200)
include/deal.II/base/mpi_compute_index_owner_internal.h
source/base/mpi_compute_index_owner_internal.cc

index 1ebed7543b85b5e9b1e5baac8a215b4e48b15036..ecba6c38b463b686b989e25c9622d87dc96f6560 100644 (file)
@@ -291,18 +291,15 @@ namespace Utilities
           Dictionary dict;
 
           /**
-           * Array to collect the indices to look up, sorted by the rank in
+           * Array to collect the indices to look up (first vector) and their
+           * local index among indices (second vector), sorted by the rank in
            * the dictionary.
            */
-          std::map<unsigned int, std::vector<types::global_dof_index>>
+          std::map<unsigned int,
+                   std::pair<std::vector<types::global_dof_index>,
+                             std::vector<unsigned int>>>
             indices_to_look_up_by_dict_rank;
 
-          /**
-           * The field where the indices for incoming data from the process
-           * are stored.
-           */
-          std::map<unsigned int, std::vector<unsigned int>> recv_indices;
-
           /**
            * Implementation of
            * Utilities::MPI::ConsensusAlgorithms::Process::answer_request(),
@@ -356,12 +353,12 @@ namespace Utilities
 
         private:
           /**
-           * Stores the index request in the `requesters` field. We first find
-           * out the owner of the index that was requested (using the guess in
-           * `owner_index`, as we typically might look up on the same rank
-           * several times in a row, which avoids the binary search in
-           * Dictionary::get_owning_rank_index()). Once we know the rank of
-           * the owner, we fill the vector entry with the rank of the
+           * Stores the index request in the `requesters` field. Given the
+           * rank of the owner, we start with a guess for the index at the
+           * owner's site. This is because we typically might look up on the
+           * same rank several times in a row, hence avoiding the binary
+           * search in Dictionary::get_owning_rank_index()). Once we know the
+           * index at the owner, we fill the vector entry with the rank of the
            * request. Here, we utilize the fact that requests are processed
            * rank-by-rank, so we can simply look at the end of the vector
            * whether there is already some data stored or not. Finally, we
@@ -369,9 +366,10 @@ namespace Utilities
            * therefore only need to append at the end.
            */
           void
-          append_index_origin(const types::global_dof_index index,
-                              unsigned int &                owner_index,
-                              const unsigned int            rank_of_request);
+          append_index_origin(const unsigned int index_within_dictionary,
+                              const unsigned int rank_of_request,
+                              const unsigned int rank_of_owner,
+                              unsigned int &     owner_index_guess);
         };
 
         /* ------------------------- inline functions ----------------------- */
index c96e46a9f05e62353ec02715ffb7a664d3fe6898..3beb029a2d609e38a35e4f663c62ae79f2475925 100644 (file)
@@ -486,7 +486,7 @@ namespace Utilities
                                       types::global_dof_index>> &buffer_recv,
           std::vector<unsigned int> &                            request_buffer)
         {
-          unsigned int owner_index = 0;
+          unsigned int owner_index_guess = 0;
           for (const auto &interval : buffer_recv)
             for (auto i = interval.first; i < interval.second; ++i)
               {
@@ -495,7 +495,10 @@ namespace Utilities
                 request_buffer.push_back(actual_owner);
 
                 if (track_index_requests)
-                  append_index_origin(i, owner_index, other_rank);
+                  append_index_origin(i - dict.local_range.first,
+                                      other_rank,
+                                      actual_owner,
+                                      owner_index_guess);
               }
         }
 
@@ -506,50 +509,34 @@ namespace Utilities
         {
           std::vector<unsigned int> targets;
 
-          // 1) collect relevant processes and process local dict entries
-          {
-            unsigned int index       = 0;
-            unsigned int owner_index = 0;
-            for (auto i : indices_to_look_up)
-              {
-                unsigned int other_rank = dict.dof_to_dict_rank(i);
-                if (other_rank == my_rank)
-                  {
-                    owning_ranks[index] =
-                      dict.actually_owning_ranks[i - dict.local_range.first];
-                    if (track_index_requests)
-                      append_index_origin(i, owner_index, my_rank);
-                  }
-                else if (targets.empty() || targets.back() != other_rank)
-                  targets.push_back(other_rank);
-                index++;
-              }
-          }
-
-
-          for (auto i : targets)
+          indices_to_look_up_by_dict_rank.clear();
+          unsigned int index             = 0;
+          unsigned int owner_index_guess = 0;
+          for (auto i : indices_to_look_up)
             {
-              recv_indices[i]                    = {};
-              indices_to_look_up_by_dict_rank[i] = {};
+              unsigned int other_rank = dict.dof_to_dict_rank(i);
+              if (other_rank == my_rank)
+                {
+                  owning_ranks[index] =
+                    dict.actually_owning_ranks[i - dict.local_range.first];
+                  if (track_index_requests)
+                    append_index_origin(i - dict.local_range.first,
+                                        my_rank,
+                                        owning_ranks[index],
+                                        owner_index_guess);
+                }
+              else
+                {
+                  if (targets.empty() || targets.back() != other_rank)
+                    targets.push_back(other_rank);
+                  auto &indices = indices_to_look_up_by_dict_rank[other_rank];
+                  indices.first.push_back(i);
+                  indices.second.push_back(index);
+                }
+              index++;
             }
 
-          // 3) collect indices for each process
-          {
-            unsigned int index = 0;
-            for (auto i : indices_to_look_up)
-              {
-                unsigned int other_rank = dict.dof_to_dict_rank(i);
-                if (other_rank != my_rank)
-                  {
-                    recv_indices[other_rank].push_back(index);
-                    indices_to_look_up_by_dict_rank[other_rank].push_back(i);
-                  }
-                index++;
-              }
-          }
-
-          Assert(targets.size() == recv_indices.size() &&
-                   targets.size() == indices_to_look_up_by_dict_rank.size(),
+          Assert(targets.size() == indices_to_look_up_by_dict_rank.size(),
                  ExcMessage("Size does not match!"));
 
           return targets;
@@ -564,7 +551,7 @@ namespace Utilities
                                 types::global_dof_index>> &send_buffer)
         {
           // create index set and compress data to be sent
-          auto &   indices_i = indices_to_look_up_by_dict_rank[other_rank];
+          auto &indices_i = indices_to_look_up_by_dict_rank[other_rank].first;
           IndexSet is(dict.size);
           is.add_indices(indices_i.begin(), indices_i.end());
           is.compress();
@@ -582,13 +569,11 @@ namespace Utilities
           const unsigned int               other_rank,
           const std::vector<unsigned int> &recv_buffer)
         {
-          Assert(recv_buffer.size() == recv_indices[other_rank].size(),
-                 ExcInternalError());
-          Assert(recv_indices[other_rank].size() == recv_buffer.size(),
-                 ExcMessage("Sizes do not match!"));
-
-          for (unsigned int j = 0; j < recv_indices[other_rank].size(); ++j)
-            owning_ranks[recv_indices[other_rank][j]] = recv_buffer[j];
+          const auto &recv_indices =
+            indices_to_look_up_by_dict_rank[other_rank].second;
+          AssertDimension(recv_indices.size(), recv_buffer.size());
+          for (unsigned int j = 0; j < recv_indices.size(); ++j)
+            owning_ranks[recv_indices[j]] = recv_buffer[j];
         }
 
 
@@ -752,30 +737,29 @@ namespace Utilities
 
         void
         ConsensusAlgorithmsPayload::append_index_origin(
-          const types::global_dof_index index,
-          unsigned int &                owner_index,
-          const unsigned int            rank_of_request)
+          const unsigned int index_within_dict,
+          const unsigned int rank_of_request,
+          const unsigned int rank_of_owner,
+          unsigned int &     owner_index_guess)
         {
           // remember who requested which index. We want to use an
           // std::vector with simple addressing, via a good guess from the
           // preceding index, rather than std::map, because this is an inner
           // loop and it avoids the map lookup in every iteration
-          const unsigned int rank_of_owner =
-            dict.actually_owning_ranks[index - dict.local_range.first];
-          owner_index = dict.get_owning_rank_index(rank_of_owner, owner_index);
-          if (requesters[owner_index].empty() ||
-              requesters[owner_index].back().first != rank_of_request)
-            requesters[owner_index].emplace_back(
+          owner_index_guess =
+            dict.get_owning_rank_index(rank_of_owner, owner_index_guess);
+
+          auto &request = requesters[owner_index_guess];
+          if (request.empty() || request.back().first != rank_of_request)
+            request.emplace_back(
               rank_of_request,
               std::vector<std::pair<unsigned int, unsigned int>>());
-          if (requesters[owner_index].back().second.empty() ||
-              requesters[owner_index].back().second.back().second !=
-                index - dict.local_range.first)
-            requesters[owner_index].back().second.emplace_back(
-              index - dict.local_range.first,
-              index - dict.local_range.first + 1);
+
+          auto &intervals = request.back().second;
+          if (intervals.empty() || intervals.back().second != index_within_dict)
+            intervals.emplace_back(index_within_dict, index_within_dict + 1);
           else
-            ++requesters[owner_index].back().second.back().second;
+            ++intervals.back().second;
         }
       } // namespace ComputeIndexOwner
     }   // namespace internal

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.