]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use the new CA interface instead of using AnonymousProcess. 13326/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Tue, 1 Feb 2022 19:01:21 +0000 (12:01 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 2 Feb 2022 19:12:24 +0000 (12:12 -0700)
include/deal.II/fe/fe_tools_extrapolate.templates.h
include/deal.II/matrix_free/matrix_free.templates.h
source/base/mpi.cc
source/grid/grid_tools.cc
source/grid/tria_description.cc

index aa58ae54114fe22d5587634e5f3046932c78ad23..8eb09d692831312aa73535b901de8091d1ab72f2 100644 (file)
@@ -1139,30 +1139,26 @@ namespace FETools
       // communication of information in cases where we do not know up
       // front which processes (and from how many processes) we have to
       // expect information from.
-      const auto get_destinations = [&destinations]() { return destinations; };
-
       const auto create_request =
-        [&cells_to_send](const types::subdomain_id other_rank,
-                         std::vector<char> &       send_buffer) {
+        [&cells_to_send](const types::subdomain_id other_rank) {
           std::vector<CellData> cells_for_this_destination;
           for (const auto &cell : cells_to_send)
             if (cell.receiver == other_rank)
               cells_for_this_destination.emplace_back(cell);
 
-          send_buffer = Utilities::pack(cells_for_this_destination, false);
+          return Utilities::pack(cells_for_this_destination, false);
         };
 
       const auto answer_request =
         [&received_cells](const unsigned int       other_rank,
-                          const std::vector<char> &buffer_recv,
-                          std::vector<char> &      request_buffer) {
+                          const std::vector<char> &request) {
           // We got a message from 'other_rank', so let us decode the
           // message in the same way as we have assembled it above.
           // Note that the cells just received do not contain
           // information where they came from, and we have to add that
           // ourselves for later use.
           for (CellData &cell_data :
-               Utilities::unpack<std::vector<CellData>>(buffer_recv, false))
+               Utilities::unpack<std::vector<CellData>>(request, false))
             {
               cell_data.receiver = other_rank;
               received_cells.emplace_back(std::move(cell_data));
@@ -1170,26 +1166,23 @@ namespace FETools
 
           // Nothing left to do here, we don't actually need to provide an
           // answer:
-          request_buffer.clear();
+          return std::vector<char>();
         };
 
       const auto read_answer = [](const unsigned int /*other_rank*/,
-                                  const std::vector<char> &recv_buffer) {
+                                  const std::vector<char> &answer) {
         // We don't put anything into the answers, so nothing should
         // have been coming out at this end either:
-        (void)recv_buffer;
-        Assert(recv_buffer.size() == 0, ExcInternalError());
+        (void)answer;
+        Assert(answer.size() == 0, ExcInternalError());
       };
 
-      Utilities::MPI::ConsensusAlgorithms::AnonymousProcess<char, char>
-        operations(get_destinations,
-                   create_request,
-                   answer_request,
-                   read_answer);
-
-      Utilities::MPI::ConsensusAlgorithms::Selector<char, char>(operations,
-                                                                communicator)
-        .run();
+      Utilities::MPI::ConsensusAlgorithms::Selector<char, char>().run(
+        destinations,
+        create_request,
+        answer_request,
+        read_answer,
+        communicator);
     }
 
 
index 253f141f7f1bcdda2419322452e43959d1d54c3c..7a341ee05e227edef5ff3d1ec9cc75c1b7382abd 100644 (file)
@@ -2094,65 +2094,66 @@ MatrixFree<dim, Number, VectorizedArrayType>::initialize_indices(
                     });
 
           // get offsets of shared cells
-          Utilities::MPI::ConsensusAlgorithms::
-            AnonymousProcess<dealii::types::global_dof_index, unsigned int>
-              process(
-                [&]() {
-                  std::vector<unsigned int> targets;
-                  for (unsigned int i = 0; i < cells_shared_ghosts.size(); ++i)
-                    if (cells_shared_ghosts[i].size() > 0)
-                      targets.push_back(i);
-                  return targets;
-                },
-                [&](const auto other_rank, auto &send_buffer) {
-                  auto &source = cells_shared_ghosts[other_rank];
-                  std::sort(source.begin(),
-                            source.end(),
-                            [](const auto &a, const auto &b) {
-                              return a.second < b.second;
-                            });
-
-                  send_buffer.reserve(source.size());
-                  for (const auto &i : source)
-                    send_buffer.push_back(i.second);
-                },
-                [&](const auto &other_rank,
-                    const auto &buffer_recv,
-                    auto &      request_buffer) {
-                  (void)other_rank;
-
-                  request_buffer.reserve(buffer_recv.size());
-
-                  unsigned int j = 0;
-
-                  for (unsigned int i = 0; i < buffer_recv.size(); ++i)
-                    {
-                      for (; j < cells_locally_owned.size(); ++j)
-                        if (cells_locally_owned[j].second == buffer_recv[i])
-                          {
-                            request_buffer.push_back(
-                              cells_locally_owned[j].first);
-                            break;
-                          }
-                    }
+          std::vector<unsigned int> targets;
+          for (unsigned int i = 0; i < cells_shared_ghosts.size(); ++i)
+            if (cells_shared_ghosts[i].size() > 0)
+              targets.push_back(i);
+
+          // Then set up the callbacks the consensus algorithm needs:
+          const auto create_request = [&](const auto other_rank) {
+            auto &source = cells_shared_ghosts[other_rank];
+            std::sort(source.begin(),
+                      source.end(),
+                      [](const auto &a, const auto &b) {
+                        return a.second < b.second;
+                      });
+
+            std::vector<dealii::types::global_dof_index> send_buffer;
+            send_buffer.reserve(source.size());
+            for (const auto &i : source)
+              send_buffer.push_back(i.second);
+
+            return send_buffer;
+          };
+
+          const auto answer_request = [&](const auto, const auto &request) {
+            std::vector<unsigned int> answer;
+            answer.reserve(request.size());
 
-                  AssertDimension(request_buffer.size(), buffer_recv.size());
-                },
-                [&](const auto other_rank, const auto &recv_buffer) {
-                  Assert(recv_buffer.size() ==
-                           cells_shared_ghosts[other_rank].size(),
-                         ExcInternalError());
-                  for (unsigned int i = 0; i < recv_buffer.size(); ++i)
+            unsigned int j = 0;
+
+            for (unsigned int i = 0; i < request.size(); ++i)
+              {
+                for (; j < cells_locally_owned.size(); ++j)
+                  if (cells_locally_owned[j].second == request[i])
                     {
-                      cells[cells_shared_ghosts[other_rank][i].first] = {
-                        other_rank, recv_buffer[i]};
+                      answer.push_back(cells_locally_owned[j].first);
+                      break;
                     }
-                });
+              }
+
+            AssertDimension(answer.size(), request.size());
+            return answer;
+          };
 
-          Utilities::MPI::ConsensusAlgorithms::Selector<
-            dealii::types::global_dof_index,
-            unsigned int>(process, communicator_sm)
-            .run();
+          const auto process_answer = [&](const auto  other_rank,
+                                          const auto &answer) {
+            Assert(answer.size() == cells_shared_ghosts[other_rank].size(),
+                   ExcInternalError());
+            for (unsigned int i = 0; i < answer.size(); ++i)
+              {
+                cells[cells_shared_ghosts[other_rank][i].first] = {other_rank,
+                                                                   answer[i]};
+              }
+          };
+
+          Utilities::MPI::ConsensusAlgorithms::
+            Selector<dealii::types::global_dof_index, unsigned int>()
+              .run(targets,
+                   create_request,
+                   answer_request,
+                   process_answer,
+                   communicator_sm);
 
           return cells;
         }();
index ed78d4f6e227d4fddf853fe04fd699f3a14fc173..481102eeb642e68513912c5ee7d2852d332f6564 100644 (file)
@@ -432,11 +432,8 @@ namespace Utilities
       if (Utilities::MPI::min((my_destinations_are_unique ? 1 : 0), mpi_comm) ==
           1)
         {
-          ConsensusAlgorithms::AnonymousProcess<char, char> process(
-            [&]() { return destinations; });
-          ConsensusAlgorithms::NBX<char, char> consensus_algorithm(process,
-                                                                   mpi_comm);
-          return consensus_algorithm.run();
+          return ConsensusAlgorithms::NBX<char, char>().run(
+            destinations, {}, {}, {}, mpi_comm);
         }
         // If that was not the case, we need to use the remainder of the code
         // below, i.e., just fall through the if condition above.
@@ -591,11 +588,9 @@ namespace Utilities
       if (Utilities::MPI::min((my_destinations_are_unique ? 1 : 0), mpi_comm) ==
           1)
         {
-          ConsensusAlgorithms::AnonymousProcess<char, char> process(
-            [&]() { return destinations; });
-          ConsensusAlgorithms::NBX<char, char> consensus_algorithm(process,
-                                                                   mpi_comm);
-          return consensus_algorithm.run().size();
+          return ConsensusAlgorithms::NBX<char, char>()
+            .run(destinations, {}, {}, {}, mpi_comm)
+            .size();
         }
       else
 #  endif
index 620ae281d1ccdc2b268a68e9e70fdb4291dc2388..9cd09727d8d16cee168ca8b868789461a0226806 100644 (file)
@@ -5933,94 +5933,98 @@ namespace GridTools
         (std::find(marked_vertices.begin(), marked_vertices.end(), true) !=
          marked_vertices.end());
 
-      Utilities::MPI::ConsensusAlgorithms::AnonymousProcess<char, char> process(
-        [&]() { return potential_owners_ranks; },
-        [&](const unsigned int other_rank, std::vector<char> &send_buffer) {
-          const auto other_rank_index = translate(other_rank);
+      const auto create_request = [&](const unsigned int other_rank) {
+        const auto other_rank_index = translate(other_rank);
 
-          std::vector<std::pair<unsigned int, Point<spacedim>>> temp;
-          temp.reserve(potential_owners_ptrs[other_rank_index + 1] -
-                       potential_owners_ptrs[other_rank_index]);
+        std::vector<std::pair<unsigned int, Point<spacedim>>> temp;
+        temp.reserve(potential_owners_ptrs[other_rank_index + 1] -
+                     potential_owners_ptrs[other_rank_index]);
 
-          for (unsigned int i = potential_owners_ptrs[other_rank_index];
-               i < potential_owners_ptrs[other_rank_index + 1];
-               ++i)
-            temp.emplace_back(potential_owners_indices[i],
-                              points[potential_owners_indices[i]]);
+        for (unsigned int i = potential_owners_ptrs[other_rank_index];
+             i < potential_owners_ptrs[other_rank_index + 1];
+             ++i)
+          temp.emplace_back(potential_owners_indices[i],
+                            points[potential_owners_indices[i]]);
 
-          send_buffer = Utilities::pack(temp, false);
-        },
-        [&](const unsigned int &     other_rank,
-            const std::vector<char> &recv_buffer,
-            std::vector<char> &      request_buffer) {
-          const auto recv_buffer_unpacked = Utilities::unpack<
-            std::vector<std::pair<unsigned int, Point<spacedim>>>>(recv_buffer,
-                                                                   false);
+        return Utilities::pack(temp, false);
+      };
 
-          std::vector<unsigned int> request_buffer_temp(
-            recv_buffer_unpacked.size(), 0);
+      const auto answer_request =
+        [&](const unsigned int &     other_rank,
+            const std::vector<char> &request) -> std::vector<char> {
+        const auto recv_buffer_unpacked = Utilities::unpack<
+          std::vector<std::pair<unsigned int, Point<spacedim>>>>(request,
+                                                                 false);
 
-          if (has_relevant_vertices)
-            {
-              cell_hint = cache.get_triangulation().begin_active();
+        std::vector<unsigned int> request_buffer_temp(
+          recv_buffer_unpacked.size(), 0);
 
-              for (unsigned int i = 0; i < recv_buffer_unpacked.size(); ++i)
-                {
-                  const auto &index_and_point = recv_buffer_unpacked[i];
+        if (has_relevant_vertices)
+          {
+            cell_hint = cache.get_triangulation().begin_active();
 
-                  const auto cells_and_reference_positions =
-                    find_all_locally_owned_active_cells_around_point(
-                      cache,
+            for (unsigned int i = 0; i < recv_buffer_unpacked.size(); ++i)
+              {
+                const auto &index_and_point = recv_buffer_unpacked[i];
+
+                const auto cells_and_reference_positions =
+                  find_all_locally_owned_active_cells_around_point(
+                    cache,
+                    index_and_point.second,
+                    cell_hint,
+                    marked_vertices,
+                    tolerance,
+                    enforce_unique_mapping);
+
+                for (const auto &cell_and_reference_position :
+                     cells_and_reference_positions)
+                  {
+                    send_components.emplace_back(
+                      std::pair<int, int>(
+                        cell_and_reference_position.first->level(),
+                        cell_and_reference_position.first->index()),
+                      other_rank,
+                      index_and_point.first,
+                      cell_and_reference_position.second,
                       index_and_point.second,
-                      cell_hint,
-                      marked_vertices,
-                      tolerance,
-                      enforce_unique_mapping);
+                      numbers::invalid_unsigned_int);
+                  }
 
-                  for (const auto &cell_and_reference_position :
-                       cells_and_reference_positions)
-                    {
-                      send_components.emplace_back(
-                        std::pair<int, int>(
-                          cell_and_reference_position.first->level(),
-                          cell_and_reference_position.first->index()),
-                        other_rank,
-                        index_and_point.first,
-                        cell_and_reference_position.second,
-                        index_and_point.second,
-                        numbers::invalid_unsigned_int);
-                    }
+                request_buffer_temp[i] = cells_and_reference_positions.size();
+              }
+          }
 
-                  request_buffer_temp[i] = cells_and_reference_positions.size();
-                }
-            }
+        if (perform_handshake)
+          return Utilities::pack(request_buffer_temp, false);
+        else
+          return {};
+      };
 
-          if (perform_handshake)
-            request_buffer = Utilities::pack(request_buffer_temp, false);
-        },
-        [&](const unsigned int       other_rank,
-            const std::vector<char> &recv_buffer) {
-          if (perform_handshake)
-            {
-              const auto recv_buffer_unpacked =
-                Utilities::unpack<std::vector<unsigned int>>(recv_buffer,
-                                                             false);
-
-              const auto other_rank_index = translate(other_rank);
-
-              for (unsigned int i = 0; i < recv_buffer_unpacked.size(); ++i)
-                for (unsigned int j = 0; j < recv_buffer_unpacked[i]; ++j)
-                  recv_components.emplace_back(
-                    other_rank,
-                    potential_owners_indices
-                      [i + potential_owners_ptrs[other_rank_index]],
-                    numbers::invalid_unsigned_int);
-            }
-        });
+      const auto process_answer = [&](const unsigned int       other_rank,
+                                      const std::vector<char> &answer) {
+        if (perform_handshake)
+          {
+            const auto recv_buffer_unpacked =
+              Utilities::unpack<std::vector<unsigned int>>(answer, false);
+
+            const auto other_rank_index = translate(other_rank);
+
+            for (unsigned int i = 0; i < recv_buffer_unpacked.size(); ++i)
+              for (unsigned int j = 0; j < recv_buffer_unpacked[i]; ++j)
+                recv_components.emplace_back(
+                  other_rank,
+                  potential_owners_indices
+                    [i + potential_owners_ptrs[other_rank_index]],
+                  numbers::invalid_unsigned_int);
+          }
+      };
 
-      Utilities::MPI::ConsensusAlgorithms::Selector<char, char>(
-        process, cache.get_triangulation().get_communicator())
-        .run();
+      Utilities::MPI::ConsensusAlgorithms::Selector<char, char>().run(
+        potential_owners_ranks,
+        create_request,
+        answer_request,
+        process_answer,
+        cache.get_triangulation().get_communicator());
 
       if (true)
         {
index 1ea4609460e1c627d4b2930a343f15feb4ae7ac0..c65388a3ae75792a0764dfb9fd8053cf3a4487a7 100644 (file)
@@ -113,37 +113,42 @@ namespace TriangulationDescription
           const MPI_Comm &                                   comm,
           const bool vertices_have_unique_ids)
         {
-          dealii::Utilities::MPI::ConsensusAlgorithms::AnonymousProcess<char,
-                                                                        char>
-            process(
-              [&]() { return relevant_processes; },
-              [&](const unsigned int other_rank,
-                  std::vector<char> &send_buffer) {
-                const auto ptr = std::find(relevant_processes.begin(),
-                                           relevant_processes.end(),
-                                           other_rank);
-
-                Assert(ptr != relevant_processes.end(), ExcInternalError());
-
-                const auto other_rank_index =
-                  std::distance(relevant_processes.begin(), ptr);
-
-                send_buffer =
-                  dealii::Utilities::pack(description_temp[other_rank_index],
-                                          false);
-              },
-              [&](const unsigned int &,
-                  const std::vector<char> &recv_buffer,
-                  std::vector<char> &) {
-                this->merge(
-                  dealii::Utilities::unpack<DescriptionTemp<dim, spacedim>>(
-                    recv_buffer, false),
-                  vertices_have_unique_ids);
-              });
-
-          dealii::Utilities::MPI::ConsensusAlgorithms::Selector<char, char>(
-            process, comm)
-            .run();
+          const auto create_request = [&](const unsigned int other_rank) {
+            const auto ptr = std::find(relevant_processes.begin(),
+                                       relevant_processes.end(),
+                                       other_rank);
+
+            Assert(ptr != relevant_processes.end(), ExcInternalError());
+
+            const auto other_rank_index =
+              std::distance(relevant_processes.begin(), ptr);
+
+            return dealii::Utilities::pack(description_temp[other_rank_index],
+                                           false);
+          };
+
+          const auto answer_request = [&](const unsigned int,
+                                          const std::vector<char> &request) {
+            this->merge(
+              dealii::Utilities::unpack<DescriptionTemp<dim, spacedim>>(request,
+                                                                        false),
+              vertices_have_unique_ids);
+            return std::vector<char>();
+          };
+
+          const auto process_answer = [](const unsigned int,
+                                         const std::vector<char> &answer) {
+            (void)answer;
+            Assert(answer.size() == 0, ExcInternalError());
+          };
+
+
+          dealii::Utilities::MPI::ConsensusAlgorithms::Selector<char, char>()
+            .run(relevant_processes,
+                 create_request,
+                 answer_request,
+                 process_answer,
+                 comm);
         }
 
         /**

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.