]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Provide a serial implementation of ConsensusAlgorithms::Interface 11230/head
authorPeter Munch <peterrmuench@gmail.com>
Mon, 23 Nov 2020 07:34:54 +0000 (08:34 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Mon, 23 Nov 2020 09:32:41 +0000 (10:32 +0100)
include/deal.II/base/mpi_consensus_algorithms.h
include/deal.II/base/mpi_consensus_algorithms.templates.h
source/base/mpi_consensus_algorithms.cc

index eb5c98e532e57738818586bba42cc5048761ddb9..707644ad002ed9c784bf2f694a68577f5edacc53 100644 (file)
@@ -190,6 +190,11 @@ namespace Utilities
          */
         const MPI_Comm &comm;
 
+        /**
+         * Cache if job supports MPI.
+         */
+        const bool job_supports_mpi;
+
         /**
          * Rank of this process.
          */
@@ -440,6 +445,29 @@ namespace Utilities
         clean_up_and_end_communication();
       };
 
+      /**
+       * A serial fall back for the above classes to allow programming
+       * independently of whether MPI is used or not.
+       */
+      template <typename T1, typename T2>
+      class Serial : public Interface<T1, T2>
+      {
+      public:
+        /**
+         * Constructor.
+         *
+         * @param process Process to be run during consensus algorithm.
+         * @param comm MPI Communicator (ignored)
+         */
+        Serial(Process<T1, T2> &process, const MPI_Comm &comm);
+
+        /**
+         * @copydoc Interface::run()
+         */
+        virtual void
+        run() override;
+      };
+
       /**
        * A class which delegates its task to other
        * ConsensusAlgorithms::Interface implementations depending on the number
index 010f0697851a3f2a4ec374f7a88645a66bcfd0e4..725c611ea0a3e738cb30d02a00e0d27f3238a529 100644 (file)
@@ -83,10 +83,9 @@ namespace Utilities
                                    const MPI_Comm & comm)
         : process(process)
         , comm(comm)
-        , my_rank(Utilities::MPI::job_supports_mpi() ? this_mpi_process(comm) :
-                                                       0)
-        , n_procs(Utilities::MPI::job_supports_mpi() ? n_mpi_processes(comm) :
-                                                       1)
+        , job_supports_mpi(Utilities::MPI::job_supports_mpi())
+        , my_rank(job_supports_mpi ? this_mpi_process(comm) : 0)
+        , n_procs(job_supports_mpi ? n_mpi_processes(comm) : 1)
       {}
 
 
@@ -555,6 +554,37 @@ namespace Utilities
 
 
 
+      template <typename T1, typename T2>
+      Serial<T1, T2>::Serial(Process<T1, T2> &process, const MPI_Comm &comm)
+        : Interface<T1, T2>(process, comm)
+      {}
+
+
+
+      template <typename T1, typename T2>
+      void
+      Serial<T1, T2>::run()
+      {
+        const auto targets = this->process.compute_targets();
+
+        if (targets.size() == 0)
+          return; // nothing to do
+
+        AssertDimension(targets[0], 0);
+
+        std::vector<T1> send_buffer;
+        std::vector<T2> recv_buffer;
+        std::vector<T2> request_buffer;
+        std::vector<T1> buffer_recv;
+
+        this->process.create_request(0, send_buffer);
+        this->process.prepare_buffer_for_answer(0, recv_buffer);
+        this->process.answer_request(0, buffer_recv, request_buffer);
+        this->process.read_answer(0, recv_buffer);
+      }
+
+
+
       template <typename T1, typename T2>
       Selector<T1, T2>::Selector(Process<T1, T2> &process, const MPI_Comm &comm)
         : Interface<T1, T2>(process, comm)
@@ -577,6 +607,8 @@ namespace Utilities
 #endif
           if (this->n_procs > 1)
           consensus_algo.reset(new PEX<T1, T2>(process, comm));
+        else
+          consensus_algo.reset(new Serial<T1, T2>(process, comm));
       }
 
 
@@ -585,8 +617,7 @@ namespace Utilities
       void
       Selector<T1, T2>::run()
       {
-        if (consensus_algo)
-          consensus_algo->run();
+        consensus_algo->run();
       }
 
 
index 1eebaa29eaa0d8b3148c69a3bba73defa9a6765c..0f1b560198744595f7534411ad48b91370b36f35 100644 (file)
@@ -29,6 +29,8 @@ namespace Utilities
 
       template class PEX<unsigned int, unsigned int>;
 
+      template class Serial<unsigned int, unsigned int>;
+
       template class Selector<unsigned int, unsigned int>;
 
 
@@ -44,6 +46,10 @@ namespace Utilities
         std::pair<types::global_dof_index, types::global_dof_index>,
         unsigned int>;
 
+      template class Serial<
+        std::pair<types::global_dof_index, types::global_dof_index>,
+        unsigned int>;
+
       template class PEX<
         std::pair<types::global_dof_index, types::global_dof_index>,
         unsigned int>;
@@ -53,6 +59,8 @@ namespace Utilities
 
       template class NBX<types::global_dof_index, unsigned int>;
 
+      template class Serial<types::global_dof_index, unsigned int>;
+
       template class PEX<types::global_dof_index, unsigned int>;
 
       template class Selector<types::global_dof_index, unsigned int>;

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.