天津大学并行计算实验三：多进程计算矩阵幂

华年似锦

已于 2024-05-17 10:28:01 修改

阅读量577

点赞数 5

文章标签：矩阵

于 2024-04-22 19:16:41 首次发布

本文链接：https://blog.csdn.net/Fun_Total/article/details/138087894

版权

往期实验的传送门：
实验二：多线程计算矩阵幂
 实验一：多线程计算正弦值

前言

大的要来咯，当初就是写这个实验写得绝望了，才激发我写下这些攻略，网上找不到一点，GPT也写不了一点。考虑到一些潜在的风险，就不提供完整的报告了，就下面给出一些我的参考代码供大家更好地学习并行计算。

免责声明：注意：代码等仅供同学们加深对课程知识点的理解，严禁抄袭，要是查重（不确定有没有）什么的被老师抓到了后果自负

环境要求

本次实验因为不是用多线程，而是用MPI多进程，所以环境上有所改变。注意参考最新的实验指导书，上面的指令跟智慧树上的截图有所不同（如果程序能跑出来结果，但是下面有跟着一堆段错误报错，大概率是这个环境问题）
环境要求，编译前先输入

module load openmpi/4.1.4-mpi-x-gcc9.3.0

就是让虚拟机配置一下mpi的环境，然后才有办法编译用mpi库的程序。

算法介绍

MPI实现矩阵乘法，左边是我参考的代码，我是看这篇博客里第二条的矩阵乘法代码，将其改编成矩阵幂计算的，同学们也可以试一下，自己从这篇博客里面的代码改写矩阵幂计算。因为这篇博客用的是c++，所以此次实验，我后面代码也是用的c++风格实现，编译的时候指令记得后缀不是.c而是.cpp

然后正式开始介绍一下此次的算法，这次跟前一次的矩阵幂在思想上类似，将要计算的矩阵，按照行数划分给若干个进程，每个进程计算一部分，然后每次矩阵乘法计算后，都需要同步一下各个进程的结果，然后再进行下一次的计算。
由于上课没听，所以不知道如何给各个进程分发任务以及同步各个进程，参考那篇博客的主要感悟也就是知道了MPI_Gather()等函数如何使用。

算法流程图：

先上代码，然后讲解在代码注释里面：

#include <iostream>
using namespace std;
#include <mpi.h>

int main(int argc, char** argv) {
    int my_rank;
    int num_procs;
    int size = atoi(argv[1]); //矩阵规模
    int cnt = atoi(argv[2]); // 计算幂次

    MPI_Init(NULL, NULL);
    MPI_Comm_rank(MPI_COMM_WORLD, &my_rank);
    MPI_Comm_size(MPI_COMM_WORLD, &num_procs);

    int line = size / num_procs;  // 每个进程计算的行数
    int* local_a = new int[line * size];
    int* b = new int[size * size];
    int* ans = new int[line * size];
    int* a = new int[size * size];
    int* c = new int[size * size];  // 用来缓冲区接受数据

    if (my_rank == 0) {
        for (int i = 0; i < size; i++) {  // a是计算矩阵，b是结果矩阵
            for (int j = 0; j < size; j++) {
                a[i * size + j] = 1; //这里可以随机化计算矩阵
                b[i * size + j] = 0;
                if (i == j)
                    b[i * size + j]  = 1; // 这个是因为0次幂是单位矩阵
            }
        }

        for (int i = 1; i <= cnt; i++) {
            MPI_Scatter(a, line * size, MPI_INT, local_a, line * size, MPI_INT,
                        0, MPI_COMM_WORLD); // 将待计算矩阵a分发个各个进程

            MPI_Bcast(b, size * size, MPI_INT, 0, MPI_COMM_WORLD); // 将矩阵b广播给各个进程

            for (int i = 0; i < line; i++) {
                for (int j = 0; j < size; j++) {
                    int temp = 0;
                    for (int k = 0; k < size; k++)
                        temp += a[i * size + k] * b[k * size + j];
                    ans[i * size + j] = temp; // 主进程自己要计算的部分
                }
            }

            MPI_Gather(ans, line * size, MPI_INT, c, line * size, MPI_INT, 0,
                       MPI_COMM_WORLD); // 主进程汇给结果
            for (int i = 0; i < size; i++)
                for (int j = 0; j < size; j++)
                    b[i * size + j] = c[i * size + j]; // 交换矩阵，保存当前轮的计算结果
        }
        cout << "ans:" << endl;
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                cout << b[i * size + j] << " ";
            }
            cout << endl;
        }
    } else { // 是其他进程
        while (cnt--) {
            int* buffer = new int[size * line]; // 接受主进程发送数据的缓冲区
            MPI_Scatter(a, line * size, MPI_INT, buffer, line * size, MPI_INT,
                        0, MPI_COMM_WORLD);
            MPI_Bcast(b, size * size, MPI_INT, 0, MPI_COMM_WORLD);

            for (int i = 0; i < line; i++) {
                for (int j = 0; j < size; j++) {
                    int temp = 0;
                    for (int k = 0; k < size; k++)
                        temp += buffer[i * size + k] * b[k * size + j];
                    ans[i * size + j] = temp;
                }
            }
            MPI_Gather(ans, line * size, MPI_INT, c + my_rank * line * size,
                       line * size, MPI_INT, 0, MPI_COMM_WORLD); // 其他进程发送给主进程结果，并且按照对应的次序发送位置 c + my_rank * line * size
            delete[] buffer;
        }
    }

    delete[] a, local_a, b, ans, c;

    MPI_Finalize();
    return 0;
}

这里面精髓的地方，感觉在MPI_Gather()参数的调整，原来的博客里面好像是乱的，直接cout出去了。但是因为矩阵幂有多轮的矩阵乘法，所以必须让每个进程发送到对应的位置，通过my_rank控制发送的位置。
然后因为总的每个进程都要计算 $c n t$ 次，所以不是主进程的else那里，也要 while(cnt–)，重复cnt次的计算才能退出。

写好程序后，注意现在的编译命令改为

mpic++ -o bing.o bing.cpp

如果前面没有输那句 module load配一下mpi环境的话编译会报错。

注意

因为这里面用的函数是MPI_Gather()，好像是这个函数必须每个进程发送的数据大小相同，所以矩阵划分的时候必须得整除，所以就规定输入的矩阵规模必须得是线程数8的整数倍，在脚本设置数据范围的时候要注意一下。

因为这次懒得写串行了，我就直接用单核单线程来替代串行程序，脚本如下：

#!/bin/bash
module load openmpi/4.1.4-mpi-x-gcc9.3.0
time yhrun -p thcp1 -N 2 -n 8 bing.o 80 500 &> run2.log
time yhrun -p thcp1 -N 2 -n 8 bing.o 240 100 &> run2.log
time yhrun -p thcp1 -N 2 -n 8 bing.o 240 200 &> run2.log
time yhrun -p thcp1 -N 2 -n 8 bing.o 496 100 &> run2.log
time yhrun -p thcp1 -N 2 -n 8 bing.o 496 200 &> run2.log
time yhrun -p thcp1 -N 2 -n 8 bing.o 960 200 &> run2.log
time yhrun -p thcp1 -N 2 -n 8 bing.o 1600 10 &> run2.log

time yhrun -p thcp1 -N 1 -n 1 bing.o 80 500 &> run2.log
time yhrun -p thcp1 -N 1 -n 1 bing.o 240 100 &> run2.log
time yhrun -p thcp1 -N 1 -n 1 bing.o 240 200 &> run2.log
time yhrun -p thcp1 -N 1 -n 1 bing.o 496 100 &> run2.log
time yhrun -p thcp1 -N 1 -n 1 bing.o 496 200 &> run2.log
time yhrun -p thcp1 -N 1 -n 1 bing.o 960 200 &> run2.log
time yhrun -p thcp1 -N 1 -n 1 bing.o 1600 10 &> run2.log