Scatter/Gather mit ungleichmäßiger Verteilung

Die Operationen MPI_Scatter und MPI_Gather setzten voraus, dass jeder Prozess gleich viele Elemente (wie beispielsweise Matrixzeilen) erhält. Was kann getan werden, wenn die Zahl der Zeilen kein Vielfaches der Zahl der Prozesse ist, über die eine Matrix verteilt wird?

MPI bietet für diese Fälle Operationen an, bei denen zusätzlich die individuellen Häppchengrößen (zusammen mit den Offsets in displs[]) angegeben werden:

int MPI_Scatterv(const void* sendbuf, int sendcounts[], const int displs[],
                 MPI_Datatype sendtype,
                 void* recvbuf, int recvcount, MPI_Datatype recvtype,
                 int root, MPI_Comm comm);

int MPI_Gatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype,
                void* recvbuf, const int recvcounts[], const int displs[],
                MPI_Datatype recvtype, int root, MPI_Comm comm);

Aufgabe

Es erscheint sinnvoll, Funktionen für das zeilenweise Verteilen von Matrizen anzubieten, die die Parameter soweit wie möglich in Abhängigkeit der übergebenen Matrizen bestimmen. Die Aufteilung kann mit Hilfe eines hpc::aux::Slices-Objekt aus der Vorlesungsbibliothek erfolgen. Diese Objekte liefern die fertigen Werte für sendcounts[]/recvcounts[] und displs[].

Die beiden folgenden Funktionen sollen jeweils von A nach B kopieren. Bei scatter_by_row ist A die große zu verteilende Matrix und B der lokale Ausschnitt. Bei gather_by_row ist A der lokale Ausschnitt und B die große Matrix, die alle Abschnitte einsammelt.

In beiden Fällen ist A jeweils const deklariert. Da die MPI-Installation auf den Debian-Maschinen im E.44 etwas veraltet ist, erwarten MPI_Scatterv und MPI_Gatherv noch jeweils void* an Stelle von const void* beim ersten Argument. Das Problem lässt sich mit einem Cast lösen. Statt &A(0,0) einfach (void*) &A(0,0) übergeben. Dann wird es vom Übersetzer akzeptiert.

#include <type_traits>
#include <vector>
#include <mpi.h>
#include <hpc/matvec/gematrix.h>
#include <hpc/matvec/isgematrix.h>
#include <hpc/mpi/fundamental.h>
#include <hpc/aux/slices.h>

template<typename MA, typename MB>
typename std::enable_if<
   hpc::matvec::IsGeMatrix<MA>::value && hpc::matvec::IsGeMatrix<MB>::value &&
      std::is_same<typename MA::ElementType, typename MB::ElementType>::value,
   int>::type
scatter_by_row(const MA& A, MB& B, int root, MPI_Comm comm) {

   /* ... */
}

template<typename MA, typename MB>
typename std::enable_if<
   hpc::matvec::IsGeMatrix<MA>::value && hpc::matvec::IsGeMatrix<MB>::value &&
      std::is_same<typename MA::ElementType, typename MB::ElementType>::value,
   int>::type
gather_by_row(const MA& A, MB& B, int root, MPI_Comm comm) {

   /* ... */
}

Testen Sie die beiden Funktionen mit Hilfe des folgenden Testprogramms:

#include <mpi.h>
#include <hpc/mpi/matrix.h>
#include <hpc/matvec/gematrix.h>
#include <hpc/matvec/print.h>
#include <hpc/matvec/apply.h>

int main(int argc, char** argv) {
   MPI_Init(&argc, &argv);

   int nof_processes; MPI_Comm_size(MPI_COMM_WORLD, &nof_processes);
   int rank; MPI_Comm_rank(MPI_COMM_WORLD, &rank);

   using namespace hpc::matvec;
   using namespace hpc::mpi;
   using namespace hpc::aux;

   using Matrix = GeMatrix<double>;
   using Index = typename Matrix::Index;
   int share = 3;
   int num_rows = nof_processes * share + 1;
   int num_cols = 5;

   Matrix A(num_rows, num_cols); /* entire matrix */
   MPI_Datatype rowtype_A = get_row_type(A);

   Slices<int> slices(nof_processes, num_rows);
   Matrix B(slices.size(rank), num_cols,
      StorageOrder::RowMajor); /* individual share */
   MPI_Datatype rowtype_B = get_row_type(B);

   if (rank == 0) {
      apply(A, [](double& val, Index i, Index j) -> void {
         val = i * 100 + j;
      });
   }

   scatter_by_row(A, B, 0, MPI_COMM_WORLD);
   apply(B, [=](double& val, Index i, Index j) -> void {
      val += 10000 * (rank + 1);
   });
   gather_by_row(B, A, 0, MPI_COMM_WORLD);

   MPI_Finalize();

   if (rank == 0) {
      print(A, "A");
   }
}