跳到主要内容

丁致宇第五周学习报告

[TOC]

矩阵的知识

矩阵与向量相乘

矩阵与向量的乘法是线性代数中的基本操作之一,它遵循特定的规则。在数学中,一个矩阵可以被看作是一个线性变换,而向量则可以被看作是空间中的一个点或者箭头。当我们将一个矩阵与一个向量相乘时,我们实际上是在将这个线性变换应用到这个向量上。

矩阵与向量乘法的规则

假设有一个 m×nm \times n 的矩阵 AA 和一个 nn 维的列向量 x\mathbf{x},它们的乘积是一个 mm 维的列向量 y\mathbf{y}。这个乘积定义为:

Ax=yA\mathbf{x} = \mathbf{y}

其中 AA 是这样一个矩阵:

A=[a11a12a1na21a22a2nam1am2amn]A = \begin{bmatrix} a_{11} & a_{12} & \cdots & a_{1n} \\ a_{21} & a_{22} & \cdots & a_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m1} & a_{m2} & \cdots & a_{mn} \end{bmatrix}

x\mathbf{x} 是这样一个向量:

x=[x1x2xn]\mathbf{x} = \begin{bmatrix} x_{1} \\ x_{2} \\ \vdots \\ x_{n} \end{bmatrix}

那么乘积 y\mathbf{y} 的每个元素 yiy_iii 从 1 到 mm)可以通过下面的方式计算:

yi=ai1x1+ai2x2++ainxny_i = a_{i1}x_1 + a_{i2}x_2 + \cdots + a_{in}x_n

或者更紧凑的表示为:

yi=j=1naijxjy_i = \sum_{j=1}^{n} a_{ij}x_j

例子

假设我们有一个 2×32 \times 3 的矩阵 AA 和一个 33 维的列向量 x\mathbf{x}

A=[201134],x=[351]A = \begin{bmatrix} 2 & 0 & 1 \\ -1 & 3 & 4 \end{bmatrix}, \quad \mathbf{x} = \begin{bmatrix} 3 \\ 5 \\ 1 \end{bmatrix}

那么它们的乘积 y\mathbf{y} 将是:

y=Ax=[201134][351]=[(23)+(05)+(11)(13)+(35)+(41)]=[714]\mathbf{y} = A\mathbf{x} = \begin{bmatrix} 2 & 0 & 1 \\ -1 & 3 & 4 \end{bmatrix} \begin{bmatrix} 3 \\ 5 \\ 1 \end{bmatrix} = \begin{bmatrix} (2 \cdot 3) + (0 \cdot 5) + (1 \cdot 1) \\ (-1 \cdot 3) + (3 \cdot 5) + (4 \cdot 1) \end{bmatrix} = \begin{bmatrix} 7 \\ 14 \end{bmatrix}

这里我们分别计算了 y1y_1y2y_2 的值:

  • 对于 y1y_1,我们计算了 23+05+11=72 \cdot 3 + 0 \cdot 5 + 1 \cdot 1 = 7
  • 对于 y2y_2,我们计算了 13+35+41=14-1 \cdot 3 + 3 \cdot 5 + 4 \cdot 1 = 14

注意事项

  • 矩阵 AA 的列数必须与向量 x\mathbf{x} 的行数相等。
  • 结果 y\mathbf{y} 的维数将与矩阵 AA 的行数相同。
  • 矩阵与向量的乘法不是交换的,即 AxxAA\mathbf{x} \neq \mathbf{x}A
  • 矩阵与向量的乘法是分配的和结合的,即 A(x+y)=Ax+AyA(\mathbf{x} + \mathbf{y}) = A\mathbf{x} + A\mathbf{y},且 (AB)x=A(Bx)(AB)\mathbf{x} = A(B\mathbf{x}),其中 AABB 是矩阵,x\mathbf{x}y\mathbf{y} 是向量。

矩阵与向量的乘法在许多领域都非常重要,包括计算机图形学、工程、物理学、统计学以及机器学习和数据科学等。

矩阵乘法

就是前面的矩阵的行乘以后面矩阵的列,每一行都与每一列相乘

结果的行数由前面矩阵决定,结果所在的列数由后面的矩阵决定

矩阵乘法是线性代数中的一个核心概念,用于结合两个矩阵的信息。假设我们有两个矩阵 AABB,它们的乘积是第三个矩阵 CC。矩阵乘法的定义如下:

  • 矩阵 AA 的大小为 m×nm \times n,即它有 mm 行和 nn 列。
  • 矩阵 BB 的大小为 n×pn \times p,即它有 nn 行和 pp 列。
  • 为了能够进行乘法,矩阵 AA 的列数必须等于矩阵 BB 的行数。

如果上述条件满足,那么矩阵 AA 和矩阵 BB 的乘积是一个 m×pm \times p 的矩阵 CC。矩阵 CC 中的每个元素 cijc_{ij} 是通过将矩阵 AA 的第 ii 行与矩阵 BB 的第 jj 列对应元素相乘然后求和得到的:

cij=ai1b1j+ai2b2j++ainbnj=k=1naikbkjc_{ij} = a_{i1}b_{1j} + a_{i2}b_{2j} + \cdots + a_{in}b_{nj} = \sum_{k=1}^{n} a_{ik}b_{kj}

其中 i=1,2,,mi = 1, 2, \ldots, mj=1,2,,pj = 1, 2, \ldots, p

例子

让我们通过一个具体的例子来说明矩阵乘法是如何进行的:

假设我们有如下两个矩阵 AABB

A=[2013],B=[1421]A = \begin{bmatrix} 2 & 0 \\ -1 & 3 \end{bmatrix}, \quad B = \begin{bmatrix} 1 & 4 \\ 2 & -1 \end{bmatrix}

矩阵 AA 是一个 2×22 \times 2 的矩阵,矩阵 BB 也是一个 2×22 \times 2 的矩阵。它们的乘积 CC 将是:

C=AB=[(21+02)(24+01)(11+32)(14+31)]=[2857]C = AB = \begin{bmatrix} (2 \cdot 1 + 0 \cdot 2) & (2 \cdot 4 + 0 \cdot -1) \\ (-1 \cdot 1 + 3 \cdot 2) & (-1 \cdot 4 + 3 \cdot -1) \end{bmatrix} = \begin{bmatrix} 2 & 8 \\ 5 & -7 \end{bmatrix}

这里,CC 的每个元素是通过以下计算得到的:

  • c11=(21)+(02)=2c_{11} = (2 \cdot 1) + (0 \cdot 2) = 2
  • c12=(24)+(01)=8c_{12} = (2 \cdot 4) + (0 \cdot -1) = 8
  • c21=(11)+(32)=5c_{21} = (-1 \cdot 1) + (3 \cdot 2) = 5
  • c22=(14)+(31)=7c_{22} = (-1 \cdot 4) + (3 \cdot -1) = -7

注意事项

  • 矩阵乘法不是交换的,即 ABBAAB \neq BA 一般来说。
  • 矩阵乘法是结合的,即 (AB)C=A(BC)(AB)C = A(BC)
  • 矩阵乘法是分配的,即 A(B+C)=AB+ACA(B + C) = AB + AC
  • 矩阵乘法通常包含大量的乘法和加法操作,因此计算上可能相当耗时,特别是对于大型矩阵。

矩阵乘法在数据处理、物理科学、工程、计算机图形学以及经济学等许多领域都有广泛的应用。

现在,我们用C语言来实现矩阵乘法。以下是一个简单的程序,它实现了两个矩阵的乘法:

#include <stdio.h>

#define MAX_SIZE 100

void matrixMultiply(int m, int n, int p, int A[][MAX_SIZE], int B[][MAX_SIZE], int C[][MAX_SIZE]) {
for (int i = 0; i < m; i++) {
for (int j = 0; j < p; j++) {
C[i][j] = 0; // 初始化结果矩阵的当前元素为0
for (int k = 0; k < n; k++) {
C[i][j] += A[i][k] * B[k][j]; // 计算结果矩阵的当前元素
}
}
}
}

int main() {
int m, n, p;
int A[MAX_SIZE][MAX_SIZE], B[MAX_SIZE][MAX_SIZE], C[MAX_SIZE][MAX_SIZE];

// 假设用户将输入矩阵的大小和元素
printf("Enter rows and columns for matrix A: ");
scanf("%d %d", &m, &n);

printf("Enter elements of matrix A:\n");
for (int i = 0; i < m; i++)
for (int j = 0; j < n; j++)
scanf("%d", &A[i][j]);

printf("Enter rows and columns for matrix B: ");
scanf("%d %d", &n, &p); // 注意:这里n应该与之前输入的n相同

printf("Enter elements of matrix B:\n");
for (int i = 0; i < n; i++)
for (int j = 0; j < p; j++)
scanf("%d", &B[i][j]);

// 执行矩阵乘法
matrixMultiply(m, n, p, A, B, C);

// 打印结果矩阵
printf("Result of matrix multiplication:\n");
for (int i = 0; i < m; i++) {
for (int j = 0; j < p; j++)
printf("%d ", C[i][j]);
printf("\n");
}

return 0;
}

这个程序首先定义了一个matrixMultiply函数,它接受两个矩阵的尺寸和元素,然后计算它们的乘积。main函数用于获取用户输入的矩阵尺寸和元素,并调用matrixMultiply函数来计算乘积,最后打印结果矩阵。

请注意,这个程序假设用户会输入有效的矩阵尺寸,其中矩阵A的列数等于矩阵B的行数,且所有矩阵的尺寸都不超过MAX_SIZE定义的大小。在实际应用中,你可能需要添加额外的错误检查来确保输入的有效性。

实践:用MPI并行化矩阵乘运算

代码在Matrix_MPI文件夹内

第一次矩阵并行化出现的问题

代码中存在的一些问题如下:

  1. 未处理矩阵大小不均匀分配给进程的情况。在A_local_row = m / size;计算中,如果m不能被size整除,将会导致矩阵行数分配不均。

  2. MPI_Gather调用可能不正确。如果m不能被size整除,那么最后一个进程可能拥有不同数量的行,需要使用MPI_Gatherv来处理不同的接收数量。

  3. MPI_ScatterMPI_Gather使用了AC数组,但是非零进程中的这些数组是未初始化的。

  4. MPI_Reduce调用的max_elapsed_time只在rank == 0的代码块中声明,这将导致编译错误,因为在其他进程中它是未声明的。

n体问题

简介

在并行编程和计算物理学中,n体问题通常指的是模拟和计算一个由n个相互作用的粒子组成的系统的动力学。在天体物理学中,这些粒子可以是星星、行星或其他天体,它们通过万有引力相互作用;在分子动力学中,粒子可以是原子或分子,它们通过电磁力相互作用。n体问题的目标是确定系统中所有粒子随时间的运动。

n体问题是一个经典的物理问题,因为它涉及到非线性的多体相互作用,使得问题的解析解通常是不可能的,除了非常简单的情况(如两体问题)。因此,科学家们通常使用数值方法来近似求解n体问题,这涉及到通过离散时间步骤来迭代计算粒子的位置和速度。

并行编程在解决n体问题中非常重要,因为:

  1. 计算量大:n体问题的计算复杂度随着粒子数目的增加而显著增长。每个时间步长,每个粒子都需要与其他所有粒子计算相互作用力,这会导致计算量随着粒子数目的增加而呈平方增长。

  2. 可分割性:n体问题的计算可以很自然地分割成多个任务,每个任务计算一部分粒子的相互作用力,因此非常适合并行处理。

  3. 实时性要求:在某些应用中,如视频游戏或现实世界模拟中,对于n体问题的解需要实时或接近实时计算,而并行处理可以提供足够的计算资源来满足这些要求。

并行解决n体问题通常涉及以下步骤:

  • 分解任务:将整个问题分解为可以并行处理的较小任务。
  • 计算相互作用:并行计算每对粒子之间的相互作用力。
  • 整合结果:将并行计算的结果合并,更新每个粒子的位置和速度。
  • 时间推进:将系统向前推进一个时间步长,并重复上述过程。

实现并行计算的方法包括使用多线程、多处理器、多核心CPU以及图形处理单元(GPU)等技术。特别是GPU由于其大量的并行处理单元,非常适合处理这种类型的计算密集型任务。

在编程实践中,解决n体问题的并行算法需要仔细设计以最小化进程或线程间的通信开销,并最大化计算与通信的重叠,以提高并行效率。

n体问题中的基本物理学知识

  1. 万有引力定律

    牛顿的万有引力定律描述了两个物体之间的引力,其大小与两个物体的质量成正比,与它们之间距离的平方成反比。公式如下:

    F=Gm1m2r2F = G \frac{m_1 m_2}{r^2}

    其中:

    • FF 是两个物体之间的引力。
    • GG 是万有引力常数,大约为 6.67430×1011m3kg1s26.67430 \times 10^{-11} \text{m}^3\text{kg}^{-1}\text{s}^{-2}
    • m1m_1m2m_2 是两个物体的质量。
    • rr 是两个物体之间的距离。
  2. 牛顿的运动定律

    牛顿的第二运动定律描述了力和物体运动状态变化之间的关系,即力等于质量乘以加速度:

    F=maF = m a

    其中:

    • FF 是作用在物体上的合力。
    • mm 是物体的质量。
    • aa 是物体的加速度。
  3. 计算引力

    代码中的 gravitational_force 函数实现了万有引力定律的计算。首先,计算出两个物体之间的距离:

    distance=(x2x1)2+(y2y1)2+(z2z1)2\text{distance} = \sqrt{(x_2 - x_1)^2 + (y_2 - y_1)^2 + (z_2 - z_1)^2}

    然后根据万有引力定律计算引力的大小,最后将这个力按照三个坐标轴方向分解,得到向量形式的力。

  4. 更新速度和位置

    update_bodies 函数中,首先计算出每个物体受到的合力,然后根据牛顿的第二运动定律更新物体的速度:

    vnew=vold+FmΔtv_{\text{new}} = v_{\text{old}} + \frac{F}{m} \Delta t

    其中 Δt\Delta t 是时间步长。

    最后,更新物体的位置,假设在很小的时间步长 Δt\Delta t 内速度是恒定的:

    xnew=xold+vΔtx_{\text{new}} = x_{\text{old}} + v \Delta t

    这个简单的更新方法称为欧拉方法,它在数值上可能不是非常稳定或精确,尤其是在较大的时间步长或复杂的动力学系统中。对于更复杂或需要更高精度的系统,可能需要使用更高级的积分方法,如Runge-Kutta方法。

  5. 模拟循环main 函数中,有一个模拟循环,它重复调用 update_bodies 函数来模拟物体随时间的运动。每次循环都代表时间向前推进了 Δt\Delta t

注意事项

在实际的数值模拟中,还需要考虑诸如数值稳定性和能量守恒等因素。例如,当物体非常接近时,上述方法可能会导致数值不稳定,因为引力会变得非常大。此外,欧拉方法在长时间积分时可能会导致能量逐渐增加或减少,这在物理上是不正确的。为了解决这些问题,可以采用更高阶的积分方法,或者引入软化长度以避免在物体非常接近时引力变得无限大。

实践:使用MPI并行化n体问题

代码在N_Body_Problem文件夹内

知识点

MPI_Allgather的特殊用法(MPI_IN_PLACE)

在MPI中,MPI_Allgather函数是一个集合通信操作,它用于在所有参与的进程之间收集数据,并将收集到的数据分发给所有进程。这个函数通常用于当每个进程都有一个数据块,并且想要每个进程都得到所有其他进程的数据块时。

函数原型如下:

int MPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm)

参数解释:

  • sendbuf: 指向发送缓冲区的起始地址,即包含该进程要发送数据的内存地址。
  • sendcount: 发送数据的元素数量。
  • sendtype: 发送数据的MPI数据类型。
  • recvbuf: 指向接收缓冲区的起始地址,即用于存储所有进程发送数据的内存地址。
  • recvcount: 每个进程接收数据的元素数量。
  • recvtype: 接收数据的MPI数据类型。
  • comm: MPI通信器,通常是MPI_COMM_WORLD,它包含了所有的进程。

在您提供的代码段中:

MPI_Allgather(MPI_IN_PLACE, local_n * sizeof(Vector3D), MPI_BYTE,
positions, local_n * sizeof(Vector3D), MPI_BYTE, MPI_COMM_WORLD);

参数解释如下:

  • MPI_IN_PLACE: 这是一个特殊的参数,用于告诉MPI函数在原地执行操作。当使用MPI_IN_PLACE时,输入数据(即发送缓冲区)和输出数据(即接收缓冲区)使用同一个缓冲区。这意味着每个进程的发送缓冲区就是它在接收缓冲区中的相应位置。使用MPI_IN_PLACE可以避免不必要的数据复制,提高效率。
  • local_n * sizeof(Vector3D): 这指定了每个进程要发送数据的大小,即每个进程的局部数组local_positionslocal_velocities包含的字节总数。这里local_n是每个进程负责的Vector3D元素的数量,sizeof(Vector3D)是一个Vector3D结构的大小(以字节为单位)。
  • MPI_BYTE: 这是一个基本的MPI数据类型,表示数据将以字节为单位发送和接收。这意味着MPI将不会解释数据的实际结构,只是简单地传输字节。
  • positions: 这是接收缓冲区,所有进程的数据将被收集并存储在这个数组中。
  • local_n * sizeof(Vector3D): 这指定了每个进程接收数据的大小,也就是每个进程在positions数组中的分段大小。
  • MPI_COMM_WORLD: 这是一个通信器,包含了所有的MPI进程。

所以,这行代码的作用是:每个进程将其负责的local_nVector3D结构的数据发送出去,并从其他进程收集同样数量和类型的数据,将结果存储在positions数组中。因为使用了MPI_IN_PLACE,所以发送和接收都是在positions数组上进行的,不需要额外的发送缓冲区。这样,每个进程在执行完这个操作之后,都将拥有整个系统的完整位置信息。

那么将MPI_IN_PLACE放在接收缓冲区是正确的吗

不,这样的代码是不正确的。在MPI_Allgather调用中,当使用MPI_IN_PLACE时,它应该只出现在发送缓冲区(sendbuf)的位置,并且接收缓冲区(recvbuf)应该是一个不同的缓冲区,用于接收所有进程的数据。MPI_IN_PLACE用于指示MPI库,发送数据就在接收缓冲区的相应位置,因此不需要提供一个单独的发送缓冲区。

正确的MPI_Allgather调用,使用MPI_IN_PLACE,应该看起来像这样:

MPI_Allgather(MPI_IN_PLACE, 0, MPI_DATATYPE_NULL,
positions, local_n * sizeof(Vector3D), MPI_BYTE, MPI_COMM_WORLD);

这里,MPI_IN_PLACE告诉MPI库在原地收集数据,即直接在positions数组中操作。由于我们在原地操作,发送数据的数量(sendcount)和发送数据的类型(sendtype)都不会被使用,所以可以将sendcount设置为0,将sendtype设置为MPI_DATATYPE_NULL

在你提供的代码中:

MPI_Allgather(positions, local_n * sizeof(Vector3D), MPI_BYTE,
MPI_IN_PLACE, local_n * sizeof(Vector3D), MPI_BYTE, MPI_COMM_WORLD);

**如果你尝试将MPI_IN_PLACE用在接收缓冲区的位置,这会导致一个错误,因为MPI_IN_PLACE是专门设计用于发送缓冲区的,以告知MPI库发送和接收使用相同的缓冲区。**此外,发送缓冲区(sendbuf)和接收缓冲区(recvbuf)的位置不能互换。正确的用法是将MPI_IN_PLACE放在sendbuf的位置,并提供一个适当的接收缓冲区。

多维数组传参

第一维可以省略,其他维度不能省略,因为编译器要根据每一维的元素数量在一维数组上进行划分

在C语言中,当你在函数参数中传递一个多维数组时,除了第一维可以省略外,其余的维度大小是必须指定的。这是因为编译器需要这些信息来计算数组元素的位置。

对于一维数组,当你将其作为参数传递给函数时,实际上传递的是数组的首元素指针。因此,一维数组的长度在函数参数中通常是可选的,例如 int vector[]int *vector

但是,对于多维数组,除了第一维之外,你必须指定数组的其他维度,以便编译器能够正确地计算偏移量。例如,如果你有一个二维数组 int matrix[3][3],并且你想将它传递给一个函数,你的函数参数应该指定第二维的大小,如 int matrix[][3]

在我之前给出的例子中,N 是一个预处理器宏,它定义了数组的大小。因此,当你看到 int matrix[N][N] 作为函数参数时,这意味着 N 必须在函数被调用之前就已经定义了。

这是一个修改后的 matrixVectorMultiply 函数,它使用了正确的参数声明方式:

void matrixVectorMultiply(int matrix[][N], int vector[N], int result[N]);

这里,我们指定了 matrix 的第二维大小为 N,而 vectorresult 作为一维数组,它们的大小是可选的。但是为了代码的清晰性和可读性,我在之前的例子中包含了大小 N。如果你在函数的实现中不需要知道一维数组的大小,你可以选择省略它,如下:

void matrixVectorMultiply(int matrix[][N], int vector[], int result[]);

在这种情况下,你需要确保在函数调用时传递正确大小的数组。

三种点对点通信方式的选用

选择使用阻塞、缓冲或非阻塞通信的最佳时机通常取决于应用程序的特定需求、通信模式和性能目标。下面是一些指导原则帮助你决定:

阻塞通信 (MPI_Send, MPI_Recv)

  • 简单性:如果你的程序逻辑简单,不需要同时进行计算和通信,标准的阻塞通信可能是最简单的选择。
  • 确定性:当你需要确保在执行后续代码之前消息已经被发送或接收时,阻塞通信提供了这种确定性。
  • 小消息:对于小消息,阻塞通信的开销可能可以忽略不计,因为小消息通常很快就能被发送出去或者接收。

缓冲发送 (MPI_Bsend)

  • 可用缓冲区:如果你的系统有足够的缓冲区资源,并且你希望避免发送操作可能的阻塞,缓冲发送可以是一个好的选择。
  • 中等大小的消息:对于中等大小的消息,使用缓冲发送可以减少发送操作的阻塞时间,因为数据会被复制到缓冲区中。
  • 计算与通信重叠:如果你希望在消息发送的同时执行一些计算,缓冲发送可以提供这种重叠的可能性,尽管它不如非阻塞通信灵活。

非阻塞通信 (MPI_Isend, MPI_Irecv)

  • 性能:当你需要最大限度地提高程序性能,特别是在需要计算和通信重叠的情况下,非阻塞通信通常是首选。
  • 大消息:对于大消息,非阻塞通信允许发送操作在数据传输的同时进行其他计算,从而提高资源利用率。
  • 复杂的通信模式:在具有复杂通信模式的程序中,非阻塞通信可以提供更好的控制,因为它允许同时启动多个通信操作,并在它们完成时进行处理。
  • 流水线操作:如果你的应用程序可以分为多个可以并行处理的阶段,非阻塞通信可以帮助你设置流水线,其中计算和通信可以在不同阶段并行执行。

总结

  • 如果你的应用程序通信模式简单,或者你刚开始使用MPI,那么从标准的阻塞通信开始是合理的。
  • 如果你的应用程序需要在通信时做一些计算,而且你不想处理非阻塞通信的复杂性,那么缓冲发送可能是一个好的中间选择。
  • 如果你需要最大化性能,尤其是在有大量并发通信和计算的情况下,那么非阻塞通信是最好的选择。

在任何情况下,最好的方法是通过实验和性能分析来确定哪种通信方式最适合你的应用程序。不同的硬件和网络架构也可能影响最佳选择。

缓冲发送MPI_Bsend函数

MPI_Bsend 函数是一个缓冲发送函数,它的原型如下:

int MPI_Bsend(const void *buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm)

下面是每个参数的详细解释:

  • const void *buf: 这是指向发送消息的起始位置的指针。这个缓冲区包含了要发送的数据。

  • int count: 这个参数指定了要发送的数据元素的数量。结合数据类型(datatype),这个参数定义了发送缓冲区中数据的总大小。

  • MPI_Datatype datatype: 这个参数指定了每个数据元素的数据类型。MPI预定义了一系列的数据类型,比如MPI_INT用于整数,MPI_FLOAT用于浮点数等。

  • int dest: 这个参数指定了消息的目的地进程的rank(标识符)。在MPI中,每个进程都有一个唯一的rank,用于标识通信的目标或来源。

  • int tag: 这个参数是一个整数标签,用于区分不同的消息。发送和接收操作中的标签必须匹配,才能正确地配对消息。

  • MPI_Comm comm: 这个参数指定了通信器(communicator),它是一个进程组的上下文,其中包含了可以进行通信的进程集合。最常用的通信器是MPI_COMM_WORLD,它包含了所有的MPI进程。

关于分配缓冲区的方法,你需要在调用MPI_Bsend之前手动分配缓冲区。缓冲区的大小应该足够大,以容纳所有的出站消息加上MPI库可能需要的额外开销。MPI库定义了一个常量MPI_BSEND_OVERHEAD,它表示每个缓冲发送操作可能需要的额外开销。

下面是分配缓冲区和附加缓冲区的示例代码:

int buffer_size = messages_count * (message_size + MPI_BSEND_OVERHEAD);
void* buffer = malloc(buffer_size);

// 附加缓冲区
MPI_Buffer_attach(buffer, buffer_size);

// ... 进行缓冲发送操作 ...

// 分离缓冲区
void* bsend_buff;
int bsend_size;
MPI_Buffer_detach(&bsend_buff, &bsend_size);

// 释放缓冲区内存
free(buffer);

在这个例子中,messages_count 是你打算发送的消息的数量,message_size 是单个消息的大小。这个大小应该根据实际应用程序中要发送的最大消息大小来确定。请注意,MPI_Buffer_attach 函数需要缓冲区的大小(以字节为单位),包括每个消息的MPI_BSEND_OVERHEAD。在不再需要缓冲区时,使用MPI_Buffer_detach来分离缓冲区,并在适当的时候释放缓冲区的内存。

非阻塞通信与缓冲发送在资源使用上的区别

非阻塞通信与缓冲发送(例如在MPI中的MPI_Bsend)在资源使用上有一些区别,但这并不意味着非阻塞通信总是会占用更多的缓冲资源。实际上,它们各自的资源使用取决于多种因素,包括通信库的实现、系统的配置、通信模式等。

缓冲区分为两个,一个是发送缓冲区,一个是用户提前分配的缓冲区,Isend是从发送缓冲区发送,而Bsend是将数据从发送缓冲区复制到用户提前分配的缓冲区,然后进行发送,函数返回后可以立即更改发送缓冲区的数据因为有一个副本在用户提前分配的缓冲区不会影响发送

非阻塞通信

  • 非阻塞通信操作(如MPI_IsendMPI_Irecv)允许程序在等待数据传输完成的同时继续执行其他操作。
  • 非阻塞通信可能需要额外的状态跟踪资源,因为通信库需要记录操作的进展,以便在稍后检查或等待时能够获取状态。
  • 在非阻塞通信中,数据可能会被复制到内部缓冲区,也可能直接从用户提供的缓冲区中发送,这取决于通信库的实现和消息的大小。

缓冲发送

  • 缓冲发送(如MPI_Bsend)需要用户提前分配一个足够大的缓冲区,MPI库使用这个缓冲区来存储即将发送的消息。
  • 一旦消息被复制到这个缓冲区,发送操作就可以返回,而实际的传输可能稍后发生。
  • 缓冲发送的优势在于,一旦消息被复制到缓冲区,发送缓冲区就可以立即重用,而无需等待消息实际被发送到接收方。
  • 缓冲发送可能会占用更多的用户空间内存,因为必须为缓冲区分配足够的空间以存储所有的出站消息。

在某些情况下,非阻塞通信可能比缓冲发送更加高效,因为它可以减少对缓冲区的需求,尤其是当通信库能够利用底层硬件的直接内存访问(DMA)功能进行数据传输时。然而,在其他情况下,如果非阻塞通信需要频繁地复制数据到内部缓冲区,那么它可能会使用更多的缓冲资源。

总的来说,非阻塞通信与缓冲发送各有利弊,它们在资源使用上的差异取决于特定的使用场景和MPI实现的细节。通常,选择哪种通信方式应基于对应用程序的性能要求和资源限制的理解。

阻塞通信示例

首先,我们来看一个使用阻塞通信的简化的矩阵向量乘法示例。在这个例子中,我们假设矩阵已经按行分割并分配给了各个进程。

#include <mpi.h>
#include <stdio.h>
#include <stdlib.h>

int main(int argc, char** argv) {
int rank, size;
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
MPI_Comm_size(MPI_COMM_WORLD, &size);

// 假设n是全局大小,local_n是每个进程的局部大小
int n = 4; // 举例,全局大小为4
int local_n = n / size; // 假设能够整除

// 分配内存
double A[local_n][n], x[n], y[local_n];

// 初始化A和x
// ...

// 分发向量x中的元素到所有进程
for (int i = 0; i < size; ++i) {
if (rank == i) {
// 主进程发送x的部分到所有其他进程
for (int j = 0; j < size; ++j) {
if (j != i) {
MPI_Send(x + i * local_n, local_n, MPI_DOUBLE, j, 0, MPI_COMM_WORLD);
}
}
} else {
// 非主进程接收x的部分
MPI_Recv(x + i * local_n, local_n, MPI_DOUBLE, i, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
}
}

// 计算y的局部部分
for (int i = 0; i < local_n; ++i) {
y[i] = 0;
for (int j = 0; j < n; ++j) {
y[i] += A[i][j] * x[j];
}
}

// 打印结果
for (int i = 0; i < local_n; ++i) {
printf("Process %d: y[%d] = %f\n", rank, i + rank * local_n, y[i]);
}

MPI_Finalize();
return 0;
}

这个例子中,每个进程计算其对应的结果向量的一部分。进程0(主进程)将全局向量x的相应部分发送给其他所有进程,而其他进程接收它们需要的x的部分。然后,每个进程计算自己的结果部分。这里,通信是阻塞的,即MPI_SendMPI_Recv会阻塞直到操作完成。

缓冲发送示例

接下来,我们将相同的逻辑修改为使用缓冲发送。这意味着发送操作将数据复制到MPI的内部缓冲区,然后立即返回,不需要等待接收方实际接收数据。

// ...(与上面相同的初始化和定义)

// 设置缓冲区
int bufsize = n * local_n * sizeof(double) + MPI_BSEND_OVERHEAD;
void* buf = malloc(bufsize);
MPI_Buffer_attach(buf, bufsize);

// 分发向量x中的元素到所有进程
// 同样的循环结构,但是使用MPI_Bsend而不是MPI_Send

// ...(计算和打印结果的代码与上面相同)

// 清理缓冲区
MPI_Buffer_detach(&buf, &bufsize);
free(buf);

// ...

在这个例子中,我们需要首先分配一个缓冲区,并通过MPI_Buffer_attach函数告诉MPI系统这个缓冲区的位置和大小。然后我们使用MPI_Bsend代替MPI_Send。在完成所有的缓冲发送操作后,我们通过MPI_Buffer_detach解除缓冲区的绑定,并释放内存。

非阻塞通信示例

最后,我们来看一个使用非阻塞通信的示例。非阻塞通信允许进程在等待数据发送和接收完成的同时执行其他操作。

// ...(与上面相同的初始化和定义)

// 分发向量x中的元素到所有进程
MPI_Request reqs[size * 2]; // 存储请求
int req_count = 0;

for (int i = 0; i < size; ++i) {
if (rank == i) {
// 非阻塞发送x的部分到所有其他进程
for (int j = 0; j < size; ++j) {
if (j != i) {
MPI_Isend(x + i * local_n, local_n, MPI_DOUBLE, j, 0, MPI_COMM_WORLD, &reqs[req_count++]);
}
}
} else {
// 非阻塞接收x的部分
MPI_Irecv(x + i * local_n, local_n, MPI_DOUBLE, i, 0, MPI_COMM_WORLD, &reqs[req_count++]);
}
}

// 在这里可以执行其他计算或者操作
// ...

// 等待所有非阻塞通信完成
MPI_Waitall(req_count, reqs, MPI_STATUSES_IGNORE);

// ...(计算和打印结果的代码与上面相同)

// ...

在这个例子中,我们使用MPI_IsendMPI_Irecv来启动非阻塞的发送和接收操作。我们存储所有的请求对象,然后使用MPI_Waitall等待所有通信操作完成。这样做可以让进程在通信完成之前执行其他有意义的工作,从而可能提高程序的整体效率。

分析

在这三个示例中,我们可以看到以下区别:

  • 阻塞通信:最简单,适合初学者,但可能不是最高效的,特别是在涉及大量数据交换时。
  • 缓冲发送:减少了发送操作的阻塞时间,但需要额外的内存作为缓冲区,并且需要管理缓冲区的生命周期。
  • 非阻塞通信:提供了最高的灵活性和潜在的性能提升,但代码更复杂,需要管理通信请求和等待操作。

在实际应用中,通信模式的选择应基于对应用程序的通信模式、数据大小和性能要求的理解。通常,非阻塞通信在需要高性能的并行计算应用中更受欢迎,尽管它增加了编程的复杂性。

MPI_Scatterv

MPI_Scatterv 是一个用于在并行计算中分发数据的 MPI (Message Passing Interface) 函数。与 MPI_Scatter 类似,它将一个数组中的数据分发到一组进程中,但与 MPI_Scatter 不同的是,它允许发送不同数量的数据到不同的进程。

MPI_Scatterv 的函数原型如下:

int MPI_Scatterv(
const void *sendbuf, // 根进程中待发送数据的起始地址
const int sendcounts[], // 数组,包含发送到每个进程的数据数量
const int displs[], // 数组,包含每个进程接收的数据在sendbuf中的偏移量
MPI_Datatype sendtype, // 发送数据的类型
void *recvbuf, // 接收数据的起始地址(对于接收进程)
int recvcount, // 接收数据的数量(对于接收进程)
MPI_Datatype recvtype, // 接收数据的类型
int root, // 发送数据的根进程的排名
MPI_Comm comm // 通信器
);

这里是一个使用 MPI_Scatterv 的例子,假设我们有一个根进程,它有一个待发送的整数数组,希望将这个数组的不同部分发送到不同的进程中。每个进程接收的元素数量是不同的。

#include <mpi.h>
#include <stdio.h>
#include <stdlib.h>

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

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

// 根进程的数据
int *sendbuf = NULL;
int sendcounts[size];
int displs[size];

// 每个进程接收的数据缓冲区
int recvbuf[10]; // 假设最大接收数量为10

if (rank == 0) {
// 根进程初始化发送缓冲区
int sendbuf_size = 0;
for (int i = 0; i < size; ++i) {
sendcounts[i] = i + 1; // 第i个进程将接收i+1个元素
sendbuf_size += sendcounts[i];
}

sendbuf = (int*)malloc(sendbuf_size * sizeof(int));

// 填充发送缓冲区
for (int i = 0; i < sendbuf_size; ++i) {
sendbuf[i] = i;
}

// 初始化偏移量数组
displs[0] = 0;
for (int i = 1; i < size; ++i) {
displs[i] = displs[i - 1] + sendcounts[i - 1];
}
}

// 分发数据
MPI_Scatterv(sendbuf, sendcounts, displs, MPI_INT, recvbuf, 10, MPI_INT, 0, MPI_COMM_WORLD);

// 打印接收到的数据
printf("Process %d received:", rank);
for (int i = 0; i < sendcounts[rank]; ++i) {
printf(" %d", recvbuf[i]);
}
printf("\n");

// 根进程需要释放发送缓冲区
if (rank == 0) {
free(sendbuf);
}

MPI_Finalize();
return 0;
}

在这个例子中,根进程(rank 0)有一个整数数组 sendbuf,它想要将这个数组分散给所有进程。每个进程将接收的元素数量由 sendcounts 数组指定,并且 displs 数组指定了每个进程接收的元素在 sendbuf 中的起始位置。每个进程都有一个接收缓冲区 recvbuf,在这个例子中,我们假设每个进程最多接收10个元素,这是一个简化的假设,实际上你会根据实际情况来分配接收缓冲区的大小。

MPI_Gatherv

MPI_Gatherv 是 MPI (Message Passing Interface) 中的一个函数,它用于从一组进程中收集不同数量的数据,并将这些数据聚集到根进程的接收缓冲区中。与 MPI_Gather 相比,MPI_Gatherv 允许每个进程发送不同数量的数据到根进程。

MPI_Gatherv 的函数原型如下:

int MPI_Gatherv(
const void *sendbuf, // 发送数据的起始地址(对于发送进程)
int sendcount, // 发送数据的数量(对于发送进程)
MPI_Datatype sendtype, // 发送数据的类型
void *recvbuf, // 接收数据的起始地址(仅对根进程有效)
const int recvcounts[], // 数组,包含每个进程将发送的数据数量
const int displs[], // 数组,包含每个进程的数据在recvbuf中的偏移量
MPI_Datatype recvtype, // 接收数据的类型(仅对根进程有效)
int root, // 接收数据的根进程的排名
MPI_Comm comm // 通信器
);

下面是一个使用 MPI_Gatherv 的例子,假设我们有一组进程,每个进程都有一个整数数组,它们希望将这个数组中的一部分数据发送到根进程中。

#include <mpi.h>
#include <stdio.h>
#include <stdlib.h>

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

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

// 每个进程的发送缓冲区
int sendbuf[10]; // 假设每个进程发送10个整数
for (int i = 0; i < 10; ++i) {
sendbuf[i] = rank * 10 + i;
}

// 根进程的接收缓冲区和相关数组
int *recvbuf = NULL;
int recvcounts[size];
int displs[size];

if (rank == 0) {
// 根进程计算总的接收数量和每个进程的偏移量
int total_count = 0;
for (int i = 0; i < size; ++i) {
recvcounts[i] = i + 10; // 假设第i个进程发送i+10个整数
displs[i] = total_count;
total_count += recvcounts[i];
}
recvbuf = (int*)malloc(total_count * sizeof(int));
}

// 收集数据
MPI_Gatherv(sendbuf, 10, MPI_INT, recvbuf, recvcounts, displs, MPI_INT, 0, MPI_COMM_WORLD);

// 根进程打印接收到的数据
if (rank == 0) {
printf("Root process has gathered the following data:\n");
for (int i = 0; i < displs[size - 1] + recvcounts[size - 1]; ++i) {
printf("%d ", recvbuf[i]);
}
printf("\n");
free(recvbuf);
}

MPI_Finalize();
return 0;
}

在这个例子中,我们假设每个进程有一个包含10个整数的发送缓冲区 sendbuf,每个整数初始化为该进程的排名乘以10加上索引值。根进程(排名为0的进程)需要准备一个足够大的接收缓冲区 recvbuf 来接收所有其他进程发送的数据。

每个进程调用 MPI_Gatherv,发送它的 sendbuf 中的数据。根进程使用 recvcounts 数组来指定它期望从每个进程接收的数据数量,displs 数组来指定每个进程的数据在接收缓冲区中的偏移量。

使用Modules进行环境变量管

新用户创建后会将系统默认的环境变量填写到用户~/.bashrc 文件下(除intel 环境变量外,其他全部用#注释,默认不生效)。用户可作为参考,根据使用情况进行修改。

(推荐)使用Environment Modules进行环境变量管理

Environment Modules 包是一个简化 shell 初始化的工具,它允许用户在使用 modulefiles 进行会话期间轻松修改其环境。每个模块文件都包含为应用程序配置 shell 所需的信息。模块文件可以由系统上的许多用户共享,并且用户可以拥有自己的集合来补充或替换共享模块文件。

Modules 指令

  • module avail:列出当前 module path 中的所有可用模块文件。
[1907160330@login02 ~]$ module avail

在这里插入图片描述

  • module load MODULEFILE:加载模块文件/类。
[1907160330@login02 ~]$ module load matlab/R2016a
  • module list: 显示已经加载的模块。
[1907160330@login02 ~]$ module list
Currently Loaded Modulefiles:
1) matlab/R2016a
  • module unload MODULEFILE:卸载模块文件/类
[1907160330@login02 ~]$ module unload matlab/R2016a
[1907160330@login02 ~]$ module list
No Modulefiles Currently Loaded.
  • module switch MODULEFILE-A MODULEFILE-B:切换模块

执行命令会卸载模块A,加载模块B

[1907160330@login02 ~]$ module load matlab/R2016a
[1907160330@login02 ~]$ module list
Currently Loaded Modulefiles:
1) matlab/R2016a
[1907160330@login02 ~]$ module switch matlab/R2016a cmake/3.8.1
[1907160330@login02 ~]$ module list
Currently Loaded Modulefiles:
1) cmake/3.8.1
[1907160330@login02 ~]$

自定义MODULEFILE

在管理员配置的MODULE中没有你需要的环境时,普通用户就可以创建属于自己的MODULEFILE文件。

  1. 进入到用户主目录下,使用mkdir命令创建一个名为privatemodules的文件夹。
[1907160330@login02 ~]$ mkdir privatemodules
[1907160330@login02 ~]$ ls
3.6.1 ai_datastore matlab-0 matlab_crash_dump.1099-1 privatemodules soft
  1. 进入到privatemodules文件夹中,创建变量,此处以ffmpeg为例。
[1907160330@login02 privatemodules]$ mkdir ffmpeg
[1907160330@login02 privatemodules]$ cd ffmpeg
[1907160330@login02 ffmpeg]$ vim 4.3.1

此时会打开vim编辑器,在里面输入或者复制如下脚本(英文状态下,按i开始编辑):

#%Module1.0
proc ModulesHelp { } {
puts stderr "\t FFmpeg \n"
}
module-whatis "\t For more information, $module help ffmpeg \n"
conflict modulefile
prepend-path PATH /gpfs/users_home/1907160330/soft/ffmpeg

编辑完成后,按Esc退出编辑模式,然后输入:wq!,保存并退出。

  1. 使用cd ~回到用户主目录,修改当前用户的环境变量文件.bashrc,为module加入我们自建的目录。
[1907160330@login02 ffmpeg]$ cd ~
[1907160330@login02 ~]$ vim .bashrc

英文状态下按i在当中插入:export MODULEPATH=/gpfs/users_home/xxx/privatemodules:$MODULEPATH(xxx为你的用户名)

# .bashrc
# Source global definitions
if [ -f /etc/bashrc ]; then
. /etc/bashrc
fi
# Uncomment the following line if you don't like systemctl's auto-paging feature:
# export SYSTEMD_PAGER=
export MODULEPATH=/gpfs/users_home/1907160330/privatemodules:$MODULEPATH
# User specific aliases and functions

编辑完成后,按Esc退出编辑模式,然后输入:wq!,保存并退出。

  1. 使用source命令,让你的配置文件生效。
[1907160330@login02 ~]$ source .bashrc
-bash: PROMPT_COMMAND: readonly variable
  1. 在使用module av命令查看,即可看到我们自行添加的MODULEFILE。
[1907160330@login02 ~]$ module avail
---------- /gpfs/users_home/1907160330/privatemodules ----------
ffmpeg/4.3.1 python3/3.7.1

请注意,Markdown不支持HTML的一些特定属性,例如id属性(例如<a id="MODULEFILE_306"></a>)。因此,这些元素在Markdown版本中被省略了。此外,Markdown不支持<code>标签中的HTML实体(如&#xff1a;),这些已被转换为其对应的字符(如)。图片和代码块已按Markdown语法进行格式化。