如何在MPI中将矩阵从一个进程转移到另一个进程

这是一个每周至少弹出一次的主题。我通常会重复回答你的问题，但由于你使用的是std::vector，而不是原始指针，我想你应该得到一个更详细的答案，涉及MPI的一个鲜为人知但非常强大的功能，即它的数据类型系统。

首先，std::vector<std::vector<T>>不是一个连续类型。想想它是如何在记忆中排列的。std::vector<T>本身通常实现为一个结构，该结构包含指向堆上分配的数组的指针和一组记账信息，例如数组容量及其当前大小。向量的向量就是一个结构，它包含一个指向结构堆数组的指针，每个结构都包含一个指针指向另一个堆数组：

p1 [ data ] ---> p1[0] [ data ] ---> [ p1[0][0] | p1[0][1] | ... ]
[ size ]            [ size ]
[ cap. ]            [ cap. ]
p1[1] [ data ] ---> [ p1[1][0] | p1[1][1] | ... ]
[ size ]
[ cap. ]
...

这是两个级别的指针间接寻址，只需要访问给定行的数据。编写p1[i][j]时，编译器读取std::vector<T>::operator[]()的代码两次，最后生成指针算术和解引用表达式，该表达式给出特定矩阵元素的地址。

MPI不是编译器扩展。它也不是某种深奥的模板库。它对C++容器对象的内部结构一无所知。它只是一个通信库，只提供了一个在C和Fortran中都能工作的抽象级别。在构思MPI的时候，Fortran甚至没有对用户定义的聚合类型(C/C++中的结构)的主流支持，因此MPIneneneba API在很大程度上是以数组为中心的。也就是说，这并不意味着MPI只能发送数组。相反，它有一个非常复杂的类型系统，如果你愿意投入额外的时间和代码来发送任意形状的对象，它允许你发送。让我们来看看可能的不同方法。

MPI将愉快地从连续内存区域发送数据或在连续内存区域中接收数据。将非连续内存布局转换为连续内存布局并不困难。不是NxN形状的std::vector<std::vector<T>>，而是创建一个大小为N²的平面std::vector<T>，然后在平面中构建一个辅助指针阵列：

vector<int> mat_data();
vector<int *> mat;
mat_data.resize(N*N);
for (int i = 0; i < N; i++)
mat.push_back(&mat_data[0] + i*N);

您可能希望将其封装在一个新的Matrix2D类中。有了这种安排，您仍然可以使用mat[i][j]来引用矩阵元素，但现在所有的行都整齐地排列在内存中。如果你想发送这样一个对象，你只需调用：

MPI_Send(mat[0], N*N, MPI_INT, ...);

如果你已经在接收器端分配了一个NxN矩阵，只需执行：

MPI_Recv(mat[0], N*N, MPI_INT, ...);

如果您还没有分配矩阵，并且希望能够接收任意大小的平方矩阵，请执行：

MPI_Status status;
// Probe for a message
MPI_Probe(..., &status);
// Get message size in number of integers
int nelems;
MPI_Get_count(&status, MPI_INT, &nelems);
N = sqrt(nelems);
// Allocate an NxN matrix mat as show above
// Receive the message
MPI_Recv(mat[0], N*N, MPI_INT, status.MPI_SOURCE, status.MPI_TAG, ...);

不幸的是，并不总是可以简单地将vector<vector<T>>交换为平面阵列类型，尤其是当您调用无法控制的外部库时。在这种情况下，您还有两个选择。

当矩阵很小时，手动打包和解包数据以进行通信并非不可行：

std::vector<int> p1_flat;
p1_flat.reserve(p1.size() * p1.size());
for (auto const &row : p1)
std::copy(row.begin(), row.end(), std::back_inserter(p1_flat));
MPI_Send(&p1_flat[0], ...);

在接收方，你会做相反的事情。

当矩阵很大时，打包和拆包就变成了耗费时间和内存的活动。幸运的是，MPI有一些条款允许您跳过该部分，让它为您打包。如前所述，由于MPI只是一个简单的通信库，它不能自动理解语言类型，并且它使用MPI数据类型形式的提示来正确处理底层语言类型。MPI数据类型类似于一个配方，它告诉MPI在哪里以及如何访问内存中的数据。它是形式为(offset, primitive type):的元组的集合

• offset告诉MPI相应的数据段相对于给MPI_Send()等函数的地址的位置
• 基元类型告诉MPI在该特定偏移处的基元语言数据类型是什么

最简单的MPI数据类型是与语言标量类型相对应的预定义数据类型。例如，MPI_INT是元组(0, int)的底层，它告诉MPI将直接位于所提供的地址处的存储器视为int的实例。当你告诉MPI你实际上正在发送一个完整的MPI_INT数组时，它知道它需要从缓冲区位置获取一个元素，然后在内存中前进到int的大小，再获取另一个，以此类推。C++的数据串行化库的工作方式不太可能。就像您可以在C++中从更简单的数据类型构建聚合类型一样，MPI允许您从更简单数据类型构建复杂的数据类型。例如，[(0, int), (16, float)]数据类型告诉MPI从缓冲区地址取int，从超过缓冲区地址的16个字节取float。

有两种方法可以构造数据类型。您可以创建一个更简单类型的数组，重复某个访问模式(这允许您在该模式中指定统一的间隙)，也可以创建任意更简单数据类型的结构。你需要的是后者。您需要能够告诉MPI以下内容："听着。我有那些N数组要发送/接收，但它们不可预测地分散在堆中。这是它们的地址。请将它们连接起来，并将它们作为一条消息发送/接收。"然后通过使用MPI_Type_create_struct构建结构数据类型来告诉它这一点。

struct-datatype构造函数接受四个输入参数：

• int count-新数据类型中的块(组件)数量，在您的情况下为p.size()(p是vector<vector<int>>之一)
• int array_of_blocklengths[]-同样，由于MPI的数组性质，结构化数据类型实际上是更简单数据类型的数组(块)的结构；在这里，您必须指定一个数组，其中的元素设置为相应行的大小
• MPI_Aint array_of_displacements[]——对于每个块，MPI需要知道它相对于数据缓冲区地址的位置；这既可以是正的也可以是负的，并且最简单的方法是在这里简单地传递所有数组的地址
• MPI_Datatype array_of_types[]——结构的每个块中的数据类型；您需要传递一个元素设置为MPI_INT的数组

在代码中：

// Block lengths
vector<int> block_lengths;
// Block displacements
vector<MPI_Aint> block_displacements;
// Block datatypes
vector<MPI_Datatype> block_dtypes(p.size(), MPI_INT);
for (auto const &row : p) {
block_lengths.push_back(row.size());
block_displacements.push_back(static_cast<MPI_Aint>(&row[0]));
}
// Create the datatype
MPI_Datatype my_matrix_type;
MPI_Type_create_struct(p.size(), block_lengths, block_displacements, block_dtypes, &my_matrix_type);
// Commit the datatatype to make it usable for communication
MPI_Type_commit(&my_matrix_type);

最后一步告诉MPI新创建的数据类型将用于通信。如果这只是构建更复杂数据类型的中间步骤，那么提交步骤可以省略。

我们现在可以使用my_matrix_type发送p:中的数据

MPI_Send(MPI_BOTTOM, 1, my_matrix_type, ...);

MPI_BOTTOM到底是什么？这是地址空间的底部，在许多平台上基本上是0。在大多数系统上，它与NULL或nullptr相同，但没有指向任何地方的指针的语义。我们在这里使用MPI_BOTTOM，因为在上一步中，我们使用了每个数组的地址作为相应块的偏移量。我们可以减去第一行的地址：

for (auto const &row : p) {
block_lengths.push_back(row.size());
block_displacements.push_back(static_cast<MPI_Aint>(&row[0] - &p[0][0]));
}

然后，我们使用以下方式发送p：

MPI_Send(&p[0][0], 1, my_matrix_type, ...);

请注意，您只能使用此数据类型来发送p的内容，而不能发送其他vector<vector<int>>实例的内容，因为偏移量会有所不同。使用绝对地址还是从第一行的地址偏移量创建my_matrix_type并不重要。因此，MPI数据类型的生存期应该与p本身的生存期相同。

不再需要时，应释放my_matrix_type：

MPI_Type_free(&my_matrix_type);

这同样适用于在vector<vector<T>>中接收数据。首先需要调整外部向量的大小，然后调整内部向量的大小以准备内存。然后构建MPI数据类型并使用它来接收数据。如果不再重用同一缓冲区，请释放MPI数据类型。

您可以将上述所有步骤整齐地打包到一个支持MPI的2D矩阵类中，该类可以释放类析构函数中的MPI数据类型。它还将确保为每个矩阵构建一个单独的MPI数据类型。

与第一种方法相比，这有多快？它比简单地使用平面阵列慢一点，而且可能比打包和拆包慢或快。这肯定比将每一行作为单独的消息发送要快。此外，一些网络适配器支持聚集读取和分散写入，这意味着MPI库只需将MPI数据类型中的偏移量直接传递给硬件，而硬件将承担将分散数组组装成单个消息的重任。这可以在沟通渠道的两侧非常有效地完成。

请注意，您不必在发送方和接收方都执行相同的操作。在发送方使用用户定义的MPI数据类型，在接收方使用简单的平面阵列，这是非常好的。反之亦然。MPI不在乎，只要发送器总共发送N²MPI_INTs，并且接收器期望N²MPI_INTs的整数倍，即，发送和接收类型是否一致。

注意：MPI数据类型是相当可移植的，可以在许多平台上工作，甚至可以在异构环境中进行通信。但它们的构造可能比看起来更棘手。例如，块偏移是类型MPI_Aint，这是指针大小的带符号整数，这意味着它可以用于在给定以地址空间的中间为中心的基的情况下可靠地对整个存储器进行寻址。但它不能代表相隔一半以上内存大小的地址之间的差异。这在大多数将虚拟地址空间1:1拆分为用户和内核部分的操作系统上都不是问题，其中包括x86上的32位Linux、x86上没有4 GB调整的32位Windows、x86和ARM上的64位版本的Linux、Windows和macOS，以及大多数其他32位和64位体系结构。但是，有些系统要么完全分离用户和内核地址空间，32位macOS就是一个显著的例子，要么可以进行1:1以外的拆分，32位Windows就是一个例子，它具有4 GB的调优功能，可以进行3:1的拆分。在这样的系统上，不应该使用具有块偏移的绝对地址的MPI_BOTTOM，也不应该使用来自第一行的相对偏移。相反，应该派生一个指向地址空间中间的指针，并计算其偏移量，然后将该指针用作MPI通信原语中的缓冲区地址。

免责声明：这是一篇很长的帖子，可能有一些错误在我的雷达下。期待编辑。此外，我声称没有能力编写惯用的C++。