【MPI】基于GEMM的Cannon算法的MPI并行编程优化

最新推荐文章于 2024-05-20 18:16:07 发布
最新推荐文章于 2024-05-20 18:16:07 发布
阅读量411 收藏 1
点赞数 1
分类专栏： 高性能计算 文章标签： 算法 c++ 矩阵
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/weixin_43730667/article/details/131551530
版权

一、并行算法程序实现

在GEMM串行程序源码基础上，分别实现基于行列分块并行算法和cannon算法的两种MPI并行程序
代码debug阶段，采用的是M = N = P =4, np =4情况下进行验证,之后换成M = N = P = 16, np = 16再次验证

|| 0 1 2 3 ||      || 1 0 0 0 ||
|| 4 5 6 7 ||      || 0 1 0 0 ||
|| 8 9 10 11 ||    || 0 0 1 0 ||
|| 12 13 14 15 ||  || 0 0 0 1 ||

1、行列分块算法的程序实现

1.1 阻塞通信版本

  • 相当于A分成4行，B分成4列，每个进程得到一个C，最后C拼起来就是答案
  • 第一步：初始化AB
  • 第二步：确定邻居进程
  • 第三步：本地进程AB进行gemm运算
  • 第四步：B向上收发子块
  • 第五步：继续gemm运算，共迭代np次
  • 第六步：收集D子块，打印
// 一维阻塞版本:sendrecv
// 一维阻塞版本:sendrecv
#include<stdio.h>
#include<mpi.h>
#include<stdlib.h>
#include<unistd.h>
#include<iostream>

#define np 1
#define M 64
#define N 64
#define P 64
#define SIZE M / np
#define UP 0
#define DOWN 1
#define TT 1

void print_myRows(int, float [][SIZE]);
void print_myRowsA(int , float myRows[][P]);
void print_myRowsB(int , float myRows[][SIZE]);

int main(int argc, char* argv[]){
    int myid, size;
	float myRowsA[SIZE][P], myRowsB[P][SIZE], myRowsC[SIZE][SIZE];
    float myRowsD[SIZE][SIZE * np], tmp[SIZE][SIZE * np], c[M][N];
	int nebor[2];
    double start, end;
    int tag_up = 0;

	MPI_Init(&argc, &argv);
	MPI_Comm_rank(MPI_COMM_WORLD, &myid);
	MPI_Comm_size(MPI_COMM_WORLD, &size);
    MPI_Status status;

    // 1. 数据初始化
    for (int i = 0; i < SIZE; i ++){
        for (int j = 0; j < P; j ++){
            myRowsA[i][j] = 0;
        }
    }

    for (int i = 0; i < P; i ++){
        for (int j = 0; j < SIZE; j ++){
            myRowsB[i][j] = 0;
        }
    }

    float cn = myid * SIZE * P;
    for (int i = 0; i < SIZE; i ++){
        for (int j = 0; j < P; j ++){
            myRowsA[i][j] = cn ++;
        }
    }
    // sleep(myid);
    // print_myRowsA(myid, myRowsA);

    int cnt = myid * SIZE;
    for (int i = 0; i < SIZE; i ++) myRowsB[cnt + i][i] = 1.0;
    // sleep(myid);
    // print_myRowsB(myid, myRowsB);

    // 2. 确定邻居进程
    nebor[UP] = (myid == 0) ? np - 1: myid - 1;
    nebor[DOWN] = (myid == np - 1) ? 0: myid + 1;
    
    // 迭代开始，计时启动
    MPI_Barrier(MPI_COMM_WORLD);
    start = MPI_Wtime();
    for (int T = 0; T < np; T ++){
        // 3. 本地进程当前子块A和B进行计算
        for (int i = 0; i < SIZE; i ++){
            for (int j = 0; j < SIZE; j ++){
                myRowsC[i][j] = 0;
            }
        }

        for (int i = 0; i < SIZE; i ++){
            for (int j = 0; j < SIZE; j ++){
                for (int k = 0; k < P; k ++){
                    myRowsC[i][j] += myRowsA[i][k] * myRowsB[k][j];
                    // 将结果放入D中
                    myRowsD[i][(T +myid) % np * SIZE + j] = myRowsC[i][j];  
                }
            }
        }

        // sleep(myid);
        // print_myRows(myid, myRowsC);

        // 4. 发送和接收B子块
        MPI_Sendrecv(&myRowsB[0][0], P * SIZE, MPI_FLOAT, nebor[UP], tag_up, &myRowsB[0][0], P * SIZE, MPI_FLOAT, nebor[DOWN], tag_up, MPI_COMM_WORLD, &status);
        // MPI_Send(&myRowsB[0][0], P * SIZE, MPI_FLOAT, nebor[UP], tag_up, MPI_COMM_WORLD);
        // MPI_Recv(&myRowsB[0][0], P * SIZE, MPI_FLOAT, nebor[DOWN], tag_up, MPI_COMM_WORLD, &status);
    }
    MPI_Barrier(MPI_COMM_WORLD); /* IMPORTANT */
    end = MPI_Wtime();

    // 5. 收集D子块到tmp中，存放到c中
    for (int i = 0; i < SIZE; i ++){
		for (int j = 0; j < SIZE * np; j ++){
			tmp[i][j] = myRowsD[i][j];
		}
	}
    MPI_Barrier(MPI_COMM_WORLD);
	MPI_Gather(tmp, SIZE * SIZE * np, MPI_FLOAT, c, SIZE * SIZE * np, MPI_FLOAT, 0, MPI_COMM_WORLD);
    
    // 6. 输出汇总后的结果
	if (0 == myid){
		fprintf(stderr, "\n收集后结果\n");
		for (int i = 0; i < M; i ++){
			for (int j = 0; j < N; j ++){
				fprintf(stderr, "%.3f\t", c[i][j]);
			}
			fprintf(stderr, "\n");
		}
		fprintf(stderr, "\n");
	}

    printf("Result in process %d, Runtime = %f\n", myid, end-start);

    MPI_Finalize();
    return 0;
}

// 打印函数
void print_myRows(int myid, float myRows[][SIZE]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < SIZE;i ++){
		for (j = 0; j < SIZE; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}
void print_myRowsA(int myid, float myRows[][P]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < SIZE;i ++){
		for (j = 0; j < P; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}
void print_myRowsB(int myid, float myRows[][SIZE]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < P;i ++){
		for (j = 0; j < SIZE; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}

1.2 非阻塞通信版本

  • 使用了MPI_Isend函数，可以让通信和gemm计算重叠
  • 如果是没有引入buff而是使用myRowsB进行传递，那么发送可以在gemm前，而接收必须在gemm后，因为gemm运算需要原始B
  • 所以引入了buffB,则在计算gemm之前就可以开启MPI_Isend``MPI_Irecv，gemm完成后MPI_Waitall接收完成
// 非阻塞通信版本：Isend 加上buff
#include<stdio.h>
#include<mpi.h>
#include<stdlib.h>
#include<unistd.h>
#include<iostream>

#define np 1
#define M 64
#define N 64
#define P 64
#define SIZE M / np
#define UP 0
#define DOWN 1
#define TT 1

void print_myRows(int, float [][SIZE]);
void print_myRowsA(int , float myRows[][P]);
void print_myRowsB(int , float myRows[][SIZE]);

int main(int argc, char* argv[]){
    int myid, size;
	float myRowsA[SIZE][P], myRowsB[P][SIZE], myRowsC[SIZE][SIZE];
    float myRowsD[SIZE][SIZE * np], tmp[SIZE][SIZE * np], c[M][N];
	int nebor[2];
    double start, end;
    float buffB[P][SIZE];

	MPI_Init(&argc, &argv);
	MPI_Comm_rank(MPI_COMM_WORLD, &myid);
	MPI_Comm_size(MPI_COMM_WORLD, &size);
    // MPI_Status status;
    MPI_Status status[2], status1;
    MPI_Request request[2];

    // 1. 数据初始化
    for (int i = 0; i < SIZE; i ++){
        for (int j = 0; j < P; j ++){
            myRowsA[i][j] = 0;
        }
    }

    for (int i = 0; i < P; i ++){
        for (int j = 0; j < SIZE; j ++){
            myRowsB[i][j] = 0;
        }
    }

    float cn = myid * SIZE * P;
    for (int i = 0; i < SIZE; i ++){
        for (int j = 0; j < P; j ++){
            myRowsA[i][j] = cn ++;
        }
    }
    // sleep(myid);
    // print_myRowsA(myid, myRowsA);

    int cnt = myid * SIZE;
    for (int i = 0; i < SIZE; i ++) myRowsB[cnt + i][i] = 1.0;
    // sleep(myid);
    // print_myRowsB(myid, myRowsB);

    // 2. 确定邻居进程
    nebor[UP] = (myid == 0) ? np - 1: myid - 1;
    nebor[DOWN] = (myid == np - 1) ? 0: myid + 1;
    
    MPI_Barrier(MPI_COMM_WORLD);
    start = MPI_Wtime();
    for (int T = 0; T < np; T ++){
        // 4. 发送和接收B子块
        int tag_up = 0, tag_down = 1;
        MPI_Isend(&myRowsB[0][0], P * SIZE, MPI_FLOAT, nebor[UP], tag_up, MPI_COMM_WORLD, &request[0]);
        MPI_Irecv(&buffB[0][0], P * SIZE, MPI_FLOAT, nebor[DOWN], tag_up, MPI_COMM_WORLD, &request[1]);
        
        // 3. 本地进程当前子块A和B进行计算
        for (int i = 0; i < SIZE; i ++){
            for (int j = 0; j < SIZE; j ++){
                myRowsC[i][j] = 0;
            }
        }

        for (int i = 0; i < SIZE; i ++){
            for (int j = 0; j < SIZE; j ++){
                for (int k = 0; k < P; k ++){
                    myRowsC[i][j] += myRowsA[i][k] * myRowsB[k][j];
                    myRowsD[i][(T +myid) % np * SIZE + j] = myRowsC[i][j];  
                }
            }
        }

        // sleep(myid);
        // print_myRows(myid, myRowsC);

        MPI_Waitall(2, &request[0], &status[0]);

        for (int i = 0; i < P; i ++){
            for (int j = 0; j < SIZE; j ++){
                myRowsB[i][j] = buffB[i][j];
            }
        }
    }
    MPI_Barrier(MPI_COMM_WORLD); /* IMPORTANT */
    end = MPI_Wtime();

    printf("Result in process %d, Runtime = %f\n", myid, end-start);

    // 5. 收集D子块到tmp中，存放到c中
    for (int i = 0; i < SIZE; i ++){
		for (int j = 0; j < SIZE * np; j ++){
			tmp[i][j] = myRowsD[i][j];
		}
	}
    MPI_Barrier(MPI_COMM_WORLD);
	MPI_Gather(tmp, SIZE * SIZE * np, MPI_FLOAT, c, SIZE * SIZE * np, MPI_FLOAT, 0, MPI_COMM_WORLD);
    
    // 6. 输出汇总后的结果
	if (0 == myid){
		fprintf(stderr, "\n收集后结果\n");
		for (int i = 0; i < M; i ++){
			for (int j = 0; j < N; j ++){
				fprintf(stderr, "%.3f\t", c[i][j]);
			}
			fprintf(stderr, "\n");
		}
		fprintf(stderr, "\n");
	}

    MPI_Finalize();
    return 0;
}
    
void print_myRows(int myid, float myRows[][SIZE]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < SIZE;i ++){
		for (j = 0; j < SIZE; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}
void print_myRowsA(int myid, float myRows[][P]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < SIZE;i ++){
		for (j = 0; j < P; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}
void print_myRowsB(int myid, float myRows[][SIZE]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < P;i ++){
		for (j = 0; j < SIZE; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}

1.3 重复非阻塞通信版本

  • 和上述非阻塞版本类似，只不过把初始化提前，可以共用初始化

代码

// 重复非阻塞通信版本：MPI_Send_Init 加上buff
#include<stdio.h>
#include<mpi.h>
#include<stdlib.h>
#include<unistd.h>
#include<iostream>

#define np 4
#define M 16
#define N 16
#define P 16
#define SIZE M / np
#define UP 0
#define DOWN 1
#define TT 4

void print_myRows(int, float [][SIZE]);
void print_myRowsA(int , float myRows[][P]);
void print_myRowsB(int , float myRows[][SIZE]);

int main(int argc, char* argv[]){
    int myid, size;
	float myRowsA[SIZE][P], myRowsB[P][SIZE], myRowsC[SIZE][SIZE];
    float myRowsD[SIZE][SIZE * np], tmp[SIZE][SIZE * np], c[M][N];
	int nebor[2];
    double start, end;
    int buffsize = P * N * sizeof(float);
    float buffB[P][SIZE];

	MPI_Init(&argc, &argv);
	MPI_Comm_rank(MPI_COMM_WORLD, &myid);
	MPI_Comm_size(MPI_COMM_WORLD, &size);
    // MPI_Status status;
    MPI_Status status[2], status1;
    MPI_Request request[2];

    // 1. 数据初始化
    for (int i = 0; i < SIZE; i ++){
        for (int j = 0; j < P; j ++){
            myRowsA[i][j] = 0;
        }
    }

    for (int i = 0; i < P; i ++){
        for (int j = 0; j < SIZE; j ++){
            myRowsB[i][j] = 0;
        }
    }

    float cn = myid * SIZE * P;
    for (int i = 0; i < SIZE; i ++){
        for (int j = 0; j < P; j ++){
            myRowsA[i][j] = cn ++;
        }
    }
    // sleep(myid);
    // print_myRowsA(myid, myRowsA);

    int cnt = myid * SIZE;
    for (int i = 0; i < SIZE; i ++) myRowsB[cnt + i][i] = 1.0;
    // sleep(myid);
    // print_myRowsB(myid, myRowsB);

    // 2. 确定邻居进程
    nebor[UP] = (myid == 0) ? np - 1: myid - 1;
    nebor[DOWN] = (myid == np - 1) ? 0: myid + 1;
    
    // 5. 初始化B的发送过程
    int tag_up = 0;
    MPI_Send_init(&myRowsB[0][0], P * SIZE, MPI_FLOAT, nebor[UP], tag_up, MPI_COMM_WORLD, &request[0]);
    MPI_Recv_init(&buffB[0][0], P * SIZE, MPI_FLOAT, nebor[DOWN], tag_up, MPI_COMM_WORLD, &request[1]);

    MPI_Barrier(MPI_COMM_WORLD);
    start = MPI_Wtime();
    for (int T = 0; T < np; T ++){
        // 4. 发送和接收B子块
        for (int i = 0; i < P; i ++){
            for (int j = 0; j < SIZE; j ++){
                buffB[i][j] = myRowsB[i][j];
            }
        }
        MPI_Startall(2, &request[0]);
        
        // 3. 本地进程当前子块A和B进行计算
        for (int i = 0; i < SIZE; i ++){
            for (int j = 0; j < SIZE; j ++){
                myRowsC[i][j] = 0;
            }
        }

        for (int i = 0; i < SIZE; i ++){
            for (int j = 0; j < SIZE; j ++){
                for (int k = 0; k < P; k ++){
                    myRowsC[i][j] += myRowsA[i][k] * myRowsB[k][j];  
                    myRowsD[i][(T +myid) % np * SIZE + j] = myRowsC[i][j];  
                }
            }
        }

        // sleep(myid);
        // print_myRows(myid, myRowsC);

        MPI_Waitall(2, &request[0], &status[0]);
        
        for (int i = 0; i < P; i ++){
            for (int j = 0; j < SIZE; j ++){
                myRowsB[i][j] = buffB[i][j];
            }
        }
        // sleep(myid);
        // print_myRowsB(myid, myRowsB);
    }
    MPI_Barrier(MPI_COMM_WORLD); /* IMPORTANT */
    end = MPI_Wtime();

    printf("Result in process %d, Runtime = %f\n", myid, end-start);

    // 5. 收集D子块到tmp中，存放到c中
    for (int i = 0; i < SIZE; i ++){
		for (int j = 0; j < SIZE * np; j ++){
			tmp[i][j] = myRowsD[i][j];
		}
	}
    MPI_Barrier(MPI_COMM_WORLD);
	MPI_Gather(tmp, SIZE * SIZE * np, MPI_FLOAT, c, SIZE * SIZE * np, MPI_FLOAT, 0, MPI_COMM_WORLD);
    
    // 6. 输出汇总后的结果
	if (0 == myid){
		fprintf(stderr, "\n收集后结果\n");
		for (int i = 0; i < M; i ++){
			for (int j = 0; j < N; j ++){
				fprintf(stderr, "%.3f\t", c[i][j]);
			}
			fprintf(stderr, "\n");
		}
		fprintf(stderr, "\n");
	}

    MPI_Finalize();
    return 0;
}


void print_myRows(int myid, float myRows[][SIZE]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < SIZE;i ++){
		for (j = 0; j < SIZE; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}
void print_myRowsA(int myid, float myRows[][P]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < SIZE;i ++){
		for (j = 0; j < P; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}
void print_myRowsB(int myid, float myRows[][SIZE]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < P;i ++){
		for (j = 0; j < SIZE; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}

2、cannon算法的程序实现

  • 第一步：初始化A、B矩阵
  • 第二步：建立笛卡尔坐标，确定邻居今称
  • 第三步：初始对准，A向左移动i，B向上移动j
  • 第四步：计算gemm
  • 第五步：累加gemm结果
  • 第六步：循环移位，A向左移动1，B向上移动1
  • 第七步：重复循环sqrt(np)次
  • 第八步：收集C子块

遇到的问题

  1. 移位的时候，邻居进程判断是正确的，但是有的进程接收或者发送不成功，会出错，因此在初始化的时候把MPI_Sendrecv换成了MPI_Send``MPI_Recv这样每个都阻塞通信一定能保证初始化没问题
  2. 收集C子块，MPI_Gather`是按照进程顺序依次排布，和我们想要的二维分块是不一致的，因此采用的是收集是按照每个进程按行排布，然后对矩阵进行重排得到结果

代码

// Cannon版本
#include<stdio.h>
#include<mpi.h>
#include<stdlib.h>
#include<unistd.h>
#include<iostream>
#include<math.h>

#define np 16
#define M 16
#define N 16
#define P 16
#define UP 0
#define DOWN 1
#define LEFT 2 
#define RIGHT 3
#define TT 4 //sqrt(np)
#define SIZE 4 // N / sqrt(np)
void print_myRows(int, float [][SIZE]);

int main(int argc, char* argv[]){
    int myid, size;
	float myRowsA[SIZE][SIZE], myRowsB[SIZE][SIZE], myRowsC[SIZE][SIZE];
    float tmp[SIZE][SIZE], c[np][SIZE * SIZE], c1[M][N];
	int nebor[2];
    double start, end;
    int buffsize = SIZE * SIZE * sizeof(float);
    float buffA[SIZE][SIZE], buffB[SIZE][SIZE], buff[buffsize];

    int dims[2] = {TT, TT};
	int periods[2] = {0, 0};
	int reorder = 0;
    int coords[2];
    int nbrs[4];
    int position;

	MPI_Init(&argc, &argv);
	MPI_Comm_size(MPI_COMM_WORLD, &size);
    // MPI_Status status;
    MPI_Status status[4], status1;
    MPI_Request request[4];
    MPI_Comm cartcomm;

    // 为笛卡尔拓扑创建一个新的通信子
    MPI_Cart_create(MPI_COMM_WORLD, 2, dims, periods, reorder, &cartcomm);	
    // 获取当前进程在通信子中的秩
    MPI_Comm_rank(cartcomm, &myid);	
    // 笛卡尔转换函数，把秩转换成坐标
    MPI_Cart_coords(cartcomm, myid, 2, coords);	

    // 1. 初始化二维进程网络
     for (int i = 0; i < SIZE; i ++){
        for (int j = 0; j < SIZE; j ++){
            myRowsA[i][j] = myRowsB[i][j] = myRowsC[i][j] = 0;
        }
    }

    int row = coords[0], col = coords[1];
    float cn = row * SIZE * P + col * SIZE;
    for (int i = 0; i < SIZE; i ++){
        for (int j = 0; j < SIZE; j ++){
            myRowsA[i][j] = cn + i * P + j;
        }
    }

    if (row == col){
        for (int i = 0; i < SIZE; i ++){
            myRowsB[i][i] = 1.0;
        }
    }

    // sleep(myid);
    // print_myRows(myid, myRowsA);

    // sleep(myid);
    // print_myRows(myid, myRowsB);

    // 2. 确定邻居进程
    // 笛卡尔移位操作,获取当前坐标位置的上下位置秩
    MPI_Cart_shift(cartcomm, 0, 1, &nbrs[UP], &nbrs[DOWN]);	
    // 笛卡尔移位操作,获取当前坐标位置的左右位置秩
    MPI_Cart_shift(cartcomm, 1, 1, &nbrs[LEFT], &nbrs[RIGHT]);

    // 3. 起始对准
    int tag_left = 0, tag_up = 1;
    int leftid = col, upid = row * TT;
    int rightid = col, downid = row * TT;
    if (row != 0){
        leftid = row * TT + (col - row + TT) % TT;
        rightid = row * TT + abs(col + row) % TT;
    }
    sleep(myid);
    std::cout << "rank: " << myid << "leftid: " << leftid << "rightid: " << rightid << std::endl;
    if (col != 0){
        upid = (row - col + TT) % TT * TT + col;
        downid = (row + col) % TT * TT + col;
    }
    sleep(myid);
    std::cout << "rank: " << myid << "upid: " << upid << "downid: " << downid << std::endl;
    
    if (row != 0){
        MPI_Send(&myRowsA[0][0], SIZE * SIZE, MPI_FLOAT, leftid, tag_left, MPI_COMM_WORLD);
        MPI_Recv(&myRowsA[0][0], SIZE * SIZE, MPI_FLOAT, rightid, tag_left, MPI_COMM_WORLD, &status1);
    }
        // MPI_Sendrecv(&myRowsA[0][0], SIZE * SIZE, MPI_FLOAT, leftid, tag_left, &myRowsA[0][0], SIZE * SIZE, MPI_FLOAT, rightid, tag_left, MPI_COMM_WORLD, &status1);
    sleep(myid);
    print_myRows(myid, myRowsA);

    if (col != 0){
        MPI_Send(&myRowsB[0][0], SIZE * SIZE, MPI_FLOAT, upid, tag_up, MPI_COMM_WORLD);
        MPI_Recv(&myRowsB[0][0], SIZE * SIZE, MPI_FLOAT, downid, tag_up, MPI_COMM_WORLD, &status1);
    }
    
        // MPI_Sendrecv(&myRowsB[0][0], SIZE * SIZE, MPI_FLOAT, upid, tag_up, &myRowsB[0][0], SIZE * SIZE, MPI_FLOAT, downid, tag_up, MPI_COMM_WORLD, &status1);
    sleep(myid);
    print_myRows(myid, myRowsB);


    leftid = (coords[1] == 0) ? (coords[0] + 1) * TT - 1: nbrs[LEFT];
    rightid = (coords[1] == TT - 1) ? coords[0] * TT : nbrs[RIGHT];
    upid = (coords[0] == 0) ? (TT - 1) * TT + coords[1]: nbrs[UP];
    downid = (coords[0] == TT - 1) ? coords[1]: nbrs[DOWN];
    MPI_Send_init(&myRowsA[0][0], SIZE * SIZE, MPI_FLOAT, leftid, tag_left, MPI_COMM_WORLD, &request[0]);
    MPI_Recv_init(&buffA[0][0], SIZE * SIZE, MPI_FLOAT, rightid, tag_left, MPI_COMM_WORLD, &request[1]);
    MPI_Send_init(&myRowsB[0][0], SIZE * SIZE, MPI_FLOAT, upid, tag_up, MPI_COMM_WORLD, &request[2]);
    MPI_Recv_init(&buffB[0][0], SIZE * SIZE, MPI_FLOAT, downid, tag_up, MPI_COMM_WORLD, &request[3]);
    
    MPI_Barrier(MPI_COMM_WORLD);
    start = MPI_Wtime();
    for (int T = 0; T < TT; T ++){
        // 4. 本地进程A与B进行gemm运算
        MPI_Startall(4, &request[0]);
        for (int i = 0; i < SIZE; i ++){
            for (int j = 0; j < SIZE; j ++){
                for (int k = 0; k < SIZE; k ++){
                    myRowsC[i][j] += myRowsA[i][k] * myRowsB[k][j];
                }
            }
        }
        // 5. 循环移位 向左1 向上1
        MPI_Waitall(4, &request[0], &status[0]);
        for (int i = 0; i < SIZE; i ++){
            for (int j = 0; j < SIZE; j ++){
                myRowsA[i][j] = buffA[i][j];
                myRowsB[i][j] = buffB[i][j];
            }
        }
        // sleep(myid);
        // print_myRows(myid, myRowsA);
        // sleep(myid);
        // print_myRows(myid, myRowsB); 
    }
    MPI_Barrier(MPI_COMM_WORLD); /* IMPORTANT */
    end = MPI_Wtime();
    printf("Result in process %d, Runtime = %f\n", myid, end-start);

    for (int i = 0; i < 4; i ++) MPI_Request_free(&request[i]);

    // sleep(myid);
	
    // print_myRows(myid, myRowsC);

    // if (myid == 1) sleep(1);
	for (int i = 0; i < SIZE; i ++){
		for (int j = 0; j < SIZE; j ++){
			tmp[i][j] = myRowsC[i][j];
		}
	}
    MPI_Barrier(MPI_COMM_WORLD);
	MPI_Gather(tmp, SIZE * SIZE, MPI_FLOAT, c, SIZE * SIZE, MPI_FLOAT, 0, MPI_COMM_WORLD);

    int a, b;
    for (int i = 0; i < np; i ++){
        for (int j = 0; j < SIZE * SIZE; j ++){
            a = i / TT * SIZE + j / SIZE;
            b = i % TT * SIZE + j % SIZE;
            c1[a][b] = c[i][j];
        }
    }
    
    // 输出结果
	if (0 == myid){
		fprintf(stderr, "\n收集后结果\n");
		for (int i = 0; i < M; i ++){
			for (int j = 0; j < N; j ++){
				fprintf(stderr, "%.3f\t", c1[i][j]);
			}
			fprintf(stderr, "\n");
		}
		fprintf(stderr, "\n");
	}

    MPI_Finalize();

    return 0;
}

void print_myRows(int myid, float myRows[][SIZE]){
	int i, j;
	printf("Result in process %d:\n", myid);
	for (i = 0; i < SIZE;i ++){
		for (j = 0; j < SIZE; j ++){
			printf("%1.3f\t", myRows[i][j]);
		}
		printf("\n");
	}
	MPI_Barrier(MPI_COMM_WORLD);
}

二、并行性能测试

对各版本算法的并行程序，进行强可拓展性测试，注意计算时间只统计内核gemm计算过程

  • 多个并行规模：进程数量由少到多，如NP为1，2，4，8，16，64，100，144等
  • 多个问题规模：小矩阵阶数，如M/N/P=512/512/512；大矩阵阶数，如M/N/P=5120/5120/5120，10240/10240/10240
  • 不同算法版本之间，对比gemm计算时间及并行效率，分析、理解各类算法、各种通信版本的优势与不足
  • 注意每次run的时候-n是不同的，每个节点大概只能算32500个数，大概都是在np的平方取得最优

Canon性能测试（gemm4_2.cpp)

np\M64128256360512720144018002160
10.0023700.022207
40.0008700.0064610.0458710.188393
160.0039970.0044100.0124300.0530790.0951600.418262
640.0999340.0819550.0767810.0947751.390502
1000.4216360.2359590.2438591.0871012.240380
1440.3182400.3762801.0801062.0601234.015687

一维阻塞通信性能测试（gemm1_1.cpp)

np\M6412825651210242048
10.0035100.0309690.3005763.179691
40.0013470.0105470.0600430.0651834
160.0045100.0026560.0150440.0139472
640.1153890.1544230.0716520.0767810.418964
1000.908158（0.980470）（1000）4.395905（4.443350）（2000）
1441.528170（1.640057）（1152）3.819908（3.938976）（2016）

一维非阻塞通信性能测试（gemm2_2.cpp)

np\M6412825651210242048
10.0037870.031790.3849214.174103
40.0012990.0080590.0595920.5337996.151261
160.0016530.0026030.0156930.1531791.32524
640.1685530.4486224.008370
1000.587287（500）0.779388（1000）5.275981（2000）
1441.592765（576）1.483951（1152）3.939960（2016）

一维重叠非阻塞通信性能测试（gemm3_1.cpp)

np\M6412825651210242048
10.0296492.2694
40.0011600.5353096.455654
160.002762（0.012021）0.132543（0.142638）1.221950（1.223053）14.979193
640.084543（0.122887）0.522625（0.552587）4.068016
1000.0799311（0.895811）（1000）5.820327（5.906116）（2000）
1441.472875（1.594119）（1152）4.490045（4.620051）（2016）

👍👍👍总结

  • 512以下，Canon < 非阻塞<= 阻塞 < 重叠，这个时候计算比较小，等待通讯的时间比较长，所以MPI_Waitall并不是很好；同时小进程表现良好，但是随着维数增大会有不断增大的np最优值出现
  • 512以上，1000以上，2000一下，大进程有效，一维时间指数增长，阻塞 < 非阻塞 < 重叠；Canon增加缓慢，进程变大后效果比一维好
  • 更大进程计算的时间一定大于通讯时间，所以重叠一定会有优势，Canon中也采用的重叠，会有更好的效果
  • 由于矩阵太大计算不了，浅浅预估一下，所以并不是非阻塞或者重叠一定好，不同的进程也要匹配不同的矩阵，要根据矩阵大小而定

存在问题

  • 一维分割必须得整除，维度和np是要能匹配上的
  • 进程数比较小的时候，负载均衡问题比较小
  • 负载均衡
    一维阻塞
    负载均衡问题
    16 512
    在这里插入图片描述
    16 64
    在这里插入图片描述
    64 512
    在这里插入图片描述
    在这里插入图片描述
    负载均衡在矩阵比较小和矩阵比较大时会出现偏差比较大的情况

GEMM 矩阵乘法

串行算法描述

  1. 矩阵太大怎么办?计算时间长 节点内存不够用
  • 基于矩阵一维分块的行列划分并行算法
  • 基于矩阵二维分块的cannon并行算法

行列划分并行算法(矩阵一维分块)

  1. 数据分块
    在这里插入图片描述
    np为处理器个, 对A和B进行分块,注意每个处理器中不一定是一行,可能是多行,所以得到C的都是subC子块而不是一个数.每个处理器从年初的是一个RowC,包含一行的subC
  2. 计算
    在这里插入图片描述
    每次计算的是相同处理器的A和B,通过消息传递,能够让处理器内部计算,完成所有A和B的计算.

cannon并行算法(矩阵二维分块)

  1. 分块方式

    相当于是一维的行块或列块,分割成更小的块了
  2. 思路推导

    subA和subB相乘,然后再求和.
  3. 基本思想
    在这里插入图片描述
  4. 流程
    在这里插入图片描述

在这里插入图片描述
在这里插入图片描述

指标

在这里插入图片描述

如何提高并行效率

在这里插入图片描述

编程步骤

  1. 行列划分并行算法
    在这里插入图片描述
  2. 阻塞接口实现版本
    在这里插入图片描述
    在这里插入图片描述
  3. 非阻塞接口实现版本
    在这里插入图片描述
  4. 非阻塞接口——重叠通信和计算
    在这里插入图片描述
  5. 非阻塞接口——重复非阻塞通信
    通讯的数据、邻居都是相同的等等
    在这里插入图片描述
  6. 矩阵二维分块
    在这里插入图片描述
张西瓜.�
关注
  • 1
    点赞
  • 踩
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
MPI 并行编程与优化
08-10
本ppt简要介绍了集群环境下的MPI 并行编程及如何进行优化
高性能计算实验——矩阵乘法基于MPI的并行实现及优化
weixin_43973089的博客
04-19 1万+
高性能计算实验——矩阵乘法基于MPI的并行实现及优化1.实验目的1.1.通过MPI实现通用矩阵乘法1.2.基于MPI的通用矩阵乘法优化1.3.改造实验1成矩阵乘法库函数2.实验过程和核心代码2.1.通过MPI实现通用矩阵乘法2.2.通用矩阵乘法优化2.3.改造实验1成矩阵乘法库函数3.实验结果3.1.通过MPI实现通用矩阵乘法3.2.基于MPI的通用矩阵乘法优化3.3.改造实验1成矩阵乘法库函数4.实验感想 1.实验目的 1.1.通过MPI实现通用矩阵乘法 熟练掌握MPI编程方法，并将通用矩阵乘法转为MPI
MPI并行计算——点对点通信（求Cos(x)dx的积分）
a429367172的博客
03-15 1116
环境 centos7 mpi 3.2 要求 求cosxdx在a到b区间的积分可以由上图方法近似，图中p是process进程的数量，n是划分小长方形的数量，在输入的时候给出。 代码 注：本代码输入的时候b和a是角度值 #include "mpi.h" #include <stdio.h> #include <math.h> #include <stri...
MPI并行编程在矩阵运算耗时分析
u010728372的博客
03-19 5279
两个矩阵相乘，例子中是1000*2000的矩阵A与2000*3000的矩阵B相乘，矩阵随机产生 运行平台是ubuntu14.04 LTS （一个物理cpu，两个逻辑cpu，每个cpu有2个core，不支持超线程） 注：若物理cpu个数*每个cpu核心数不等于逻辑cpu个数，则cpu支持超线程。 普通c语言程序计算矩阵耗时：109s top命令查看（single是该程序）
并行算法优化（1）
weixin_51942493的博客
06-25 797
通常支持超线程的多核处理器能够使用的线程数最多是物理核心数的2倍X86 流加载/流存储： SSE中的 prefetch指令 可以实现软件预取技术提高多路系统中多核处理器之间通信的带宽。 原理：访问存储器的速度与距离处理器的距离有关，为了满足分配的从“近端”起。 线程亲和性 gcc的环境变量：GOMP_CPU_AFFINITY Linux系统：numactl工具 OpenMP、pthread 性能分析工具： gprof ：通过在编译时插入代码来分析程序 nvprof ：NVIDIA开发，用于运行在GPU上
MPI.rar_MPI_并行 mpi 算法_并行计算
09-24
编写一个MPI并行程序，计算积分：要求写出算法思想、计算步骤、程序代码、重要变量和代码语句的注释。
mpi.rar_MPI_mpi并行_pdf_并行编程
09-20
消息传递并行编程环境mpi pdf
mpi-programming.rar_MPI_并行 mpi 算法
09-24
mpi编程，并行计算实验实现的十四个小算法，效率，算法简单易懂，希望对初学mpi编程的人有所帮助，有详细注释。
MPI.tar.gz_MPI_mpi并行计算_mpi矩阵_并行算法
09-19
MPI并行计算的一些实例，包括向量矩阵计算，广播等并行算法～
高性能计算之并行编程技术—— MPI并行程序设计.zip
12-28
高性能计算之并行编程技术—— MPI并行程序设计
高性能计算机部件mpi,MPI并行代码生成与通信优化技术研究
weixin_36135382的博客
07-28 341
MPI并行代码生成与通信优化技术研究科学与工程计算中重大挑战性问题的涌现和微电子技术的革新，扎实推进高性能计算机的发展。随着高性能计算机的普及和应用，MPI并行模式已成为高性能计算平台普遍支持的重要编程手段之一。因此，MPI并行化编译是有效解决开发高效并行程序和继承经典串行软件问题的重要途径。本文对MPI并行化编译中并行代码自动生成的相关技术进行了深入研究，主要内容包括代码生成和通信优化。本文的主...
MPI学习二 MPI并行程序的两种基本模式
s_sunnyy的博客
02-27 8762
MPI的两种最基本的并行程序设计模式 即对等模式和主从模式。  　　对等模式：各个部分地位相同，功能和代码基本一致，只不过是处理的数据或对象不同，也容易用同样的程序来实现。  　　主从模式：分为主进程和从进程，程序通信进程之间的一种主从或依赖关系 。MPI程序包括两套代码，主进程运行其中一套代码，从进程运行另一套代码。 1. 对等模式MPI程序设计 1.1 Jaco
MPI系列: 并行IO性能优化究竟是怎么玩的呢？
架构师技术联盟
03-02 5180
浓缩干货：MPI-IO对单个进程使用“数据筛选”技术优化，对聚合IO使用“两阶段IO”技术优化。并行存储系统的文件条块长度与严重影响性能。1997年，MPI 2.0标准...
MPI编程及性能优化
热门推荐
GarfieldEr007的专栏
04-15 1万+
第1节 MPI简介 1.1 MPI及其历史 与OpenMP相似，消息传递接口（Message Passing Interface，简称MPI）是一种编程接口标准，而不是一种具体的编程语言。该标准是由消息传递接口论坛（Message Passing Interface Forum，简称MPIF）发起讨论并进行规范化的。 MPI标准从1992年开始起草，1994年发布第一个版本MPI-1（MPI
素数筛法以及其并行优化(MPI)
qq_35571432的博客
10-07 3562
1.素数筛法1.创建一个自然数列表2,3,4……n 2.令k=2，把在列表内的所有k的倍数删去 3.然后找出最小的那个比k大的数 重复2.3步骤，知道k^2>=n2.并行思想1.交叉数据分解 即把数组元素进行分解，例如： 进程0负责2+k,2+2k,2+3k 进程1负责3+k,3+2k,3+3k 这个很容易确定哪个进程负责哪个下标，但是他的缺点也是显而易见的，最后还是需要规约化和广播操
【高性能计算】OpenMP/MPI优化技巧
weixin_45570710的博客
08-12 2644
一些关于OpenMP/MPI使用的经验
MPI并行优化高斯朴素贝叶斯算法 - 通过身高、体重和肺活量推测性别（机器学习）
admiz的博客
11-19 1284
MPI并行设计 在读取数据时使用串行，其余for循环都采用MPI并行。具体步骤如下： 1. 0进程读取数据文件，生成dataSet数组。 2.通过MPI_Scatterv()函数将dataSet数组分为多个长度不等的数据块，向各进程散射数据块。 3.各进程进行求和，通过MPI_Gather()函数将结果返回0进程累加，从而得出平均值。 4.0进程通过MPI_Bcast()函数将平均值(r)广播到所有进程。 5.各进程计算标准差的部分内容：。 6.0进程完成剩余的标准差计算得出标准差SD。 7.得到高斯分布
【算法】【二叉树，DFS，哈希集合，分类讨论】力扣1110. 删点成林
最新发布
未来就在脚下。
05-20 492
本题通过DFS遍历二叉树，结合分类讨论的方法，逐步删除指定节点并生成新的森林。该算法有效地处理了节点删除后树结构的调整问题，并通过哈希集合优化删除判断的时间复杂度，最终实现了高效的解决方案。
请使用C语言编写完整的PageRank算法，使用MPI并行编程
05-27
以下是使用C语言编写完整的PageRank算法，使用MPI并行编程的示例代码： ```c #include <stdio.h> #include <stdlib.h> #include <mpi.h> #define N 100 // 网页总数 #define DAMPING_FACTOR 0.85 // 阻尼系数 #define EPSILON 0.0001 // 收敛精度 int main(int argc, char** argv) { int rank, size; // 进程编号和总进程数 int i, j; // 网页编号 int iter = 0; // 迭代次数 double diff, diff_all; // 差值和全局差值 double pr[N], pr_new[N], pr_all[N]; // 网页PageRank值 double link[N][N]; // 网页链接矩阵 // 初始化MPI MPI_Init(&argc, &argv); MPI_Comm_rank(MPI_COMM_WORLD, &rank); MPI_Comm_size(MPI_COMM_WORLD, &size); // 初始化网页链接矩阵 for (i = 0; i < N; i++) { for (j = 0; j < N; j++) { link[i][j] = (double)rand() / RAND_MAX; } } // 初始化网页PageRank值 for (i = 0; i < N; i++) { pr[i] = 1.0 / N; } // 迭代计算PageRank值 do { // 计算本地网页PageRank值 for (i = rank; i < N; i += size) { pr_new[i] = 0.0; for (j = 0; j < N; j++) { if (link[j][i] > 0.0) { pr_new[i] += pr[j] / link[j][i]; } } pr_new[i] *= DAMPING_FACTOR; pr_new[i] += (1.0 - DAMPING_FACTOR) / N; } // 计算本地网页PageRank值的差值 diff = 0.0; for (i = rank; i < N; i += size) { diff += fabs(pr_new[i] - pr[i]); } // 求所有进程的差值和 MPI_Allreduce(&diff, &diff_all, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); // 将本地网页PageRank值的差值发送给其他进程 for (i = 0; i < size; i++) { if (i != rank) { MPI_Send(&diff, 1, MPI_DOUBLE, i, 0, MPI_COMM_WORLD); } } // 接收其他进程的网页PageRank值 for (i = 0; i < size; i++) { if (i != rank) { MPI_Recv(pr_all, N, MPI_DOUBLE, i, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE); for (j = 0; j < N; j++) { pr_new[j] += pr_all[j]; } } } // 更新网页PageRank值 for (i = rank; i < N; i += size) { pr[i] = pr_new[i]; } // 求所有进程的网页PageRank值 MPI_Allgather(pr, N, MPI_DOUBLE, pr_all, N, MPI_DOUBLE, MPI_COMM_WORLD); // 增加迭代次数 iter++; } while (diff_all > EPSILON); // 输出网页PageRank值 if (rank == 0) { printf("Iterations: %d\n", iter); for (i = 0; i < N; i++) { printf("PR[%d] = %f\n", i, pr_all[i]); } } // 结束MPI MPI_Finalize(); return 0; } ``` 以上代码使用MPI实现了PageRank算法的并行计算，其中每个进程负责计算部分网页的PageRank值，并将计算结果传递给其他进程进行汇总。同时，代码中还实现了收敛精度的控制和迭代次数的统计，以及网页链接矩阵的随机初始化。

“相关推荐”对你有帮助么？

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值