]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move ConsensusAlgorithm related functions to mpi_consensus_algorithm.h 9600/head
authorPeter Munch <peterrmuench@gmail.com>
Tue, 3 Mar 2020 10:10:10 +0000 (11:10 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Wed, 4 Mar 2020 19:40:47 +0000 (20:40 +0100)
include/deal.II/base/mpi.h
include/deal.II/base/mpi.templates.h
include/deal.II/base/mpi_compute_index_owner_internal.h
include/deal.II/base/mpi_consensus_algorithm.h
include/deal.II/base/mpi_consensus_algorithm.templates.h [new file with mode: 0644]
source/base/CMakeLists.txt
source/base/mpi.cc
source/base/mpi_consensus_algorithm.cc [new file with mode: 0644]

index ecb6ff11da11b30912864ea9c6bf9b6c3f8328c9..6431e56ea15de8c1d590be5c03a74a264790d2ab 100644 (file)
@@ -975,461 +975,6 @@ namespace Utilities
            const T &          object_to_send,
            const unsigned int root_process = 0);
 
-    /**
-     * An interface to be able to use the ConsensusAlgorithm classes. The main
-     * functionality of the implementations is to return a list of process ranks
-     * this process wants data from and to deal with the optional payload of the
-     * messages sent/received by the ConsensusAlgorithm classes.
-     *
-     * There are two kinds of messages:
-     * - send/request message: A message consisting of a data request
-     *   which should be answered by another process. This message is
-     *   considered as a request message by the receiving rank.
-     * - recv message: The answer to a send/request message.
-     *
-     * @tparam T1 the type of the elements of the vector to sent
-     * @tparam T2 the type of the elements of the vector to received
-     *
-     * @note Since the payloads of the messages are optional, users have
-     *       to deal with buffers themselves. The ConsensusAlgorithm classes 1)
-     *       deliver only references to empty vectors (of size 0)
-     *       the data to be sent can be inserted to or read from, and
-     *       2) communicate these vectors blindly.
-     *
-     * @author Peter Munch, 2019
-     */
-    template <typename T1, typename T2>
-    class ConsensusAlgorithmProcess
-    {
-    public:
-      /**
-       * Destructor.
-       */
-      virtual ~ConsensusAlgorithmProcess() = default;
-
-      /**
-       * @return A vector of ranks this process wants to send a request to.
-       *
-       * @note This is the only method which has to be implemented since the
-       *       payloads of the messages are optional.
-       */
-      virtual std::vector<unsigned int>
-      compute_targets() = 0;
-
-      /**
-       * Add to the request to the process with the specified rank a payload.
-       *
-       * @param[in]  other_rank rank of the process
-       * @param[out] send_buffer data to be sent part of the request (optional)
-       *
-       * @note The buffer is empty. Before using it, you have to set its size.
-       */
-      virtual void
-      create_request(const unsigned int other_rank,
-                     std::vector<T1> &  send_buffer);
-
-      /**
-       * Prepare the buffer where the payload of the answer of the request to
-       * the process with the specified rank is saved in. The most obvious task
-       * is to resize the buffer, since it is empty when the function is called.
-       *
-       * @param[in]  other_rank rank of the process
-       * @param[out] recv_buffer data to be sent part of the request (optional)
-       */
-      virtual void
-      prepare_buffer_for_answer(const unsigned int other_rank,
-                                std::vector<T2> &  recv_buffer);
-
-      /**
-       * Prepare the buffer where the payload of the answer of the request to
-       * the process with the specified rank is saved in.
-       *
-       * @param[in]  other_rank rank of the process
-       * @param[in]  buffer_recv received payload (optional)
-       * @param[out] request_buffer payload to be sent as part of the request
-       *             (optional)
-       *
-       * @note The request_buffer is empty. Before using it, you have to set
-       *       its size.
-       */
-      virtual void
-      answer_request(const unsigned int     other_rank,
-                     const std::vector<T1> &buffer_recv,
-                     std::vector<T2> &      request_buffer);
-
-      /**
-       * Process the payload of the answer of the request to the process with
-       * the specified rank.
-       *
-       * @param[in] other_rank rank of the process
-       * @param[in] recv_buffer data to be sent part of the request (optional)
-       */
-      virtual void
-      read_answer(const unsigned int     other_rank,
-                  const std::vector<T2> &recv_buffer);
-    };
-
-
-
-    /**
-     * A base class for algorithms that implement the task of coming up with
-     * communication patterns to retrieve data from other processes in a
-     * dynamic-sparse way. In computer science, this is often called a
-     * <a href="https://en.wikipedia.org/wiki/Consensus_algorithm">consensus
-     * problem</a>.
-     *
-     * Dynamic-sparse means in this context:
-     * - By the time this function is called, the other processes do
-     *   not know yet that they have to answer requests.
-     * - Each process only has to communicate with a small subset of
-     *   processes of the MPI communicator.
-     *
-     * Naturally, the user has to provide:
-     * - A communicator.
-     * - For each rank a list of ranks of processes this process should
-     *   communicate to.
-     * - Functionality to pack/unpack data to be sent/received.
-     *
-     * This base class only introduces a basic interface to achieve
-     * these goals, while derived classes implement different algorithms
-     * to actually compute such communication patterns.
-     * The last two features of the list above this paragraph are implemented
-     * in classes derived from ConsensusAlgorithmProcess.
-     *
-     * @tparam T1 The type of the elements of the vector to be sent.
-     * @tparam T2 The type of the elements of the vector to be received.
-     *
-     * @author Peter Munch, 2019
-     */
-    template <typename T1, typename T2>
-    class ConsensusAlgorithm
-    {
-    public:
-      ConsensusAlgorithm(ConsensusAlgorithmProcess<T1, T2> &process,
-                         const MPI_Comm &                   comm);
-
-      /**
-       * Destructor.
-       */
-      virtual ~ConsensusAlgorithm() = default;
-
-      /**
-       * Run consensus algorithm.
-       */
-      virtual void
-      run() = 0;
-
-    protected:
-      /**
-       * Reference to the process provided by the user.
-       */
-      ConsensusAlgorithmProcess<T1, T2> &process;
-
-      /**
-       * MPI communicator.
-       */
-      const MPI_Comm &comm;
-
-      /**
-       * Rank of this process.
-       */
-      const unsigned int my_rank;
-
-      /**
-       * Number of processes in the communicator.
-       */
-      const unsigned int n_procs;
-    };
-
-
-    /**
-     * This class implements a concrete algorithm for the ConsensusAlgorithm
-     * base class, using only point-to-point communications and a single
-     * IBarrier.
-     *
-     * @note This class closely follows the paper Hoefler, Siebert, Lumsdaine
-     *       "Scalable Communication Protocols for Dynamic Sparse Data
-     *       Exchange". Since the algorithm shown there is not considering
-     *       payloads, the algorithm has been modified here in such a way that
-     *       synchronous sends (Issend) have been replaced by equivalent
-     *       Isend/Irecv, where Irecv receives the answer to a request (with
-     *       payload).
-     *
-     * @tparam T1 The type of the elements of the vector to be sent.
-     * @tparam T2 The type of the elements of the vector to be received.
-     *
-     * @author Peter Munch, 2019
-     */
-    template <typename T1, typename T2>
-    class ConsensusAlgorithm_NBX : public ConsensusAlgorithm<T1, T2>
-    {
-    public:
-      /**
-       * Constructor.
-       *
-       * @param process Process to be run during consensus algorithm.
-       * @param comm MPI Communicator
-       */
-      ConsensusAlgorithm_NBX(ConsensusAlgorithmProcess<T1, T2> &process,
-                             const MPI_Comm &                   comm);
-
-      /**
-       * Destructor.
-       */
-      virtual ~ConsensusAlgorithm_NBX() = default;
-
-      /**
-       * @copydoc ConsensusAlgorithm::run()
-       */
-      virtual void
-      run() override;
-
-    private:
-#ifdef DEAL_II_WITH_MPI
-      /**
-       * List of processes this process wants to send requests to.
-       */
-      std::vector<unsigned int> targets;
-
-      /**
-       * Buffers for sending requests.
-       */
-      std::vector<std::vector<T1>> send_buffers;
-
-      /**
-       * Requests for sending requests.
-       */
-      std::vector<MPI_Request> send_requests;
-
-      /**
-       * Buffers for receiving answers to requests.
-       */
-      std::vector<std::vector<T2>> recv_buffers;
-
-
-      /**
-       * Requests for receiving answers to requests.
-       */
-      std::vector<MPI_Request> recv_requests;
-
-      /**
-       * Buffers for sending answers to requests.
-       */
-      std::vector<std::unique_ptr<std::vector<T2>>> request_buffers;
-
-      /**
-       * Requests for sending answers to requests.
-       */
-      std::vector<std::unique_ptr<MPI_Request>> request_requests;
-
-      // request for barrier
-      MPI_Request barrier_request;
-#endif
-
-#ifdef DEBUG
-      /**
-       * List of processes who have made a request to this process.
-       */
-      std::set<unsigned int> requesting_processes;
-#endif
-
-      /**
-       * Check if all request answers have been received by this rank.
-       */
-      bool
-      check_own_state();
-
-      /**
-       * Signal to all other ranks that this rank has received all request
-       * answers via entering IBarrier.
-       */
-      void
-      signal_finish();
-
-      /**
-       * Check if all ranks have received all their request answers, i.e.
-       * all ranks have reached the IBarrier.
-       */
-      bool
-      check_global_state();
-
-      /**
-       * A request message from another rank has been received: process the
-       * request and send an answer.
-       */
-      void
-      answer_requests();
-
-      /**
-       * Start to send all requests via ISend and post IRecvs for the incoming
-       * answer messages.
-       */
-      void
-      start_communication();
-
-      /**
-       * After all rank has received all answers, the MPI data structures can be
-       * freed and the received answers can be processed.
-       */
-      void
-      clean_up_and_end_communication();
-    };
-
-    /**
-     * This class implements a concrete algorithm for the ConsensusAlgorithm
-     * base class, using a two step approach. In
-     * the first step the source ranks are determined and in the second step
-     * a static sparse data exchange is performed.
-     *
-     * @note In contrast to ConsensusAlgorithm_NBX, this class splits the same
-     *       task into two distinct steps. In the first step, all processes are
-     *       identified who want to send a request to this process. In the
-     *       second step, the data is exchanged. However, since - in the second
-     *       step - now it is clear how many requests have to be answered, i.e.
-     *       when this process can stop waiting for requests, no IBarrier is
-     *       needed.
-     *
-     * @note The function
-     *       Utilities::MPI::compute_point_to_point_communication_pattern() is
-     *       used to determine the source processes, which implements a
-     *       PEX-algorithm from Hoefner et al., "Scalable Communication
-     *       Protocols for Dynamic Sparse Data Exchange".
-     *
-     * @tparam T1 The type of the elements of the vector to be sent.
-     * @tparam T2 The type of the elements of the vector to be received.
-     *
-     * @author Peter Munch, 2019
-     */
-    template <typename T1, typename T2>
-    class ConsensusAlgorithm_PEX : public ConsensusAlgorithm<T1, T2>
-    {
-    public:
-      /**
-       * Constructor.
-       *
-       * @param process Process to be run during consensus algorithm.
-       * @param comm MPI Communicator
-       */
-      ConsensusAlgorithm_PEX(ConsensusAlgorithmProcess<T1, T2> &process,
-                             const MPI_Comm &                   comm);
-
-      /**
-       * Destructor.
-       */
-      virtual ~ConsensusAlgorithm_PEX() = default;
-
-      /**
-       * @copydoc ConsensusAlgorithm::run()
-       */
-      virtual void
-      run() override;
-
-    private:
-#ifdef DEAL_II_WITH_MPI
-      /**
-       * List of ranks of processes this processes wants to send a request to.
-       */
-      std::vector<unsigned int> targets;
-
-      /**
-       * List of ranks of processes wanting to send a request to this process.
-       */
-      std::vector<unsigned int> sources;
-
-      // data structures to send and receive requests
-
-      /**
-       * Buffers for sending requests.
-       */
-      std::vector<std::vector<T1>> send_buffers;
-
-      /**
-       * Buffers for receiving answers to requests.
-       */
-      std::vector<std::vector<T2>> recv_buffers;
-
-      /**
-       * Requests for sending requests and receiving answers to requests.
-       */
-      std::vector<MPI_Request> send_and_recv_buffers;
-
-      /**
-       * Buffers for sending answers to requests.
-       */
-      std::vector<std::vector<T2>> requests_buffers;
-
-      /**
-       * Requests for sending answers to requests.
-       */
-      std::vector<MPI_Request> requests_answers;
-#endif
-
-      /**
-       * The ith request message from another rank has been received: process
-       * the request and send an answer.
-       */
-      void
-      answer_requests(int index);
-
-      /**
-       * Start to send all requests via ISend and post IRecvs for the incoming
-       * answer messages.
-       */
-      unsigned int
-      start_communication();
-
-      /**
-       * After all answers have been exchanged, the MPI data structures can be
-       * freed and the received answers can be processed.
-       */
-      void
-      clean_up_and_end_communication();
-    };
-
-    /**
-     * A class which delegates its task to other ConsensusAlgorithm
-     * implementations depending on the number of processes in the
-     * MPI communicator. For a small number of processes it uses
-     * ConsensusAlgorithm_PEX and for a large number of processes
-     * ConsensusAlgorithm_NBX. The threshold depends if the program is
-     * compiled in debug or release mode.
-     *
-     * @tparam T1 The type of the elements of the vector to be sent.
-     * @tparam T2 The type of the elements of the vector to be received.
-     *
-     * @author Peter Munch, 2019
-     */
-    template <typename T1, typename T2>
-    class ConsensusAlgorithmSelector : public ConsensusAlgorithm<T1, T2>
-    {
-    public:
-      /**
-       * Constructor.
-       *
-       * @param process Process to be run during consensus algorithm.
-       * @param comm MPI Communicator.
-       */
-      ConsensusAlgorithmSelector(ConsensusAlgorithmProcess<T1, T2> &process,
-                                 const MPI_Comm &                   comm);
-
-      /**
-       * Destructor.
-       */
-      virtual ~ConsensusAlgorithmSelector() = default;
-
-      /**
-       * @copydoc ConsensusAlgorithm::run()
-       *
-       * @note The function call is delegated to another ConsensusAlgorithm implementation.
-       */
-      virtual void
-      run() override;
-
-    private:
-      // Pointer to the actual ConsensusAlgorithm implementation.
-      std::shared_ptr<ConsensusAlgorithm<T1, T2>> consensus_algo;
-    };
-
     /**
      * Given a partitioned index set space, compute the owning MPI process rank
      * of each element of a second index set according to the partitioned index
index 5cc190f05f27f1022b72aeff114fea728986b70c..10888963672b91070792212b45c55557adfc3162 100644 (file)
@@ -413,559 +413,6 @@ namespace Utilities
       internal::all_reduce(MPI_MIN, values, mpi_communicator, minima);
     }
 
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithmProcess<T1, T2>::answer_request(const unsigned int,
-                                                      const std::vector<T1> &,
-                                                      std::vector<T2> &)
-    {
-      // nothing to do
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithmProcess<T1, T2>::create_request(const unsigned int,
-                                                      std::vector<T1> &)
-    {
-      // nothing to do
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithmProcess<T1, T2>::prepare_buffer_for_answer(
-      const unsigned int,
-      std::vector<T2> &)
-    {
-      // nothing to do
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithmProcess<T1, T2>::read_answer(const unsigned int,
-                                                   const std::vector<T2> &)
-    {
-      // nothing to do
-    }
-
-
-
-    template <typename T1, typename T2>
-    ConsensusAlgorithm<T1, T2>::ConsensusAlgorithm(
-      ConsensusAlgorithmProcess<T1, T2> &process,
-      const MPI_Comm &                   comm)
-      : process(process)
-      , comm(comm)
-      , my_rank(this_mpi_process(comm))
-      , n_procs(n_mpi_processes(comm))
-    {}
-
-
-
-    template <typename T1, typename T2>
-    ConsensusAlgorithm_NBX<T1, T2>::ConsensusAlgorithm_NBX(
-      ConsensusAlgorithmProcess<T1, T2> &process,
-      const MPI_Comm &                   comm)
-      : ConsensusAlgorithm<T1, T2>(process, comm)
-    {}
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithm_NBX<T1, T2>::run()
-    {
-      static CollectiveMutex      mutex;
-      CollectiveMutex::ScopedLock lock(mutex, this->comm);
-
-      // 1) send requests and start receiving the answers
-      start_communication();
-
-      // 2) answer requests and check if all requests of this process have been
-      //    answered
-      while (!check_own_state())
-        answer_requests();
-
-      // 3) signal to all other processes that all requests of this process have
-      //    been answered
-      signal_finish();
-
-      // 4) nevertheless, this process has to keep on answering (potential)
-      //    incoming requests until all processes have received the
-      //    answer to all requests
-      while (!check_global_state())
-        answer_requests();
-
-      // 5) process the answer to all requests
-      clean_up_and_end_communication();
-    }
-
-
-
-    template <typename T1, typename T2>
-    bool
-    ConsensusAlgorithm_NBX<T1, T2>::check_own_state()
-    {
-#ifdef DEAL_II_WITH_MPI
-      int        all_receive_requests_are_done;
-      const auto ierr = MPI_Testall(recv_requests.size(),
-                                    recv_requests.data(),
-                                    &all_receive_requests_are_done,
-                                    MPI_STATUSES_IGNORE);
-      AssertThrowMPI(ierr);
-
-      return all_receive_requests_are_done;
-#else
-      return true;
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithm_NBX<T1, T2>::signal_finish()
-    {
-#ifdef DEAL_II_WITH_MPI
-#  if DEAL_II_MPI_VERSION_GTE(3, 0)
-      const auto ierr = MPI_Ibarrier(this->comm, &barrier_request);
-      AssertThrowMPI(ierr);
-#  else
-      AssertThrow(
-        false,
-        ExcMessage(
-          "ConsensusAlgorithm_NBX uses MPI 3.0 features. You should compile with at least MPI 3.0."));
-#  endif
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    bool
-    ConsensusAlgorithm_NBX<T1, T2>::check_global_state()
-    {
-#ifdef DEAL_II_WITH_MPI
-      int        all_ranks_reached_barrier;
-      const auto ierr = MPI_Test(&barrier_request,
-                                 &all_ranks_reached_barrier,
-                                 MPI_STATUSES_IGNORE);
-      AssertThrowMPI(ierr);
-      return all_ranks_reached_barrier;
-#else
-      return true;
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithm_NBX<T1, T2>::answer_requests()
-    {
-#ifdef DEAL_II_WITH_MPI
-
-      const int tag_request =
-        Utilities::MPI::internal::Tags::consensus_algorithm_nbx_answer_request;
-      const int tag_deliver =
-        Utilities::MPI::internal::Tags::consensus_algorithm_nbx_process_deliver;
-
-      // check if there is a request pending
-      MPI_Status status;
-      int        request_is_pending;
-      const auto ierr = MPI_Iprobe(
-        MPI_ANY_SOURCE, tag_request, this->comm, &request_is_pending, &status);
-      AssertThrowMPI(ierr);
-
-      if (request_is_pending) // request is pending
-        {
-          // get rank of requesting process
-          const auto other_rank = status.MPI_SOURCE;
-
-#  ifdef DEBUG
-          Assert(requesting_processes.find(other_rank) ==
-                   requesting_processes.end(),
-                 ExcMessage("Process is requesting a second time!"));
-          requesting_processes.insert(other_rank);
-#  endif
-
-          std::vector<T1> buffer_recv;
-          // get size of of incoming message
-          int  number_amount;
-          auto ierr = MPI_Get_count(&status, MPI_BYTE, &number_amount);
-          AssertThrowMPI(ierr);
-
-          // allocate memory for incoming message
-          Assert(number_amount % sizeof(T1) == 0, ExcInternalError());
-          buffer_recv.resize(number_amount / sizeof(T1));
-          ierr = MPI_Recv(buffer_recv.data(),
-                          number_amount,
-                          MPI_BYTE,
-                          other_rank,
-                          tag_request,
-                          this->comm,
-                          &status);
-          AssertThrowMPI(ierr);
-
-          // allocate memory for answer message
-          request_buffers.emplace_back(
-            std_cxx14::make_unique<std::vector<T2>>());
-          request_requests.emplace_back(std_cxx14::make_unique<MPI_Request>());
-
-          // process request
-          auto &request_buffer = *request_buffers.back();
-          this->process.answer_request(other_rank, buffer_recv, request_buffer);
-
-          // start to send answer back
-          ierr = MPI_Isend(request_buffer.data(),
-                           request_buffer.size() * sizeof(T2),
-                           MPI_BYTE,
-                           other_rank,
-                           tag_deliver,
-                           this->comm,
-                           request_requests.back().get());
-          AssertThrowMPI(ierr);
-        }
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithm_NBX<T1, T2>::start_communication()
-    {
-#ifdef DEAL_II_WITH_MPI
-      // 1)
-      targets              = this->process.compute_targets();
-      const auto n_targets = targets.size();
-
-      const int tag_request =
-        Utilities::MPI::internal::Tags::consensus_algorithm_nbx_answer_request;
-      const int tag_deliver =
-        Utilities::MPI::internal::Tags::consensus_algorithm_nbx_process_deliver;
-
-      // 2) allocate memory
-      recv_buffers.resize(n_targets);
-      recv_requests.resize(n_targets);
-      send_requests.resize(n_targets);
-      send_buffers.resize(n_targets);
-
-      {
-        // 4) send and receive
-        for (unsigned int i = 0; i < n_targets; i++)
-          {
-            const unsigned int rank  = targets[i];
-            const unsigned int index = i;
-
-            // translate index set to a list of pairs
-            auto &send_buffer = send_buffers[index];
-            this->process.create_request(rank, send_buffer);
-
-            // start to send data
-            auto ierr = MPI_Isend(send_buffer.data(),
-                                  send_buffer.size() * sizeof(T1),
-                                  MPI_BYTE,
-                                  rank,
-                                  tag_request,
-                                  this->comm,
-                                  &send_requests[index]);
-            AssertThrowMPI(ierr);
-
-            // start to receive data
-            auto &recv_buffer = recv_buffers[index];
-            this->process.prepare_buffer_for_answer(rank, recv_buffer);
-            ierr = MPI_Irecv(recv_buffer.data(),
-                             recv_buffer.size() * sizeof(T2),
-                             MPI_BYTE,
-                             rank,
-                             tag_deliver,
-                             this->comm,
-                             &recv_requests[index]);
-            AssertThrowMPI(ierr);
-          }
-      }
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithm_NBX<T1, T2>::clean_up_and_end_communication()
-    {
-#ifdef DEAL_II_WITH_MPI
-      // clean up
-      {
-        if (send_requests.size() > 0)
-          {
-            const int ierr = MPI_Waitall(send_requests.size(),
-                                         send_requests.data(),
-                                         MPI_STATUSES_IGNORE);
-            AssertThrowMPI(ierr);
-          }
-
-        if (recv_requests.size() > 0)
-          {
-            const int ierr = MPI_Waitall(recv_requests.size(),
-                                         recv_requests.data(),
-                                         MPI_STATUSES_IGNORE);
-            AssertThrowMPI(ierr);
-          }
-
-
-        const int ierr = MPI_Wait(&barrier_request, MPI_STATUS_IGNORE);
-        AssertThrowMPI(ierr);
-
-        for (auto &i : request_requests)
-          {
-            const auto ierr = MPI_Wait(i.get(), MPI_STATUS_IGNORE);
-            AssertThrowMPI(ierr);
-          }
-
-#  ifdef DEBUG
-        // note: IBarrier seems to make problem during testing, this additional
-        // Barrier seems to help
-        MPI_Barrier(this->comm);
-#  endif
-      }
-
-      // unpack data
-      {
-        for (unsigned int i = 0; i < targets.size(); i++)
-          this->process.read_answer(targets[i], recv_buffers[i]);
-      }
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    ConsensusAlgorithm_PEX<T1, T2>::ConsensusAlgorithm_PEX(
-      ConsensusAlgorithmProcess<T1, T2> &process,
-      const MPI_Comm &                   comm)
-      : ConsensusAlgorithm<T1, T2>(process, comm)
-    {}
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithm_PEX<T1, T2>::run()
-    {
-      static CollectiveMutex      mutex;
-      CollectiveMutex::ScopedLock lock(mutex, this->comm);
-
-      // 1) send requests and start receiving the answers
-      //    especially determine how many requests are expected
-      const unsigned int n_requests = start_communication();
-
-      // 2) answer requests
-      for (unsigned int request = 0; request < n_requests; request++)
-        answer_requests(request);
-
-      // 3) process answers
-      clean_up_and_end_communication();
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithm_PEX<T1, T2>::answer_requests(int index)
-    {
-#ifdef DEAL_II_WITH_MPI
-      const int tag_request =
-        Utilities::MPI::internal::Tags::consensus_algorithm_pex_answer_request;
-      const int tag_deliver =
-        Utilities::MPI::internal::Tags::consensus_algorithm_pex_process_deliver;
-
-      MPI_Status status;
-      auto ierr = MPI_Probe(MPI_ANY_SOURCE, tag_request, this->comm, &status);
-      AssertThrowMPI(ierr);
-
-      // get rank of incoming message
-      const auto other_rank = status.MPI_SOURCE;
-
-      std::vector<T1> buffer_recv;
-
-      // get size of incoming message
-      int number_amount;
-      ierr = MPI_Get_count(&status, MPI_BYTE, &number_amount);
-      AssertThrowMPI(ierr);
-
-      // allocate memory for incoming message
-      Assert(number_amount % sizeof(T1) == 0, ExcInternalError());
-      buffer_recv.resize(number_amount / sizeof(T1));
-      ierr = MPI_Recv(buffer_recv.data(),
-                      number_amount,
-                      MPI_BYTE,
-                      other_rank,
-                      tag_request,
-                      this->comm,
-                      &status);
-      AssertThrowMPI(ierr);
-
-      // process request
-      auto &request_buffer = requests_buffers[index];
-      this->process.answer_request(other_rank, buffer_recv, request_buffer);
-
-      // start to send answer back
-      ierr = MPI_Isend(request_buffer.data(),
-                       request_buffer.size() * sizeof(T2),
-                       MPI_BYTE,
-                       other_rank,
-                       tag_deliver,
-                       this->comm,
-                       &requests_answers[index]);
-      AssertThrowMPI(ierr);
-#else
-      (void)index;
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    unsigned int
-    ConsensusAlgorithm_PEX<T1, T2>::start_communication()
-    {
-#ifdef DEAL_II_WITH_MPI
-      // 1) determine with which processes this process wants to communicate
-      targets = this->process.compute_targets();
-
-      const int tag_request =
-        Utilities::MPI::internal::Tags::consensus_algorithm_pex_answer_request;
-      const int tag_deliver =
-        Utilities::MPI::internal::Tags::consensus_algorithm_pex_process_deliver;
-
-      // 2) determine who wants to communicate with this process
-      sources =
-        compute_point_to_point_communication_pattern(this->comm, targets);
-
-      const auto n_targets = targets.size();
-      const auto n_sources = sources.size();
-
-      // 2) allocate memory
-      recv_buffers.resize(n_targets);
-      send_buffers.resize(n_targets);
-      send_and_recv_buffers.resize(2 * n_targets);
-
-      requests_answers.resize(n_sources);
-      requests_buffers.resize(n_sources);
-
-      // 4) send and receive
-      for (unsigned int i = 0; i < n_targets; i++)
-        {
-          const unsigned int rank = targets[i];
-
-          // pack data which should be sent
-          auto &send_buffer = send_buffers[i];
-          this->process.create_request(rank, send_buffer);
-
-          // start to send data
-          auto ierr = MPI_Isend(send_buffer.data(),
-                                send_buffer.size() * sizeof(T1),
-                                MPI_BYTE,
-                                rank,
-                                tag_request,
-                                this->comm,
-                                &send_and_recv_buffers[n_targets + i]);
-          AssertThrowMPI(ierr);
-
-          // start to receive data
-          auto &recv_buffer = recv_buffers[i];
-          this->process.prepare_buffer_for_answer(rank, recv_buffer);
-          ierr = MPI_Irecv(recv_buffer.data(),
-                           recv_buffer.size() * sizeof(T2),
-                           MPI_BYTE,
-                           rank,
-                           tag_deliver,
-                           this->comm,
-                           &send_and_recv_buffers[i]);
-          AssertThrowMPI(ierr);
-        }
-
-      return sources.size();
-#else
-      return 0;
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithm_PEX<T1, T2>::clean_up_and_end_communication()
-    {
-#ifdef DEAL_II_WITH_MPI
-      // finalize all MPI_Requests
-      if (send_and_recv_buffers.size() > 0)
-        {
-          auto ierr = MPI_Waitall(send_and_recv_buffers.size(),
-                                  send_and_recv_buffers.data(),
-                                  MPI_STATUSES_IGNORE);
-          AssertThrowMPI(ierr);
-        }
-
-      if (requests_answers.size() > 0)
-        {
-          auto ierr = MPI_Waitall(requests_answers.size(),
-                                  requests_answers.data(),
-                                  MPI_STATUSES_IGNORE);
-          AssertThrowMPI(ierr);
-        }
-
-      // unpack received data
-      for (unsigned int i = 0; i < targets.size(); i++)
-        this->process.read_answer(targets[i], recv_buffers[i]);
-#endif
-    }
-
-
-
-    template <typename T1, typename T2>
-    ConsensusAlgorithmSelector<T1, T2>::ConsensusAlgorithmSelector(
-      ConsensusAlgorithmProcess<T1, T2> &process,
-      const MPI_Comm &                   comm)
-      : ConsensusAlgorithm<T1, T2>(process, comm)
-    {
-      // Depending on the number of processes we switch between implementations.
-      // We reduce the threshold for debug mode to be able to test also the
-      // non-blocking implementation. This feature is tested by:
-      // tests/multigrid/transfer_matrix_free_06.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=15.output
-#ifdef DEAL_II_WITH_MPI
-#  if DEAL_II_MPI_VERSION_GTE(3, 0)
-#    ifdef DEBUG
-      if (Utilities::MPI::n_mpi_processes(comm) > 14)
-#    else
-      if (Utilities::MPI::n_mpi_processes(comm) > 99)
-#    endif
-        consensus_algo.reset(new ConsensusAlgorithm_NBX<T1, T2>(process, comm));
-      else
-#  endif
-#endif
-        consensus_algo.reset(new ConsensusAlgorithm_PEX<T1, T2>(process, comm));
-    }
-
-
-
-    template <typename T1, typename T2>
-    void
-    ConsensusAlgorithmSelector<T1, T2>::run()
-    {
-      consensus_algo->run();
-    }
-
-
-
   } // end of namespace MPI
 } // end of namespace Utilities
 
index 7ecce658dccf1f2b3be483d3cf6ed2465b64780d..d853cd8739c4524ef809c67f81c53b0c08e94c17 100644 (file)
@@ -19,6 +19,7 @@
 #include <deal.II/base/config.h>
 
 #include <deal.II/base/mpi.h>
+#include <deal.II/base/mpi_consensus_algorithm.h>
 
 DEAL_II_NAMESPACE_OPEN
 
index e9c13eb7990e18d00833a2e22b02618dcc7d605f..bfd7b3a274cfdc1f42f8e0f44c99ee3fc8fb9a47 100644 (file)
@@ -28,6 +28,461 @@ namespace Utilities
 {
   namespace MPI
   {
+    /**
+     * An interface to be able to use the ConsensusAlgorithm classes. The main
+     * functionality of the implementations is to return a list of process ranks
+     * this process wants data from and to deal with the optional payload of the
+     * messages sent/received by the ConsensusAlgorithm classes.
+     *
+     * There are two kinds of messages:
+     * - send/request message: A message consisting of a data request
+     *   which should be answered by another process. This message is
+     *   considered as a request message by the receiving rank.
+     * - recv message: The answer to a send/request message.
+     *
+     * @tparam T1 the type of the elements of the vector to sent
+     * @tparam T2 the type of the elements of the vector to received
+     *
+     * @note Since the payloads of the messages are optional, users have
+     *       to deal with buffers themselves. The ConsensusAlgorithm classes 1)
+     *       deliver only references to empty vectors (of size 0)
+     *       the data to be sent can be inserted to or read from, and
+     *       2) communicate these vectors blindly.
+     *
+     * @author Peter Munch, 2019
+     */
+    template <typename T1, typename T2>
+    class ConsensusAlgorithmProcess
+    {
+    public:
+      /**
+       * Destructor.
+       */
+      virtual ~ConsensusAlgorithmProcess() = default;
+
+      /**
+       * @return A vector of ranks this process wants to send a request to.
+       *
+       * @note This is the only method which has to be implemented since the
+       *       payloads of the messages are optional.
+       */
+      virtual std::vector<unsigned int>
+      compute_targets() = 0;
+
+      /**
+       * Add to the request to the process with the specified rank a payload.
+       *
+       * @param[in]  other_rank rank of the process
+       * @param[out] send_buffer data to be sent part of the request (optional)
+       *
+       * @note The buffer is empty. Before using it, you have to set its size.
+       */
+      virtual void
+      create_request(const unsigned int other_rank,
+                     std::vector<T1> &  send_buffer);
+
+      /**
+       * Prepare the buffer where the payload of the answer of the request to
+       * the process with the specified rank is saved in. The most obvious task
+       * is to resize the buffer, since it is empty when the function is called.
+       *
+       * @param[in]  other_rank rank of the process
+       * @param[out] recv_buffer data to be sent part of the request (optional)
+       */
+      virtual void
+      prepare_buffer_for_answer(const unsigned int other_rank,
+                                std::vector<T2> &  recv_buffer);
+
+      /**
+       * Prepare the buffer where the payload of the answer of the request to
+       * the process with the specified rank is saved in.
+       *
+       * @param[in]  other_rank rank of the process
+       * @param[in]  buffer_recv received payload (optional)
+       * @param[out] request_buffer payload to be sent as part of the request
+       *             (optional)
+       *
+       * @note The request_buffer is empty. Before using it, you have to set
+       *       its size.
+       */
+      virtual void
+      answer_request(const unsigned int     other_rank,
+                     const std::vector<T1> &buffer_recv,
+                     std::vector<T2> &      request_buffer);
+
+      /**
+       * Process the payload of the answer of the request to the process with
+       * the specified rank.
+       *
+       * @param[in] other_rank rank of the process
+       * @param[in] recv_buffer data to be sent part of the request (optional)
+       */
+      virtual void
+      read_answer(const unsigned int     other_rank,
+                  const std::vector<T2> &recv_buffer);
+    };
+
+
+
+    /**
+     * A base class for algorithms that implement the task of coming up with
+     * communication patterns to retrieve data from other processes in a
+     * dynamic-sparse way. In computer science, this is often called a
+     * <a href="https://en.wikipedia.org/wiki/Consensus_algorithm">consensus
+     * problem</a>.
+     *
+     * Dynamic-sparse means in this context:
+     * - By the time this function is called, the other processes do
+     *   not know yet that they have to answer requests.
+     * - Each process only has to communicate with a small subset of
+     *   processes of the MPI communicator.
+     *
+     * Naturally, the user has to provide:
+     * - A communicator.
+     * - For each rank a list of ranks of processes this process should
+     *   communicate to.
+     * - Functionality to pack/unpack data to be sent/received.
+     *
+     * This base class only introduces a basic interface to achieve
+     * these goals, while derived classes implement different algorithms
+     * to actually compute such communication patterns.
+     * The last two features of the list above this paragraph are implemented
+     * in classes derived from ConsensusAlgorithmProcess.
+     *
+     * @tparam T1 The type of the elements of the vector to be sent.
+     * @tparam T2 The type of the elements of the vector to be received.
+     *
+     * @author Peter Munch, 2019
+     */
+    template <typename T1, typename T2>
+    class ConsensusAlgorithm
+    {
+    public:
+      ConsensusAlgorithm(ConsensusAlgorithmProcess<T1, T2> &process,
+                         const MPI_Comm &                   comm);
+
+      /**
+       * Destructor.
+       */
+      virtual ~ConsensusAlgorithm() = default;
+
+      /**
+       * Run consensus algorithm.
+       */
+      virtual void
+      run() = 0;
+
+    protected:
+      /**
+       * Reference to the process provided by the user.
+       */
+      ConsensusAlgorithmProcess<T1, T2> &process;
+
+      /**
+       * MPI communicator.
+       */
+      const MPI_Comm &comm;
+
+      /**
+       * Rank of this process.
+       */
+      const unsigned int my_rank;
+
+      /**
+       * Number of processes in the communicator.
+       */
+      const unsigned int n_procs;
+    };
+
+
+    /**
+     * This class implements a concrete algorithm for the ConsensusAlgorithm
+     * base class, using only point-to-point communications and a single
+     * IBarrier.
+     *
+     * @note This class closely follows the paper Hoefler, Siebert, Lumsdaine
+     *       "Scalable Communication Protocols for Dynamic Sparse Data
+     *       Exchange". Since the algorithm shown there is not considering
+     *       payloads, the algorithm has been modified here in such a way that
+     *       synchronous sends (Issend) have been replaced by equivalent
+     *       Isend/Irecv, where Irecv receives the answer to a request (with
+     *       payload).
+     *
+     * @tparam T1 The type of the elements of the vector to be sent.
+     * @tparam T2 The type of the elements of the vector to be received.
+     *
+     * @author Peter Munch, 2019
+     */
+    template <typename T1, typename T2>
+    class ConsensusAlgorithm_NBX : public ConsensusAlgorithm<T1, T2>
+    {
+    public:
+      /**
+       * Constructor.
+       *
+       * @param process Process to be run during consensus algorithm.
+       * @param comm MPI Communicator
+       */
+      ConsensusAlgorithm_NBX(ConsensusAlgorithmProcess<T1, T2> &process,
+                             const MPI_Comm &                   comm);
+
+      /**
+       * Destructor.
+       */
+      virtual ~ConsensusAlgorithm_NBX() = default;
+
+      /**
+       * @copydoc ConsensusAlgorithm::run()
+       */
+      virtual void
+      run() override;
+
+    private:
+#ifdef DEAL_II_WITH_MPI
+      /**
+       * List of processes this process wants to send requests to.
+       */
+      std::vector<unsigned int> targets;
+
+      /**
+       * Buffers for sending requests.
+       */
+      std::vector<std::vector<T1>> send_buffers;
+
+      /**
+       * Requests for sending requests.
+       */
+      std::vector<MPI_Request> send_requests;
+
+      /**
+       * Buffers for receiving answers to requests.
+       */
+      std::vector<std::vector<T2>> recv_buffers;
+
+
+      /**
+       * Requests for receiving answers to requests.
+       */
+      std::vector<MPI_Request> recv_requests;
+
+      /**
+       * Buffers for sending answers to requests.
+       */
+      std::vector<std::unique_ptr<std::vector<T2>>> request_buffers;
+
+      /**
+       * Requests for sending answers to requests.
+       */
+      std::vector<std::unique_ptr<MPI_Request>> request_requests;
+
+      // request for barrier
+      MPI_Request barrier_request;
+#endif
+
+#ifdef DEBUG
+      /**
+       * List of processes who have made a request to this process.
+       */
+      std::set<unsigned int> requesting_processes;
+#endif
+
+      /**
+       * Check if all request answers have been received by this rank.
+       */
+      bool
+      check_own_state();
+
+      /**
+       * Signal to all other ranks that this rank has received all request
+       * answers via entering IBarrier.
+       */
+      void
+      signal_finish();
+
+      /**
+       * Check if all ranks have received all their request answers, i.e.
+       * all ranks have reached the IBarrier.
+       */
+      bool
+      check_global_state();
+
+      /**
+       * A request message from another rank has been received: process the
+       * request and send an answer.
+       */
+      void
+      answer_requests();
+
+      /**
+       * Start to send all requests via ISend and post IRecvs for the incoming
+       * answer messages.
+       */
+      void
+      start_communication();
+
+      /**
+       * After all rank has received all answers, the MPI data structures can be
+       * freed and the received answers can be processed.
+       */
+      void
+      clean_up_and_end_communication();
+    };
+
+    /**
+     * This class implements a concrete algorithm for the ConsensusAlgorithm
+     * base class, using a two step approach. In
+     * the first step the source ranks are determined and in the second step
+     * a static sparse data exchange is performed.
+     *
+     * @note In contrast to ConsensusAlgorithm_NBX, this class splits the same
+     *       task into two distinct steps. In the first step, all processes are
+     *       identified who want to send a request to this process. In the
+     *       second step, the data is exchanged. However, since - in the second
+     *       step - now it is clear how many requests have to be answered, i.e.
+     *       when this process can stop waiting for requests, no IBarrier is
+     *       needed.
+     *
+     * @note The function
+     *       Utilities::MPI::compute_point_to_point_communication_pattern() is
+     *       used to determine the source processes, which implements a
+     *       PEX-algorithm from Hoefner et al., "Scalable Communication
+     *       Protocols for Dynamic Sparse Data Exchange".
+     *
+     * @tparam T1 The type of the elements of the vector to be sent.
+     * @tparam T2 The type of the elements of the vector to be received.
+     *
+     * @author Peter Munch, 2019
+     */
+    template <typename T1, typename T2>
+    class ConsensusAlgorithm_PEX : public ConsensusAlgorithm<T1, T2>
+    {
+    public:
+      /**
+       * Constructor.
+       *
+       * @param process Process to be run during consensus algorithm.
+       * @param comm MPI Communicator
+       */
+      ConsensusAlgorithm_PEX(ConsensusAlgorithmProcess<T1, T2> &process,
+                             const MPI_Comm &                   comm);
+
+      /**
+       * Destructor.
+       */
+      virtual ~ConsensusAlgorithm_PEX() = default;
+
+      /**
+       * @copydoc ConsensusAlgorithm::run()
+       */
+      virtual void
+      run() override;
+
+    private:
+#ifdef DEAL_II_WITH_MPI
+      /**
+       * List of ranks of processes this processes wants to send a request to.
+       */
+      std::vector<unsigned int> targets;
+
+      /**
+       * List of ranks of processes wanting to send a request to this process.
+       */
+      std::vector<unsigned int> sources;
+
+      // data structures to send and receive requests
+
+      /**
+       * Buffers for sending requests.
+       */
+      std::vector<std::vector<T1>> send_buffers;
+
+      /**
+       * Buffers for receiving answers to requests.
+       */
+      std::vector<std::vector<T2>> recv_buffers;
+
+      /**
+       * Requests for sending requests and receiving answers to requests.
+       */
+      std::vector<MPI_Request> send_and_recv_buffers;
+
+      /**
+       * Buffers for sending answers to requests.
+       */
+      std::vector<std::vector<T2>> requests_buffers;
+
+      /**
+       * Requests for sending answers to requests.
+       */
+      std::vector<MPI_Request> requests_answers;
+#endif
+
+      /**
+       * The ith request message from another rank has been received: process
+       * the request and send an answer.
+       */
+      void
+      answer_requests(int index);
+
+      /**
+       * Start to send all requests via ISend and post IRecvs for the incoming
+       * answer messages.
+       */
+      unsigned int
+      start_communication();
+
+      /**
+       * After all answers have been exchanged, the MPI data structures can be
+       * freed and the received answers can be processed.
+       */
+      void
+      clean_up_and_end_communication();
+    };
+
+    /**
+     * A class which delegates its task to other ConsensusAlgorithm
+     * implementations depending on the number of processes in the
+     * MPI communicator. For a small number of processes it uses
+     * ConsensusAlgorithm_PEX and for a large number of processes
+     * ConsensusAlgorithm_NBX. The threshold depends if the program is
+     * compiled in debug or release mode.
+     *
+     * @tparam T1 The type of the elements of the vector to be sent.
+     * @tparam T2 The type of the elements of the vector to be received.
+     *
+     * @author Peter Munch, 2019
+     */
+    template <typename T1, typename T2>
+    class ConsensusAlgorithmSelector : public ConsensusAlgorithm<T1, T2>
+    {
+    public:
+      /**
+       * Constructor.
+       *
+       * @param process Process to be run during consensus algorithm.
+       * @param comm MPI Communicator.
+       */
+      ConsensusAlgorithmSelector(ConsensusAlgorithmProcess<T1, T2> &process,
+                                 const MPI_Comm &                   comm);
+
+      /**
+       * Destructor.
+       */
+      virtual ~ConsensusAlgorithmSelector() = default;
+
+      /**
+       * @copydoc ConsensusAlgorithm::run()
+       *
+       * @note The function call is delegated to another ConsensusAlgorithm implementation.
+       */
+      virtual void
+      run() override;
+
+    private:
+      // Pointer to the actual ConsensusAlgorithm implementation.
+      std::shared_ptr<ConsensusAlgorithm<T1, T2>> consensus_algo;
+    };
+
     /**
      * This class implements Utilities::MPI::ConsensusAlgorithmProcess,
      * using user-provided function wrappers.
diff --git a/include/deal.II/base/mpi_consensus_algorithm.templates.h b/include/deal.II/base/mpi_consensus_algorithm.templates.h
new file mode 100644 (file)
index 0000000..3b39939
--- /dev/null
@@ -0,0 +1,597 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2011 - 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_mpi_consensus_algorithm_templates_h
+#define dealii_mpi_consensus_algorithm_templates_h
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/exceptions.h>
+#include <deal.II/base/mpi.h>
+#include <deal.II/base/mpi_consensus_algorithm.h>
+#include <deal.II/base/symmetric_tensor.h>
+#include <deal.II/base/tensor.h>
+
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/lapack_full_matrix.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <deal.II/lac/vector.h>
+
+#include <vector>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Utilities
+{
+  namespace MPI
+  {
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithmProcess<T1, T2>::answer_request(const unsigned int,
+                                                      const std::vector<T1> &,
+                                                      std::vector<T2> &)
+    {
+      // nothing to do
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithmProcess<T1, T2>::create_request(const unsigned int,
+                                                      std::vector<T1> &)
+    {
+      // nothing to do
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithmProcess<T1, T2>::prepare_buffer_for_answer(
+      const unsigned int,
+      std::vector<T2> &)
+    {
+      // nothing to do
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithmProcess<T1, T2>::read_answer(const unsigned int,
+                                                   const std::vector<T2> &)
+    {
+      // nothing to do
+    }
+
+
+
+    template <typename T1, typename T2>
+    ConsensusAlgorithm<T1, T2>::ConsensusAlgorithm(
+      ConsensusAlgorithmProcess<T1, T2> &process,
+      const MPI_Comm &                   comm)
+      : process(process)
+      , comm(comm)
+      , my_rank(this_mpi_process(comm))
+      , n_procs(n_mpi_processes(comm))
+    {}
+
+
+
+    template <typename T1, typename T2>
+    ConsensusAlgorithm_NBX<T1, T2>::ConsensusAlgorithm_NBX(
+      ConsensusAlgorithmProcess<T1, T2> &process,
+      const MPI_Comm &                   comm)
+      : ConsensusAlgorithm<T1, T2>(process, comm)
+    {}
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithm_NBX<T1, T2>::run()
+    {
+      static CollectiveMutex      mutex;
+      CollectiveMutex::ScopedLock lock(mutex, this->comm);
+
+      // 1) send requests and start receiving the answers
+      start_communication();
+
+      // 2) answer requests and check if all requests of this process have been
+      //    answered
+      while (!check_own_state())
+        answer_requests();
+
+      // 3) signal to all other processes that all requests of this process have
+      //    been answered
+      signal_finish();
+
+      // 4) nevertheless, this process has to keep on answering (potential)
+      //    incoming requests until all processes have received the
+      //    answer to all requests
+      while (!check_global_state())
+        answer_requests();
+
+      // 5) process the answer to all requests
+      clean_up_and_end_communication();
+    }
+
+
+
+    template <typename T1, typename T2>
+    bool
+    ConsensusAlgorithm_NBX<T1, T2>::check_own_state()
+    {
+#ifdef DEAL_II_WITH_MPI
+      int        all_receive_requests_are_done;
+      const auto ierr = MPI_Testall(recv_requests.size(),
+                                    recv_requests.data(),
+                                    &all_receive_requests_are_done,
+                                    MPI_STATUSES_IGNORE);
+      AssertThrowMPI(ierr);
+
+      return all_receive_requests_are_done;
+#else
+      return true;
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithm_NBX<T1, T2>::signal_finish()
+    {
+#ifdef DEAL_II_WITH_MPI
+#  if DEAL_II_MPI_VERSION_GTE(3, 0)
+      const auto ierr = MPI_Ibarrier(this->comm, &barrier_request);
+      AssertThrowMPI(ierr);
+#  else
+      AssertThrow(
+        false,
+        ExcMessage(
+          "ConsensusAlgorithm_NBX uses MPI 3.0 features. You should compile with at least MPI 3.0."));
+#  endif
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    bool
+    ConsensusAlgorithm_NBX<T1, T2>::check_global_state()
+    {
+#ifdef DEAL_II_WITH_MPI
+      int        all_ranks_reached_barrier;
+      const auto ierr = MPI_Test(&barrier_request,
+                                 &all_ranks_reached_barrier,
+                                 MPI_STATUSES_IGNORE);
+      AssertThrowMPI(ierr);
+      return all_ranks_reached_barrier;
+#else
+      return true;
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithm_NBX<T1, T2>::answer_requests()
+    {
+#ifdef DEAL_II_WITH_MPI
+
+      const int tag_request =
+        Utilities::MPI::internal::Tags::consensus_algorithm_nbx_answer_request;
+      const int tag_deliver =
+        Utilities::MPI::internal::Tags::consensus_algorithm_nbx_process_deliver;
+
+      // check if there is a request pending
+      MPI_Status status;
+      int        request_is_pending;
+      const auto ierr = MPI_Iprobe(
+        MPI_ANY_SOURCE, tag_request, this->comm, &request_is_pending, &status);
+      AssertThrowMPI(ierr);
+
+      if (request_is_pending) // request is pending
+        {
+          // get rank of requesting process
+          const auto other_rank = status.MPI_SOURCE;
+
+#  ifdef DEBUG
+          Assert(requesting_processes.find(other_rank) ==
+                   requesting_processes.end(),
+                 ExcMessage("Process is requesting a second time!"));
+          requesting_processes.insert(other_rank);
+#  endif
+
+          std::vector<T1> buffer_recv;
+          // get size of of incoming message
+          int  number_amount;
+          auto ierr = MPI_Get_count(&status, MPI_BYTE, &number_amount);
+          AssertThrowMPI(ierr);
+
+          // allocate memory for incoming message
+          Assert(number_amount % sizeof(T1) == 0, ExcInternalError());
+          buffer_recv.resize(number_amount / sizeof(T1));
+          ierr = MPI_Recv(buffer_recv.data(),
+                          number_amount,
+                          MPI_BYTE,
+                          other_rank,
+                          tag_request,
+                          this->comm,
+                          &status);
+          AssertThrowMPI(ierr);
+
+          // allocate memory for answer message
+          request_buffers.emplace_back(
+            std_cxx14::make_unique<std::vector<T2>>());
+          request_requests.emplace_back(std_cxx14::make_unique<MPI_Request>());
+
+          // process request
+          auto &request_buffer = *request_buffers.back();
+          this->process.answer_request(other_rank, buffer_recv, request_buffer);
+
+          // start to send answer back
+          ierr = MPI_Isend(request_buffer.data(),
+                           request_buffer.size() * sizeof(T2),
+                           MPI_BYTE,
+                           other_rank,
+                           tag_deliver,
+                           this->comm,
+                           request_requests.back().get());
+          AssertThrowMPI(ierr);
+        }
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithm_NBX<T1, T2>::start_communication()
+    {
+#ifdef DEAL_II_WITH_MPI
+      // 1)
+      targets              = this->process.compute_targets();
+      const auto n_targets = targets.size();
+
+      const int tag_request =
+        Utilities::MPI::internal::Tags::consensus_algorithm_nbx_answer_request;
+      const int tag_deliver =
+        Utilities::MPI::internal::Tags::consensus_algorithm_nbx_process_deliver;
+
+      // 2) allocate memory
+      recv_buffers.resize(n_targets);
+      recv_requests.resize(n_targets);
+      send_requests.resize(n_targets);
+      send_buffers.resize(n_targets);
+
+      {
+        // 4) send and receive
+        for (unsigned int i = 0; i < n_targets; i++)
+          {
+            const unsigned int rank  = targets[i];
+            const unsigned int index = i;
+
+            // translate index set to a list of pairs
+            auto &send_buffer = send_buffers[index];
+            this->process.create_request(rank, send_buffer);
+
+            // start to send data
+            auto ierr = MPI_Isend(send_buffer.data(),
+                                  send_buffer.size() * sizeof(T1),
+                                  MPI_BYTE,
+                                  rank,
+                                  tag_request,
+                                  this->comm,
+                                  &send_requests[index]);
+            AssertThrowMPI(ierr);
+
+            // start to receive data
+            auto &recv_buffer = recv_buffers[index];
+            this->process.prepare_buffer_for_answer(rank, recv_buffer);
+            ierr = MPI_Irecv(recv_buffer.data(),
+                             recv_buffer.size() * sizeof(T2),
+                             MPI_BYTE,
+                             rank,
+                             tag_deliver,
+                             this->comm,
+                             &recv_requests[index]);
+            AssertThrowMPI(ierr);
+          }
+      }
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithm_NBX<T1, T2>::clean_up_and_end_communication()
+    {
+#ifdef DEAL_II_WITH_MPI
+      // clean up
+      {
+        if (send_requests.size() > 0)
+          {
+            const int ierr = MPI_Waitall(send_requests.size(),
+                                         send_requests.data(),
+                                         MPI_STATUSES_IGNORE);
+            AssertThrowMPI(ierr);
+          }
+
+        if (recv_requests.size() > 0)
+          {
+            const int ierr = MPI_Waitall(recv_requests.size(),
+                                         recv_requests.data(),
+                                         MPI_STATUSES_IGNORE);
+            AssertThrowMPI(ierr);
+          }
+
+
+        const int ierr = MPI_Wait(&barrier_request, MPI_STATUS_IGNORE);
+        AssertThrowMPI(ierr);
+
+        for (auto &i : request_requests)
+          {
+            const auto ierr = MPI_Wait(i.get(), MPI_STATUS_IGNORE);
+            AssertThrowMPI(ierr);
+          }
+
+#  ifdef DEBUG
+        // note: IBarrier seems to make problem during testing, this additional
+        // Barrier seems to help
+        MPI_Barrier(this->comm);
+#  endif
+      }
+
+      // unpack data
+      {
+        for (unsigned int i = 0; i < targets.size(); i++)
+          this->process.read_answer(targets[i], recv_buffers[i]);
+      }
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    ConsensusAlgorithm_PEX<T1, T2>::ConsensusAlgorithm_PEX(
+      ConsensusAlgorithmProcess<T1, T2> &process,
+      const MPI_Comm &                   comm)
+      : ConsensusAlgorithm<T1, T2>(process, comm)
+    {}
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithm_PEX<T1, T2>::run()
+    {
+      static CollectiveMutex      mutex;
+      CollectiveMutex::ScopedLock lock(mutex, this->comm);
+
+      // 1) send requests and start receiving the answers
+      //    especially determine how many requests are expected
+      const unsigned int n_requests = start_communication();
+
+      // 2) answer requests
+      for (unsigned int request = 0; request < n_requests; request++)
+        answer_requests(request);
+
+      // 3) process answers
+      clean_up_and_end_communication();
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithm_PEX<T1, T2>::answer_requests(int index)
+    {
+#ifdef DEAL_II_WITH_MPI
+      const int tag_request =
+        Utilities::MPI::internal::Tags::consensus_algorithm_pex_answer_request;
+      const int tag_deliver =
+        Utilities::MPI::internal::Tags::consensus_algorithm_pex_process_deliver;
+
+      MPI_Status status;
+      auto ierr = MPI_Probe(MPI_ANY_SOURCE, tag_request, this->comm, &status);
+      AssertThrowMPI(ierr);
+
+      // get rank of incoming message
+      const auto other_rank = status.MPI_SOURCE;
+
+      std::vector<T1> buffer_recv;
+
+      // get size of incoming message
+      int number_amount;
+      ierr = MPI_Get_count(&status, MPI_BYTE, &number_amount);
+      AssertThrowMPI(ierr);
+
+      // allocate memory for incoming message
+      Assert(number_amount % sizeof(T1) == 0, ExcInternalError());
+      buffer_recv.resize(number_amount / sizeof(T1));
+      ierr = MPI_Recv(buffer_recv.data(),
+                      number_amount,
+                      MPI_BYTE,
+                      other_rank,
+                      tag_request,
+                      this->comm,
+                      &status);
+      AssertThrowMPI(ierr);
+
+      // process request
+      auto &request_buffer = requests_buffers[index];
+      this->process.answer_request(other_rank, buffer_recv, request_buffer);
+
+      // start to send answer back
+      ierr = MPI_Isend(request_buffer.data(),
+                       request_buffer.size() * sizeof(T2),
+                       MPI_BYTE,
+                       other_rank,
+                       tag_deliver,
+                       this->comm,
+                       &requests_answers[index]);
+      AssertThrowMPI(ierr);
+#else
+      (void)index;
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    unsigned int
+    ConsensusAlgorithm_PEX<T1, T2>::start_communication()
+    {
+#ifdef DEAL_II_WITH_MPI
+      // 1) determine with which processes this process wants to communicate
+      targets = this->process.compute_targets();
+
+      const int tag_request =
+        Utilities::MPI::internal::Tags::consensus_algorithm_pex_answer_request;
+      const int tag_deliver =
+        Utilities::MPI::internal::Tags::consensus_algorithm_pex_process_deliver;
+
+      // 2) determine who wants to communicate with this process
+      sources =
+        compute_point_to_point_communication_pattern(this->comm, targets);
+
+      const auto n_targets = targets.size();
+      const auto n_sources = sources.size();
+
+      // 2) allocate memory
+      recv_buffers.resize(n_targets);
+      send_buffers.resize(n_targets);
+      send_and_recv_buffers.resize(2 * n_targets);
+
+      requests_answers.resize(n_sources);
+      requests_buffers.resize(n_sources);
+
+      // 4) send and receive
+      for (unsigned int i = 0; i < n_targets; i++)
+        {
+          const unsigned int rank = targets[i];
+
+          // pack data which should be sent
+          auto &send_buffer = send_buffers[i];
+          this->process.create_request(rank, send_buffer);
+
+          // start to send data
+          auto ierr = MPI_Isend(send_buffer.data(),
+                                send_buffer.size() * sizeof(T1),
+                                MPI_BYTE,
+                                rank,
+                                tag_request,
+                                this->comm,
+                                &send_and_recv_buffers[n_targets + i]);
+          AssertThrowMPI(ierr);
+
+          // start to receive data
+          auto &recv_buffer = recv_buffers[i];
+          this->process.prepare_buffer_for_answer(rank, recv_buffer);
+          ierr = MPI_Irecv(recv_buffer.data(),
+                           recv_buffer.size() * sizeof(T2),
+                           MPI_BYTE,
+                           rank,
+                           tag_deliver,
+                           this->comm,
+                           &send_and_recv_buffers[i]);
+          AssertThrowMPI(ierr);
+        }
+
+      return sources.size();
+#else
+      return 0;
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithm_PEX<T1, T2>::clean_up_and_end_communication()
+    {
+#ifdef DEAL_II_WITH_MPI
+      // finalize all MPI_Requests
+      if (send_and_recv_buffers.size() > 0)
+        {
+          auto ierr = MPI_Waitall(send_and_recv_buffers.size(),
+                                  send_and_recv_buffers.data(),
+                                  MPI_STATUSES_IGNORE);
+          AssertThrowMPI(ierr);
+        }
+
+      if (requests_answers.size() > 0)
+        {
+          auto ierr = MPI_Waitall(requests_answers.size(),
+                                  requests_answers.data(),
+                                  MPI_STATUSES_IGNORE);
+          AssertThrowMPI(ierr);
+        }
+
+      // unpack received data
+      for (unsigned int i = 0; i < targets.size(); i++)
+        this->process.read_answer(targets[i], recv_buffers[i]);
+#endif
+    }
+
+
+
+    template <typename T1, typename T2>
+    ConsensusAlgorithmSelector<T1, T2>::ConsensusAlgorithmSelector(
+      ConsensusAlgorithmProcess<T1, T2> &process,
+      const MPI_Comm &                   comm)
+      : ConsensusAlgorithm<T1, T2>(process, comm)
+    {
+      // Depending on the number of processes we switch between implementations.
+      // We reduce the threshold for debug mode to be able to test also the
+      // non-blocking implementation. This feature is tested by:
+      // tests/multigrid/transfer_matrix_free_06.with_mpi=true.with_p4est=true.with_trilinos=true.mpirun=15.output
+#ifdef DEAL_II_WITH_MPI
+#  if DEAL_II_MPI_VERSION_GTE(3, 0)
+#    ifdef DEBUG
+      if (Utilities::MPI::n_mpi_processes(comm) > 14)
+#    else
+      if (Utilities::MPI::n_mpi_processes(comm) > 99)
+#    endif
+        consensus_algo.reset(new ConsensusAlgorithm_NBX<T1, T2>(process, comm));
+      else
+#  endif
+#endif
+        consensus_algo.reset(new ConsensusAlgorithm_PEX<T1, T2>(process, comm));
+    }
+
+
+
+    template <typename T1, typename T2>
+    void
+    ConsensusAlgorithmSelector<T1, T2>::run()
+    {
+      consensus_algo->run();
+    }
+
+
+
+  } // end of namespace MPI
+} // end of namespace Utilities
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index be79de63f7b9ab51bb08c8d1cd81c7a443602dcd..ed2f0f819fd3c777fd84aa6a506fdd5cc422cf76 100644 (file)
@@ -45,6 +45,7 @@ SET(_unity_include_src
   logstream.cc
   hdf5.cc
   mpi.cc
+  mpi_consensus_algorithm.cc
   mu_parser_internal.cc
   multithread_info.cc
   named_selection.cc
index 487461cdbe2608f04dafd4730f8c246318ee52f4..799cbe1f814a18820063001ed41717048e36474c 100644 (file)
@@ -1045,10 +1045,6 @@ namespace Utilities
       return owning_ranks;
     }
 
-    template class ConsensusAlgorithmSelector<
-      std::pair<types::global_dof_index, types::global_dof_index>,
-      unsigned int>;
-
 
 
     CollectiveMutex::CollectiveMutex()
diff --git a/source/base/mpi_consensus_algorithm.cc b/source/base/mpi_consensus_algorithm.cc
new file mode 100644 (file)
index 0000000..b4ecfb9
--- /dev/null
@@ -0,0 +1,41 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2019 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/mpi_consensus_algorithm.templates.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+
+namespace Utilities
+{
+  namespace MPI
+  {
+    template class ConsensusAlgorithmProcess<unsigned int, unsigned int>;
+
+    template class ConsensusAlgorithmSelector<unsigned int, unsigned int>;
+
+
+    template class ConsensusAlgorithmProcess<
+      std::pair<types::global_dof_index, types::global_dof_index>,
+      unsigned int>;
+
+    template class ConsensusAlgorithmSelector<
+      std::pair<types::global_dof_index, types::global_dof_index>,
+      unsigned int>;
+
+  } // end of namespace MPI
+} // end of namespace Utilities
+
+DEAL_II_NAMESPACE_CLOSE

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.