From d3c741494c5f62b4856bd6dd80c1f02344d7bc22 Mon Sep 17 00:00:00 2001
From: Denis Davydov <davydden@gmail.com>
Date: Tue, 19 Mar 2019 20:39:35 +0100
Subject: [PATCH] add Utilities::MPI::create_ascending_partitioning()

---
 doc/news/changes/minor/20190319DenisDavydov   |  4 +
 include/deal.II/base/mpi.h                    | 13 +++
 source/base/mpi.cc                            | 33 ++++++++
 tests/base/mpi_create_partitioning_01.cc      | 80 +++++++++++++++++++
 ...mpi_create_partitioning_01.mpirun=2.output | 13 +++
 ...mpi_create_partitioning_01.mpirun=4.output | 43 ++++++++++
 tests/base/mpi_create_partitioning_01.output  |  4 +
 7 files changed, 190 insertions(+)
 create mode 100644 doc/news/changes/minor/20190319DenisDavydov
 create mode 100644 tests/base/mpi_create_partitioning_01.cc
 create mode 100644 tests/base/mpi_create_partitioning_01.mpirun=2.output
 create mode 100644 tests/base/mpi_create_partitioning_01.mpirun=4.output
 create mode 100644 tests/base/mpi_create_partitioning_01.output

diff --git a/doc/news/changes/minor/20190319DenisDavydov b/doc/news/changes/minor/20190319DenisDavydov
new file mode 100644
index 0000000000..ab8f90eec8
--- /dev/null
+++ b/doc/news/changes/minor/20190319DenisDavydov
@@ -0,0 +1,4 @@
+New: Add Utilities::MPI::create_ascending_partitioning() to create a one-to-one
+ascending partitioning from locally owned sizes.
+<br>
+(Denis Davydov, 2019/03/19)
diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h
index 4688b3c17f..74bef20899 100644
--- a/include/deal.II/base/mpi.h
+++ b/include/deal.II/base/mpi.h
@@ -90,6 +90,7 @@ template <int rank, int dim, typename Number>
 class SymmetricTensor;
 template <typename Number>
 class SparseMatrix;
+class IndexSet;
 
 namespace Utilities
 {
@@ -225,6 +226,18 @@ namespace Utilities
                  MPI_Comm *       new_comm);
 #endif
 
+    /**
+     * Given the number of locally owned elements @p local_size,
+     * create a 1:1 partitioning of the of elements across the MPI communicator @p comm.
+     * The total size of elements is the sum of @p local_size across the MPI communicator.
+     * Each process will store contiguous subset of indices, and the index set
+     * on process p+1 starts at the index one larger than the last one stored on
+     * process p.
+     */
+    std::vector<IndexSet>
+    create_ascending_partitioning(const MPI_Comm &           comm,
+                                  const IndexSet::size_type &local_size);
+
     /**
      * Return the sum over all processors of the value @p t. This function is
      * collective over all processors given in the
diff --git a/source/base/mpi.cc b/source/base/mpi.cc
index b8cd45e7c2..35928aa9cc 100644
--- a/source/base/mpi.cc
+++ b/source/base/mpi.cc
@@ -15,6 +15,7 @@
 
 
 #include <deal.II/base/exceptions.h>
+#include <deal.II/base/index_set.h>
 #include <deal.II/base/mpi.h>
 #include <deal.II/base/mpi.templates.h>
 #include <deal.II/base/multithread_info.h>
@@ -181,6 +182,30 @@ namespace Utilities
 
 
 
+    std::vector<IndexSet>
+    create_ascending_partitioning(const MPI_Comm &           comm,
+                                  const IndexSet::size_type &local_size)
+    {
+      const unsigned int                     n_proc = n_mpi_processes(comm);
+      const std::vector<IndexSet::size_type> sizes =
+        all_gather(comm, local_size);
+      const auto total_size =
+        std::accumulate(sizes.begin(), sizes.end(), IndexSet::size_type(0));
+
+      std::vector<IndexSet> res(n_proc, IndexSet(total_size));
+
+      IndexSet::size_type begin = 0;
+      for (unsigned int i = 0; i < n_proc; ++i)
+        {
+          res[i].add_range(begin, begin + sizes[i]);
+          begin = begin + sizes[i];
+        }
+
+      return res;
+    }
+
+
+
     std::vector<unsigned int>
     compute_point_to_point_communication_pattern(
       const MPI_Comm &                 mpi_comm,
@@ -484,6 +509,14 @@ namespace Utilities
     }
 
 
+    std::vector<IndexSet>
+    create_ascending_partitioning(const MPI_Comm & /*comm*/,
+                                  const IndexSet::size_type &local_size)
+    {
+      return std::vector<IndexSet>(1, complete_index_set(local_size));
+    }
+
+
     MPI_Comm
     duplicate_communicator(const MPI_Comm &mpi_communicator)
     {
diff --git a/tests/base/mpi_create_partitioning_01.cc b/tests/base/mpi_create_partitioning_01.cc
new file mode 100644
index 0000000000..6ca9f6de44
--- /dev/null
+++ b/tests/base/mpi_create_partitioning_01.cc
@@ -0,0 +1,80 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 - 2018 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.
+//
+// ---------------------------------------------------------------------
+
+// Test Utilities::MPI::create_ascending_partitioning()
+
+#include <deal.II/base/mpi.h>
+#include <deal.II/base/point.h>
+
+#include <vector>
+
+#include "../tests.h"
+
+void
+test()
+{
+  MPI_Comm   comm(MPI_COMM_WORLD);
+  const auto my_proc = Utilities::MPI::this_mpi_process(comm);
+  const auto n_proc  = Utilities::MPI::n_mpi_processes(comm);
+
+  const auto my_size = (my_proc == 2 ? 0 : my_proc * 2 + 5);
+  const auto vec = Utilities::MPI::create_ascending_partitioning(comm, my_size);
+
+  // correct number:
+  AssertThrow(vec.size() == n_proc, ExcInternalError());
+
+  // ascending one-to-one
+  AssertThrow(vec[my_proc].is_ascending_and_one_to_one(comm),
+              ExcInternalError());
+
+  // same size:
+  const auto size = vec[0].size();
+  for (unsigned int i = 1; i < vec.size(); ++i)
+    AssertThrow(size == vec[i].size(), ExcInternalError());
+
+  // 1:1
+  for (unsigned int i = 0; i < vec.size(); ++i)
+    for (unsigned int j = i + 1; j < vec.size(); ++j)
+      {
+        const IndexSet intersection = vec[i] & vec[j];
+        AssertThrow(intersection.is_empty(), ExcInternalError());
+      }
+
+  // union of all
+  IndexSet all(size);
+  for (unsigned int i = 0; i < vec.size(); ++i)
+    all.add_indices(vec[i]);
+
+  AssertThrow(all == complete_index_set(size), ExcInternalError());
+
+  // finally print:
+  deallog << "Total size=" << size << std::endl;
+  for (unsigned int p = 0; p < vec.size(); ++p)
+    {
+      deallog << "p=" << p << ":" << std::endl;
+      vec[p].print(deallog.get_file_stream());
+    }
+}
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  MPILogInitAll log;
+
+  test();
+  return 0;
+}
diff --git a/tests/base/mpi_create_partitioning_01.mpirun=2.output b/tests/base/mpi_create_partitioning_01.mpirun=2.output
new file mode 100644
index 0000000000..9075cba86e
--- /dev/null
+++ b/tests/base/mpi_create_partitioning_01.mpirun=2.output
@@ -0,0 +1,13 @@
+
+DEAL:0::Total size=12
+DEAL:0::p=0:
+{[0,4]}
+DEAL:0::p=1:
+{[5,11]}
+
+DEAL:1::Total size=12
+DEAL:1::p=0:
+{[0,4]}
+DEAL:1::p=1:
+{[5,11]}
+
diff --git a/tests/base/mpi_create_partitioning_01.mpirun=4.output b/tests/base/mpi_create_partitioning_01.mpirun=4.output
new file mode 100644
index 0000000000..7ac455243e
--- /dev/null
+++ b/tests/base/mpi_create_partitioning_01.mpirun=4.output
@@ -0,0 +1,43 @@
+
+DEAL:0::Total size=23
+DEAL:0::p=0:
+{[0,4]}
+DEAL:0::p=1:
+{[5,11]}
+DEAL:0::p=2:
+{}
+DEAL:0::p=3:
+{[12,22]}
+
+DEAL:1::Total size=23
+DEAL:1::p=0:
+{[0,4]}
+DEAL:1::p=1:
+{[5,11]}
+DEAL:1::p=2:
+{}
+DEAL:1::p=3:
+{[12,22]}
+
+
+DEAL:2::Total size=23
+DEAL:2::p=0:
+{[0,4]}
+DEAL:2::p=1:
+{[5,11]}
+DEAL:2::p=2:
+{}
+DEAL:2::p=3:
+{[12,22]}
+
+
+DEAL:3::Total size=23
+DEAL:3::p=0:
+{[0,4]}
+DEAL:3::p=1:
+{[5,11]}
+DEAL:3::p=2:
+{}
+DEAL:3::p=3:
+{[12,22]}
+
diff --git a/tests/base/mpi_create_partitioning_01.output b/tests/base/mpi_create_partitioning_01.output
new file mode 100644
index 0000000000..4d10ab83b3
--- /dev/null
+++ b/tests/base/mpi_create_partitioning_01.output
@@ -0,0 +1,4 @@
+
+DEAL:0::Total size=5
+DEAL:0::p=0:
+{[0,4]}
-- 
2.39.5