# Reconstructing Matrix scattered over many MPI ranks

I have a program that does some matrix computations using MPI (MPICH). Each rank has a slice of the matrix and does the computations on their slice to get a new slice of the matrix. Sometimes I need to synchronize the matrix on the manager for snapshotting, output, etc. so I need all ranks to send their current matrix slice back to the manager.

I've come up with a solution to do this, and it seems to work OK; but I wanted to know whether this is a good approach or if I'm doing something stupid since I've never worked with MPI before. One thing I'm particularly curious about is whether I'm deleting my requests vector correctly. Thanks in advance!

Here's the code:

#define __WORLD MPI_COMM_WORLD

// reconstructs the matrix on the manager node
void reconstructMatrix(double** current, const int COLS, const int ROWS,
const int MY_RANK, const int MPI_SIZE)
{
auto requests = new std::vector<MPI_Request*>();
requests->reserve(COLS);

// Manager receives data from all ranks
if (MY_RANK == 0)
{
// number of columns assigned to all ranks
const auto SLICE_SIZE = COLS / MPI_SIZE;

// receive from rank 1 to rank n-1
for (auto src = 1; src < MPI_SIZE; src++)
{
for (auto i = (SLICE_SIZE * src); (i < COLS) && (i < SLICE_SIZE * (src + 1)); i++)
{
auto temp = new MPI_Request();
MPI_Irecv(current[i], ROWS, MPI_DOUBLE, src, 1, __WORLD, temp);
requests->push_back(temp);
}
}
}
else // Workers send data to manager
{
// index 0 and COLS + 1 are ghost columns
for (auto i = 1; i <= COLS; i++)
{
auto temp = new MPI_Request();
MPI_Isend(current[i], ROWS, MPI_DOUBLE, 0, 1, __WORLD, temp);
requests->push_back(temp);
}
}
// wait for sync to finish
waitAll(requests);

requests->clear();
delete requests;
requests = nullptr;
}

// waits for all MPI requests and deletes the vector
void waitAll(std::vector<MPI_Request*>*& requests)
{
for (auto& req : *requests)
{
MPI_Status stat;
MPI_Wait(req, &stat);
if (stat.MPI_ERROR != MPI_SUCCESS)
{
std::cerr << "MPI Request failed!\n\tCode: " << stat.MPI_ERROR
<< "\n\tFrom: " << stat.MPI_SOURCE;
}
delete req;
}
}

• One thing I should mention is, I know sending a bunch of messages like this is inefficient. I saw that most of the time people store matrices in a 1D array with MPI to send many columns at once, but I usually only exchange boundary data (1 column) so I figured it would be nicer to store data in a 2D array for that. – nick Oct 26 '18 at 7:41

Note: I've written a matrix multiplication benchmark program in C and MPI:https://github.com/thefangbear/matrix-mpi. It uses a row-major format to store the two matrices in 1D arrays for better cache coherency. It also supports both synchronous and asynchronous send/receive modes, tunable via a preprocessor macro, so you can take a look.

With that being said, let's turn to a few improvable parts in your code:

1. You don't have to use an array of pointers to represent the matrix format. You can use an 1D array with row-major/column-major format. If you use a "2D array" (array of pointers), each row/colum slice of your matrix are not guaranteed to be stored in consecutive locations.

2. Isend/Irecv lead to degraded performance. You can verify this via benchmarking. MPI sends/receives extra messages when asynchronous mode is used, which makes basic operations more expensive. And since you have blocking anyways somewhere in your code, ISend/IRecv do not make a huge difference.

3. It is common for worker nodes to send back processed slices of data to be merged and outputted in the master node. You can technically just invoke printf in a loop on each worker node to output its own slice of data (since mpirun redirects all stdout streams on each node to your terminal screen), but MPI does not guarantee the order of output.

4. If you want a more elegant semantics, you should checkout MPI_Scatterv/MPI_Gatherv scatter/gather functions: https://www.mpich.org/static/docs/v3.1/www3/MPI_Scatterv.html. (But I don't think you'll get significant performance gains by using them)

• Thank you for your comments bearaqua. I'll have a look at the code you posted and see if putting my matrix into a 1D array makes it any faster. From my own benchmarking, the non-blocking versions of send and recv are faster though, do you have any references that go into more detail on when the non-blocking versions are slower? I think it could be faster in my case because of the massive number of messages this sends (size of the matrix is 10000x10000 scattered over 16 nodes). – nick Oct 31 '18 at 14:12
• @nickI believe when I benchmarked my matrix multiplication code earlier this year, the nonblocking version was slower (on 20+ nodes each with 24 cores). Feel free to re-run the test though. I think the difference might be the actual organization of code, i.e. where you put the Send/Recv statements. I've seen this mentioned elsewhere; I'm trying to find a reference for you. – BearAqua Oct 31 '18 at 20:15
• @nick I finally hovered across a useful link today! Sorry for posting this comment so late, but this link might be of help: mcs.anl.gov/research/projects/mpi/sendmode.html – BearAqua Feb 18 at 1:36
• Quote: The best performance is likely if you can write your program so that you could use just MPI_Ssend; in that case, an MPI implementation can completely avoid buffering data. Use MPI_Send instead; this allows the MPI implementation the maximum flexibility in choosing how to deliver your data. (Unfortunately, one vendor has chosen to have MPI_Send emphasize buffering over performance; on that system, MPI_Ssend may perform better.) If nonblocking routines are necessary, then try to use MPI_Isend or MPI_Irecv. Use MPI_Bsend only when it is too inconvienent to use MPI_Isend. – BearAqua Feb 18 at 1:36
• Thanks for the reference! This is very helpful, I'll have to try out that Ssend function, wasn't even aware it existed. – nick Feb 18 at 22:22