利用华为鲲鹏云服务器进行MPI编程

1. 梯形积分法

1.1 问题描述

实现第 5 章课件中的梯形积分法的 MPI 编程熟悉并掌握 MPI 编程方法，规模自行设定，可探讨 不同规模对不同实现方式的影响。

1.2 算法设计与实现

    MPI_Init(&argc, &argv);

    int rank;
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);

    int size;
    MPI_Comm_size(MPI_COMM_WORLD, &size);

    double a = 0.0;  // 积分区间的左端点
    double b = 8.0;  // 积分区间的右端点
    int n = 1000;    // 小区间的数量

    double h = (b - a) / n;  // 小区间的长度

    int local_n = n / size;  // 每个进程计算的小区间数量
    double local_a = a + rank * local_n * h;  // 每个进程计算的积分区间的左端点
    double local_b = local_a + local_n * h;   // 每个进程计算的积分区间的右端点

    double local_sum = (f(local_a) + f(local_b)) / 2.0;
    for (int i = 1; i < local_n; ++i) {
        double x = local_a + i * h;
        local_sum += f(x);
    }
    local_sum *= h;

    double total_sum;
    MPI_Reduce(&local_sum, &total_sum, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD);

    if (rank == 0) {
        printf("With n = %d trapezoids, our estimate\n", n);
        printf("The integral from %f to %f is %f", a, b, total_sum);
    }

    MPI_Finalize();
    return 0;

1.3 实验结果

本实验中，我使用了 MPI 方法来实现梯形积分法，梯形的积分区间为 [3.0, 8.0]，被积函数为 x 2 ， 积分区间分别被分为为 1000￿5000￿10000 个小区间，线程数为 6。实验结果如下图所示：

1.4 实验分析

从上图可以看出，随着划分的小区间规模增大，积分的值越来越接近真实值，但是当划分的小区间 规模增大到一定程度后，积分的值就不再变化，这是因为当划分的小区间规模增大到一定程度后，计 算机的精度就不够了，所以积分的值就不再变化了。

2 数组排序

2.1 问题描述

对于课件中“多个数组排序”的任务不均衡案例进行 MPI 编程实现，规模可自己设定、调整。

2.2 算法设计与实现

首先在数组初始化时，对于中“多个数组排序”的任务不均衡案例进行复现

    int ratio;
    srand(unsigned(time(NULL)));
    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 if (i < seg * 3) {
            ratio = 64;
        } 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;
            }
        }
    }

 init();
    //记录开始时间
    clock_t start, finish;
    double duration;
    start = clock();
    int current_task = 0;//当前的任务
    MPI_Init(NULL, NULL);
    int rank, thread_num;
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
    MPI_Comm_size(MPI_COMM_WORLD, &thread_num);
    MPI_Status status;
    int ready;
    bool done = false;
    if (rank == 0) {
        while (current_task < ARR_NUM) {
            //接受任何一个线程的状态
            MPI_Recv(&ready, 1, MPI_INT, MPI_ANY_SOURCE, MPI_ANY_TAG, MPI_COMM_WORLD, &status);
            MPI_Send(&current_task, 1, MPI_INT, status.MPI_SOURCE, 0, MPI_COMM_WORLD);
            current_task += seg;
        }
        done = true;
        //记录结束时间
        finish = clock();
        duration = (double) (finish - start) / CLOCKS_PER_SEC;
        cout << "ARR_NUM:" << ARR_NUM << endl;
        cout << "ARR_LEN:" << ARR_LEN << endl;
        cout << "time: " << duration << "ms" << endl;
        //唤醒所有线程
        for (int i = 1; i < thread_num; ++i) {
            MPI_Recv(&ready, 1, MPI_INT, MPI_ANY_SOURCE, MPI_ANY_TAG, MPI_COMM_WORLD, &status);
            MPI_Send(&current_task, 1, MPI_INT, status.MPI_SOURCE, 1, MPI_COMM_WORLD);
        }

    } else {
        while (!done) {
            int begin;
            MPI_Send(&rank, 1, MPI_INT, 0, 0, MPI_COMM_WORLD);
            MPI_Recv(&begin, 1, MPI_INT, 0, MPI_ANY_TAG, MPI_COMM_WORLD, &status);
            if (status.MPI_TAG == 1) {
                break;
            }
            doTask(begin);
        }
    }
    MPI_Finalize();
    return 0;

2.3 实验结果

2.4 实验分析

从图中可以看出，当数组长度为 8000 时，数组个数越多，排序时间越长，这是因为数组长度固定， 数组个数越多，每个数组的排序时间越长，而每个数组的排序时间是并行的，所以数组个数越多，排 序时间越长。可以看出，排序时间的增长趋势是线性的。

3 高斯消元法

3.1 问题描述

实现高斯消去法解线性方程组的 MPI 编程，与 SSE（或 AVX）编程结合，并与 Pthread、OpenMP （结合 SSE 或 AVX）版本对比，规模自己设定。

3.2 算法设计与实现

float **matrix_const = new float *[N];
    float *x_const = new float[N];
    matrix = new float *[N];
    for (int i = 0; i < N; i++) {
        matrix_const[i] = new float[N];
        matrix[i] = new float[N];
    }
    x = new float[N];

    init(N, matrix, x);
    int rank, thread_num;
    gettimeofday(&startTime, NULL);
    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
    MPI_Comm_size(MPI_COMM_WORLD, &thread_num);
    MPI_Status status;
    int ready;
    bool done = false;
    if (rank == 0) {
        printf("size : %d\n", N);
        for (line = 0; line < N - 1; ++line) {
            next_task = line + 1;
            seg = (N - next_task) / (thread_num - 1) + 1;
            for (int i = 1; i < thread_num; i++) {
                int task = (i - 1) * seg + next_task;
                MPI_Send(&task, 1, MPI_INT, i , 0, MPI_COMM_WORLD);
            }
            //等待所有线程
            for (int j = 1; j < thread_num; ++j) {
                MPI_Recv(&ready, 1, MPI_INT, MPI_ANY_SOURCE, MPI_ANY_TAG, MPI_COMM_WORLD, &status);
            }
        }
        printf("all work done!!!!!\n");
        done = true;
        gettimeofday(&stopTime, NULL);
        double trans_mul_time =
                (stopTime.tv_sec - startTime.tv_sec) * 1000 + (stopTime.tv_usec - startTime.tv_usec) * 0.001;
        printf("time: %lf ms\n", trans_mul_time);
        //唤醒所有线程
        for (int i = 1; i < thread_num; i++) {
//            MPI_Recv(&ready, 1, MPI_INT, MPI_ANY_SOURCE, MPI_ANY_TAG, MPI_COMM_WORLD, &status);
            MPI_Send(&done, 1, MPI_INT, status.MPI_SOURCE, 1, MPI_COMM_WORLD);
        }
    } else {
        while (!done) {
            int task;
            MPI_Recv(&task, 1, MPI_INT, 0, MPI_ANY_TAG, MPI_COMM_WORLD, &status);
            if(done) break;
            int min = task + seg < N ? task + seg : N;
            for (int i = task; i < min; ++i) {
                GAUSS_EXITNCTION(line, i);
            }
            MPI_Send(&ready, 1, MPI_INT, 0, 0, MPI_COMM_WORLD);
        }
    }
    MPI_Finalize();
    return 0;

3.3 实验结果