并行程序设计-实验5.MPI编程练习

并行程序设计 S5课上 实验报告
并行程序设计
发布日期:   2022-11-28
文章字数:   3.5k

MPI 编程练习实验报告

实验内容

  • 实现第5章课件中的梯形积分法的MPI编程熟悉并掌握MPI编程方法，规模自行设定，可探讨不同规模对不同实现方式的影响。
  • 对于课件中“多个数组排序”的任务不均衡案例进行MPI编程实现，规模可自己设定、调整。
  • 附加：实现高斯消去法解线性方程组的MPI编程，与SSE（或AVX）编程结合，并与Pthread、OpenMP（结合SSE或AVX）版本对比，规模自己设定。

实验一：梯形积分

问题描述

本实验分别使用MPI、Pthread和OpenMP三种方法，实现了梯形积分法。并通过调整梯形积分法划分成小梯形的个数规模，来比较不同编程方式的异同。

具体案例为：对于下图中给出的梯形积分法，实现并行编程。每个线程或进程计算a到b区间中的某一段的梯形面积，最后求取全局和得到结果。

![image-20230102165554905](/Users/zhangxiaoni/Library/Application Support/typora-user-images/image-20230102165554905.png)

程序设计与实现

MPI实现

实验共用到了3台主机的进程进行计算。每个进程根据自己的my_rank进程号得到自身的计算任务，完成局部和的计算后，（除0号进程外）使用MPI_Send将结果发送至0号进程；0号进程使用MPI_Recv阻塞式地接受其他进程传回的结果，并计算全局和。同时，0号进程还完成计时、输出结果的任务。具体的实现代码如下：

#include <stdio.h>
#include <algorithm>
#include<cstdlib>
#include<ctime>
#include<mpi.h>

//待积分的f(x)
double f(double x) {
    return 1;
}

//计算局部和
double Trap(double a, double b, int count, double h) {
    double my_result, x;
    int i;
    my_result = (f(a) + f(b)) / 2.0;
    //局部求和，避免过多通信
    for (i = 1; i < count; i++) {
        x = a + i * h;
        my_result += f(x);
    }
    my_result *= h;
    return my_result;
}

//主函数
int main()
{
    int my_rank, comm_sz, n = 1000;
    double a = 0, b = 5000;
    double h = (b - a) / n;
    double global_result;
    clock_t start, end;
    double time;

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

    //0号进程开始计时
    if (my_rank == 0) {
        start = clock();
    }
    int local_n = n / comm_sz;
    int local_a = a + my_rank * local_n * h;
    int local_b = local_a + local_n * h;
    //将未除尽的任务数分配给最后一个进程
    if (my_rank == comm_sz - 1) {
        local_n = n - (comm_sz - 1) * local_n;
        local_b = b;
    }
    //每个进程计算自身分配的任务的局部和
    double local_result = Trap(local_a, local_b, local_n, h);
    //其他进程将局部和发送给0号进程
    if (my_rank != 0) {
        MPI_Send(&local_result, 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
    }
    //0号进程收集局部和并求全局和
    else {
        global_result = local_result;
        for (int source = 1; source < comm_sz; source++)
        {
            MPI_Recv(&local_result, 1, MPI_DOUBLE, source, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
            global_result += local_result;
        }
        //完成计算，0号进程停止计时
        end = clock();
        time = (end - start) / CLOCKS_PER_SEC;
    }
    //0号进程打印出结果
    if (my_rank == 0) {
        printf("划分成小梯形的块数: %d\n", n);
        printf("计算结果是: %.15e\n", global_result);
        printf("总共耗时: %f\n\nms", time);
    }

    MPI_Finalize();
    return 0;
}

Pthread实现

对于每个线程，先求局部和最后全局求和，避免过多通信。同时注意给全局和加锁保证结果的正确性。代码实现如下：

void pthread_trap(double a,double b,int n,double*global_result_p,int thread_count){
    pthread_param params[thread_count];
    pthread_t threads[thread_count];
    pthread_mutex_t amutex;
    pthread_mutex_init(&amutex,NULL);
    for(int i=0;i<thread_count;i++){
        //为每个线程参数赋值
        params[i].a=a;
        params[i].b=b;
        params[i].n=n;
        params[i].global_result_p=global_result_p;
        params[i].thread_count=thread_count;
        params[i].my_rank=i;
        params[i].amutex=&amutex;
        //创建线程
        pthread_create(&threads[i],NULL,pthread_trap,(void*)&params[i]);
    }
    //销毁线程
    for(int i=0;i<thread_count;i++){
        pthread_join(threads[i],NULL);
    }
    //销毁互斥量锁
    pthread_mutex_destroy(&amutex);
}

其中每个线程执行的函数如下：

void* pthread_trap(void* p){
    //恢复获取参数
    pthread_param* params=(pthread_param*) p;
    double h,x,my_result;
    double local_a,local_b;
    int local_n;
    int my_rank=params->my_rank;
    int thread_count=params->thread_count;
    h=(params->b-params->a)/params->n;
    local_n=params->n/thread_count;
    local_a=params->a+my_rank*local_n*h;
    local_b=local_a+local_n*h;
    my_result=(f(local_a)+f(local_b))/2.0;
    //先局部求和，避免过多通信
    for(int i=1;i<=local_n-1;i++){
        x=local_a+i*h;
        my_result+=f(x);
    }
    my_result=my_result*h;
    //使用互斥量 amutex 上锁保证全局求和正确性
    pthread_mutex_lock(params->amutex);
    *params->global_result_p+=my_result;
    pthread_mutex_unlock(params->amutex);
}

OpenMP实现

利用OpenMP的相关方法实现多线程并行编程，对于每个线程，先求局部和最后全局求和，避免过多通信。同时注意给全局和加锁保证结果的正确性。代码实现如下：

#pragma omp parallel num_threads(thread_count)
trap(a,b,n,&global_result);

void trap(double a,double b,int n,double*global_result_p){
    double h,x,my_result;
    double local_a,local_b;
    int local_n;
    int my_rank=omp_get_thread_num();
    int thread_count=omp_get_num_threads();
    h=(b-a)/n;
    local_n=n/thread_count;
    local_a=a+my_rank*local_n*h;
    local_b=local_a+local_n*h;
    my_result=(f(local_a)+f(local_b))/2.0;
    //先局部求和，避免过多通信
    for(int i=1;i<=local_n-1;i++){
        x=local_a+i*h;
        my_result+=f(x);
    }
    my_result=my_result*h;
    //用临界区保证全局求和正确性
    #pragma omp critical
        *global_result_p+=my_result;
}

实验设计

不同任务规模下三种方法执行时间对比

待计算函数为 f(x) = 1；待计算区间为0到5000。将积分区域分别划分为100，1000，10000，100000个小梯形，比较3种方法在这些不同的任务规模下执行效率的差异。

计时方法

  • MPI版本使用的计时方法如下：

    start = clock();
end = clock();
time = (end - start) / CLOCKS_PER_SEC;

  • Pthread和OpenMP版本继续沿用之前实验中用到的计时方法

    QueryPerformanceFrequency((LARGE_INTEGER *)&freq);
QueryPerformanceCounter((LARGE_INTEGER *)&head);          QueryPerformanceCounter((LARGE_INTEGER *)&tail);

不同规模下3种方法的结果展示与对比

程序输出结果

待计算函数为 f(x) = 1；待计算区间为0到5000

  1. 将MPI的运行环境及程序部署到华为鲲鹏云服务器上，得到示例输出如下：

    ![image-20230102211323833](/Users/zhangxiaoni/Library/Application Support/typora-user-images/image-20230102211323833.png)

  2. 编辑器为codeblocks，编译器为mingw64，得到Pthread版本和OpenMP版本的示例输出如下：

    ![image-20230102193907704](/Users/zhangxiaoni/Library/Application Support/typora-user-images/image-20230102193907704.png)

    从上述截图中看出，三种实现方式均按照预期效果输出，程序正确。

结果分析

按照上述的实验设计执行程序，得到经过分析处理后的图表如下：

100 1000 10000 100000
MPI 4ms 5ms 11ms 13ms
Pthread 13ms 15ms 16ms 19ms
OpenMP 7ms 9ms 12ms 14ms

![image-20230102211927911](/Users/zhangxiaoni/Library/Application Support/typora-user-images/image-20230102211927911.png)

根据上述图表可以看出，几种实现方式的性能差异并不明显，可能是因为该计算任务本身复杂度较低，三种方式均能快速完成，故无法很好地展现出不同方法的性能差异。单从本实验的表现来看，MPI的执行效率相对来说更高一点，可能是因为多台主机多个进程多个CPU确实会更快一些。当然，这种性能上的优秀表现仅在平均值中体现（上述图表中的数据是多次实验取平均后得到的结果），在某些情况下，MPI方法的执行时间也会出现较长的现象，猜想可能是因为使用云服务器传输数据的过程耗时较长导致。

实验总结

通过本实验，我熟悉并掌握MPI编程方法，在探讨不同规模对不同实现方式的影响，我了解到了多台主机共同完成同一计算任务的思想以及在本实验中的优秀表现。不过本实验中MPI、Pthread、OpenMP三种方法的性能差异并不是特别明显，这也是本实验可以继续改进优化的一个方向。

实验二：数组排序的任务分配不均衡案例复现

问题描述

对ARR_NUM个长度为ARR_LEN的一维数组进行排序，使用MPI编程，每个进程处理一部分数组排序的任务。若数组分布不均衡，如前二分之一的数组全部升序，后二分之一的数组全部逆序，则每个进程分得的任务负载可能也会不均衡。

程序设计与实现

实验共用到了3台主机的进程完成数组排序。根据进程号连续划分任务，每个进程完成排序之后，将自己执行任务的耗时传送给0号进程，最终0号进程接受各个进程完成排序任务的时间并输出。具体的实现代码如下：

#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <vector>
#include <algorithm>
#include <mpi.h>
using namespace std;

//数组数量
const int ARR_NUM = 3000;
//每个数组的长度
const int ARR_LEN = 300;
//进程数量/主机数量
const int PROCESS_NUM = 3;
//每个进程处理的数组数量
const int seg = ARR_NUM / PROCESS_NUM;
vector<int> arr[ARR_NUM];

// 初始化待排序数组，使得
// 第一段：完全逆序
// 第二段：1/2逆序，1/2升序
// 第三段：1/4逆序，3/4升序
// 第四段：完全升序
void init(void) {
    int ratio;
    srand(unsigned(time(nullptr)));
    for (int i = 0; i < ARR_NUM; i++) {
        arr[i].resize(ARR_LEN);
        if (i < seg)ratio = 0;
        else if (i < seg * 2)ratio = 32;
        else ratio = 128;
        if ((rand() & 127) < ratio) {
            for (int j = 0; j < ARR_LEN; j++) {
                arr[i][j] = ARR_LEN - j;
            }
        }
        else {
            for (int j = 0; j < ARR_LEN; j++) {
                arr[i][j] = j;
            }
        }
    }
}

//主函数
int main() {
    int my_rank, comm_sz;
    clock_t start, end;
    double time;
    //初始化数组
    init();

    MPI_Init(NULL, NULL);
    MPI_Comm_rank(MPI_COMM_WORLD, &my_rank);
    MPI_Comm_size(MPI_COMM_WORLD, &comm_sz);
    //根据进程号分配任务
    int startTask = my_rank * seg;
    int endTask = startTask + seg;
    //排序
    start = clock();
    for (int i = startTask; i < endTask; i++) {
        sort(arr[i].begin(), arr[i].end());
    }
    end = clock();
    time = (end - start);
    //其他进程将本进程的运行时间发送给0号进程
    if (my_rank != 0) {
        MPI_Send(&time, 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
    }
    //0号进程打印各个进程的运行时间
    else if (my_rank == 0) {
        printf("各进程运行时间如下:\n");
        printf("0号进程: %fms\n", time);
        for (int source = 1; source < comm_sz; source++) {
            MPI_Recv(&time, 1, MPI_DOUBLE, source, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
            printf("%d号进程: %fms\n", source, time);
        }
    }
    MPI_Finalize();
    return 0;
}

实验设计

  • 使用init函数生成分布不均的数组；
  • 使用与实验一（梯形积分法）中相同的计时方法。
  • 改变数组规模，分别处理500，1000，1500，2000，2500，3000个长度为300的数组，比较在不同规模下各个进程的时间差异，从而理解任务负载不均衡对效率的影响。

结果展示与对比分析

程序输出结果

将MPI的运行环境及程序部署到华为鲲鹏云服务器上，得到示例输出如下：

注：本实验给每台主机分配了1个进程（在config文件中设置），总共有3个进程共同执行排序任务，故执行命令时为“3”

![截屏2023-01-02 21.32.52](/Users/zhangxiaoni/Desktop/截屏2023-01-02 21.32.52.png)

图表结果分析

经过分析处理后得到的图表结果展示如下：

500 1000 1500 2000 2500 3000
0号进程 5240ms 10564ms 15789ms 21170ms 26619ms 31875ms
1号进程 5095ms 10224ms 15539ms 20340ms 25586ms 30774ms
2号进程 4712ms 9346ms 13974ms 18574ms 23425ms 28067ms

![image-20230102214239527](/Users/zhangxiaoni/Library/Application Support/typora-user-images/image-20230102214239527.png)

从上述图表中，我们大致可以看出，由于数据分布不均衡，各个进程的执行时间存在差异，具体表现为：0号、1号、2号进程的执行效率依次增大，这与他们各自的任务负载是保持一致的。在实际情况中，这些负载任务复杂度较高、执行时间较长的进程，会拖慢任务解决的整体时长，这是我们不希望看到的。

实验总结

​ 通过本实验，我了解到数据分布均衡对于多进程、多线程程序的重要性。如果任务分布不均，将会出现执行效率较低的进程拖慢任务解决的整体效率的情况。因此，并行编程时控制任务规模、控制任务分配的粗细粒度、控制进程或线程的数量，使得程序在数据分布不均衡的情况下也能有较为优秀的表现就显得很重要。当然，也可以在程序执行任务之前，对数据进行处理，达到各个进程分配的数据分布都较为均衡的效果。

附加实验：高斯消去法解线性方程组

问题描述

对于给定的线性方程组Ax=b，使用MPI编程，完成对于线性方程组的高斯消元。

程序设计与实现

每个进程完成部分行的消元计算任务，然后使用barrier功能，使得所有进程均完成第k行的消元后再计算下一行的消元。程序实现的具体代码如下：

#include <stdio.h>
#include<iostream>
#include <algorithm>
#include<cstdlib>
#include<ctime>
#include<mpi.h>
using namespace std;

const int n = 900;
const int NUM_PROC = 3;
float a[n][n], b[n], c[n];

void init() {
    for (int i = 0; i < n; i++) {
        for (int j = 0; j < n; j++) {
            a[i][j] = rand() % (1000 - 1) + 1;
        }
        b[i] = rand() % (1000 - 1) + 1;
    }
}

int main() {
    int my_rank, comm_sz, left, right;
    int source;
    clock_t start, finish;
    double sec;
    init();

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

    start = clock();
    for (int k = 0; k < n - 1; k++) {
        int seg = (n - (k + 1)) / NUM_PROC;
        left = my_rank * seg;
        right = (my_rank + 1) * seg;
        for (int i = left; i < right; i++) {
            c[i] = a[i][k] / a[k][k];
        }
        for (int i = left; i < right; i++) {
            for (int j = 0; j < n; j++) {
                a[i][j] = a[i][j] - c[i] * a[k][j];
            }
            b[i] = b[i] - c[i] * b[k];
        }
        MPI_Barrier(MPI_COMM_WORLD);
    }
    finish = clock();
    sec = finish - start;

    if (my_rank != 0) {
        MPI_Send(&sec, 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
    }
    else {
        printf("0号进程: %f\n", sec);
        for (source = 1; source < comm_sz; source++) {
            MPI_Recv(&sec, 1, MPI_DOUBLE, source, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
            printf("%d号进程: %f\n", source, sec);
        }
    }
    MPI_Finalize();
    return 0;
}

实验总结

该部分尝试对高斯消元法进行了的MPI的编程实现，其中也遇到了很多困难，比如数据依赖等问题。经过反复尝试，我实现了高斯消元MPI的解法，使得我对于MPI编程中需要注意的内存共享的问题有了更深刻的理解。

文章作者: zxn
文章链接: http://zxnstc.github.io/2022/11/28/bing-xing/mpi-bian-cheng-lian-xi-shi-yan-bao-gao/
版权声明: 本博客所有文章除特別声明外，均采用 CC BY 4.0 许可协议。转载请注明来源 zxn !
并行程序设计 S5课上 实验报告
赏
  目录