MPI——通信集合

通信集合

树形结构通信

一开始1、3、5、7号进程都将它们的值发送给0、2、4、6。然后0、2、4、6进程都将接收到的值加到它们自己原有的值，整个过程重复两次。

但是这个解决的方案不是最理想的，因为其中一半的进程的工作量与原来方案相同。原来的方案需要0号进程接收comm_sz-1=7次。而新的方案中0号进程只需要3次接收和3次加法，并且其他进程所有的接收与加法操作不超过2次。所以，如果所有进程几乎是同时启动的，那么全局求和需要的总时间将是0号进程所需要的时间，即3次接收和3次加法操作，因为减少了超过50%的总时间。如果comm_sz=1024，则原方案需要0号进程执行1023次接收和加法操作，而采用新的方案只需要10次接收和10次加法操作，则比原方案提高了100倍。

另外一种方式

MPI_Reduce

在MPI里，涉及通信子中所有进程的通信函数称为集合通信（collective communication）.

int MPI_Reduce(
	void * 	input_data_p,		/*in*/
    void * 	output_data_p,		/*out*/
    int 	count,				/*in*/
    MPI_Datatype	datatype,	/*in*/
    MPI_Op 	operator,			/*in*/
    int 	dest_process,		/*in*/
    MPI_Comm	comm,			/*in*/
)

第5个参数operator，它的类型为MPI_Op 是一个预定义的MPI类型。

用于一组N维向量的加法，每个进程上有一个向量

double local_x[N], sum[N];
MPI_Reduce(local_x, sum, N, MPI_DOUBLE, MPI_SUM,0,MPI_COMM_WORLD);

MPI中预定于的归约操作符

运算符值	含义
MPI_MAX	求最大值
MPI_MIN	求最小值
MPI_SUM	求累加和
MPI_PROD	求累乘积
MPI_LAND	逻辑与
MPI_BAND	按位与
MPI_LOR	逻辑或
MPI_BOR	按位或
MPI_LXOR	逻辑异或
MPI_BOXR	按位异或
MPI_MAXLOC	求最大值和最小大的位置
MPI_MINLOC	求最小值和最小值的位置

集合通信与点对点通信

集合通信与点对点通信在多个方向是不同的

1. 在通信子中的所有进程必须调用相同的集合通信函数
2. 每个进程传递给MPI集合通信函数的参数必须是相容的
3. 参数output_data_p只用在destprocess上人，所有进程仍需要传递一个与output_data_p相对应的实际参数，即使它的值NULL。
4. 点对点通信函数是通过标签和通信子来匹配的。

MPI_Allreduce

如果想要所有的进程都可以得到全局总和的结果，以便完成一件大规模的计算。

一个简单的方式可以通过颠倒整棵树来发布全局综合，另外一种简单的方法就是可以让进程之间相互交换部分结果，而不是单向通信。这种通信模式成为蝶形。

int MPI_Allreduce(
	void *	intput_data_p,  //in
    void * 	output_data_p,  //out
    int		count, 			//in
    MPI_Datatype	datatype, 	//in
    MPI_Op	operator,		//in
    MPI_Comm	comm
);

使用颠倒的方式发布结果

蝶形结果的全局求和

广播

在一个集合通信中，如果属于一个进程的数据被发送到通信子中的所有进程，这样的集合通信就叫做广播（broadcast）

int MPI_Bcast(
	void * 	data_p, 	//in/out
    int 	count, 		//in
    MPI_Datatype	datatype,	//in
    int source_proc ,	//in
    MPI_Comm	comm,	//in
);

进程号为source_proc的进程将data_p所引用的内存内容发送给通信子comm中的所有进程。

void Get_input(int my_rank, int comm_sz, double * ap, double *b_p, int * n_p) {
    if(my_rank == 0) {
        printf("Enter a, b and n\n");
        scanf("%lf %lf %d", a_p, b_p, n_p);
    }
    MPI_Bcast(a_p, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
    MPI_Bcast(b_p, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
    MPI_Bcast(n_p, 1, MPI_INT, 0, MPI_COMM_WORLD);
}

数据分发

编写一个程序，用于计算向量和
$$\begin{array}{rl}x\text{x}x+y\text{y}y& =(x_0,x_1,\cdots ,x_{n-1})+(y_0,y_1,\cdots ,y_{n-1})\\ & =(x_0++y_0,x_1+y_1,\cdots ,x_{n-1}+y_{n-1})\\ & =(z_0,z_1,\cdots ,z_{n-1})\\ & =z\text{z}z\end{array}$$
向量的求和的串行实现

void Vector_sum(double x[], double y[], double z[], int n) {
    int i;
    for(i = 0;i < n; i++) {
        z[i] = x[i] + y[i];
    }
}

计算工作由向量的各个分量分别求和组成，我们可能指定各个任务求和的对应分量。各个任务间没有通信，向量的并行加法问题就归结为聚合任务以及将它们分配给核上，如果分量的个数为n，并且我们有comm_sz个核或者进程，那么我们可以简单将连续local_n个向量分量所构成的块，分配给每个进程。称为块划分。

循环划分，用轮转的方式分配变量分量。

第三种划分方法叫做：块-循环划分，用一个循环来分发向量分量所构成的块，而不是分发单个向量分量。

一旦决定如何划分分量，就能很容易地编写向量的并行加法函数：每个进程只要简单地将它所分配到向量分量加起来。

void Parallel_vector_sum() {
    double local_x[],
    double local_y[],
    double local_z[],
    int local_n
} {
   	int local_i;
    for(local_i = 0;local_i < local_n;local_i++) {
        local_z[local_i] - local_x[local_i] + local_y[local_i];
    }
}

散射

例如我们使用块划分法，那么如果0号进程只需要将$1000\sim 1999$号分量发送给1号进程，将$2000\sim 2999$分量分配给2号进程。用这用方法，$1\sim 9$号进程将只需要他们实际使用的向量分量分配存储空间。

int MPI_Scatter(
	void * 	send_buf_p, //in
    int 	send_count,	//in
    MPI_Datatype 	send_type, 	//in
    void * 	recv_buf_p,	//out
    int 	recv_count,	//in
    MPI_Datatype	recv_type, //in
    int		src_proc,	//in
    MPI_Comm	comm 	//in
)

如果通信子comm包含comm_sz个进程，那么MPI_Scatter函数会将send_buf_p所引用的数据分成comm_sz份，第一份给0进程，第二份给1进程，以此类推。send_count参数表示的是发送到每个进程的数据量，而不是send_buf_p所引用的内存的数据量。

一个读取并分发向量的函数

void Read_vector(double local_a[], int local_n, int n ,char vec_name[], int my_rank, MPI_Comm comm) {
    double* a = NULL;
    int i;
    if(my_rank == 0) {
        a = malloc(n*size(double));
        printf("Enter the vector %s\n", vec_name);
        for(int i = 0;i<n;i++) {
            scanf("%lf", &a[i]);
        }
        MPI_Scatter(a, local_n, MPI_DOUBLE, local_a, local_n, MPI_DOUBLE, 0,comm);
        free(a);
    }
    else {
        MPI_Scatter(a, local_n, MPI_DOUBLE, local_a, local_n,MPI_DOUBLE, 0, comm);
    }
}

聚集

这个函数将向量的所有分量都收集到0号进程中，然后由0号进程将所有分量打印出来。

int MPI_Gather(
	void *	send_buf_p, 	//in
    int 	send_count,		//in
    MPI_Datatype 	send_type, 	//in
    void *	recv_buf_p,		//out	
    int 	recv_count,		//in
    MPI_Datatype	recv_type,	//in
    int dest_proc,		//in
    MPI_Comm	comm,	//in
);

在0号进程中，由send_buf_p所引用的内存的数据存储在recv_buf_p的第一个块中，在1号进程中，由send_buf_p所引用的内存的数据存储在recv_buf_p的第二个块中。recv_count指的是每个进程接收到的数据量，而不是所有的接收到的数据量的总和。只有在使用块划分法，并且每个块的大小相同的情况下。

void Print_vector(double local_b[], int local_n, int n, char title[], int my_rank, MPI_Comm comm) {
    double * b = NULL;
    int i;
    if(my_rank == 0) {
        b = malloc(n * sizeof(double));
        MPI_Gather(local_b, local_n, MPI_DOUBLE, b, local_n, MPI_DOUBLE, 0, comm);
        printf("%s\n",title);
        for(i = 1;i<n;i++) {
            printf("%f ",b[i]);
        }
        printf("\n");
        free(b);
    }
    else {
        MPI_Gather(local_b, local_n, MPI_DOUBLE, b, local_n, MPI_DOUBLE, 0, comm);
    }
}

全局聚集

完成矩阵和向量的相乘，如果$A=a_{ij}$是一个$m\times n$的矩阵，$x$是一个具有$n$个分量的向量，那么$y=Ax$就是一个有$m$个分量的向量。
$$y_i=a_{io}{x}_0+a_{i1}{x}_1+\cdots +a_{in-1}{x}_{n-1}$$
串行伪代码

for(i = 0;i<m;i++) {
    y[i] = 0.0;
    for(j = 0;i<n;j++) 
        y[i] += A[i][j] * x[j];
}

在C语言中，经常用一维数据来模拟二维数组，最常见的方法就是将一行内容存储在另一行的后面。
$$\left[\begin{array}{cccc}0& 1& 2& 3\\ 4& 5& 6& 7\\ 8& 9& 10& 11\end{array}\right]$$
这样的数组作为一维数组存储
$$01234567891011$$

void Mat_vect_mult(double A[], double x[], double y[], int m, int n){
    int i,j;
    for(i = 0;i<m;i++) {
        y[i]= 0.0;
        for(j = 0;j<n;j++) {
           	y[i] += A[i * n+j]*x[j];
        }
    }
}

对A进行划分可以以使用$y[i]$的计算包含所需要的A中的元素，所以对y也应该采用块划分。如果A第i行分配给q号进程，y的第i个分量也分配给q号进程。

int MPI_Allgather(
	void * send_buf_p, 	//in
    int send_count, 	//in
    MPI_Datatype 	send_type, 	//in
    void * 	recv_buf_p,	//out
    int 	recv_count,	//in
    MPI_Datatype	recv_type,	//in
    MPI_Comm comm	//in
)

这个函数将每个进程的send_buf_p内容串起来，存储到每个进程的recv_buf_p参数中，通常，recv_count指每个进程接收的数据量。大部分情况下，recv_count的值与send_count的值相同。

MPI矩阵-向量乘法函数

void Mat_vect_mult(double local_A[], double local_x[], double local_y[], int local_m, int n, int local_n, MPI_Comm comm) {
    double *x;
    int local_i, j;
    int local_ok =1;
    x = malloc(n * sizeof(double));
    MPI_Allgather(local_x, local_n, MPI_DOUBLE, x, local_n, MPI_DOUBLE, comm);
    for(local_i = 0; local_i < local_m;local_i ++ ) {
        local_y[local_i] = 0.0;
        for(j = 0;j<n;j++) {
            local_y[local_i] += local_A[local_i*n+j] * x[j];
        }
    }
}