项目地址：https://github.com/snjl/c.openmp.git

开发环境：Clion或visual studio，vs需要在项目里开启openmp，具体可以搜索。

时间统计使用omp_get_wtime，是omp包里的函数，如下：

1	double start = omp_get_wtime();

如果使用time的函数，clock_t是记录cpu的滴答数的，并行时多个进程同时计算，cpu计时数成倍增加。

题目

编写完整OpenMP程序：要求在程序中指定用4个线程执行；定义A、B数组，实现初始化整数数组A[1000][1000]，然后计算B[i][j]=（A[i-1][j]+A[i+1][j]+A[i][j-1]+A[i][j+1])/4（其中0<i,j<1000）。

代码

头文件：

#include <iostream>
#include "time.h"
#include<omp.h>

using namespace std;

初始化：

#define MATRIXLENGTH 10000
#define TIME 10

double matrix[MATRIXLENGTH][MATRIXLENGTH] = { 0.0 };

double newMatrix[MATRIXLENGTH][MATRIXLENGTH] = { 0.0 };

核心计算函数：

void calculate(double matrix[MATRIXLENGTH][MATRIXLENGTH], double newMatrix[MATRIXLENGTH][MATRIXLENGTH], int height, int length, int ompNum) {
  omp_set_num_threads(ompNum);
#pragma omp parallel
  {
# pragma omp for
    for (int i = 0; i < height; i++) {
      for (int j = 0; j < length; j++) {
        if (i > 0) {
          newMatrix[i][j] += matrix[i - 1][j];
        }
        if (j > 0) {
          newMatrix[i][j] += matrix[i][j - 1];
        }
        if (i < height - 1) {
          newMatrix[i][j] += matrix[i + 1][j];
        }
        if (j < length - 1) {
          newMatrix[i][j] += matrix[i][j + 1];
        }
        newMatrix[i][j] /= 4.0;

      }
    }
  }
};

关键点：

没有考虑角、边的2、3邻域情况（如果考虑，需要一个变量值，需要每个线程固定，不然会混乱）；
用了很多if，可以考虑做0的填充，就不需要考虑if了。

优化

方法1：如果不考虑边角，则不需要考虑2邻域或3邻域，可以去除私有变量num，这样能减少一部分计算量，由此考虑去除私有变量num来进行计算并且对比效率

void calculate(double matrix[MATRIXLENGTH][MATRIXLENGTH], double newMatrix[MATRIXLENGTH][MATRIXLENGTH], int height, int length, int ompNum) {
    omp_set_num_threads(ompNum);
#pragma omp parallel
    {
#pragma omp for
        for (int i = 0; i < height; i++) {
            for (int j = 0; j < length; j++) {
                if (i > 0) {newMatrix[i][j] += matrix[i - 1][j];}
                if (j > 0) {newMatrix[i][j] += matrix[i][j - 1];}
                if (i < height - 1) {newMatrix[i][j] += matrix[i + 1][j];}
                if (j < length - 1) {newMatrix[i][j] += matrix[i][j + 1];}
                newMatrix[i][j] /= 4;
            }
        }
    }
};

方法2：由于会判断是否出界，所以使用了4个if，如果加上一圈padding，全部填充0，则可以去掉所有if语句，同时也去除了计算2邻域或3邻域的问题，但是会占用较高内存

#define MATRIXLENGTH 1002
void calculate(double matrix[MATRIXLENGTH][MATRIXLENGTH], double newMatrix[MATRIXLENGTH][MATRIXLENGTH], int height, int length, int ompNum) {
    omp_set_num_threads(ompNum);
#pragma omp parallel
    {
#pragma omp for
        for (int i = 1; i < height - 1; i++) {
            for (int j = 1; j < length - 1; j++) {
                newMatrix[i][j] = (matrix[i - 1][j] + matrix[i + 1][j] + matrix[i][j - 1] + matrix[i][j + 1]) / 4;
            }
        }
    }
};

平均计算时间函数：

double getAverage(const double duration[TIME]) {
  double average = 0.0;
  for (int i = 0; i < TIME; i++) {
    average += duration[i];
  }
  return average / TIME;
}

主函数：

int main() {
    double duration[TIME] = { 0.0 };
    for (int d = 0; d < TIME; d++) {
        double start = omp_get_wtime();
        calculate(matrix, newMatrix, MATRIXLENGTH, MATRIXLENGTH, 4);
        double end = omp_get_wtime();
        duration[d] = end - start;
        printf("%f s  ", duration[d]);
    }
    cout << "average time is : " << getAverage(duration) << endl;
    return 0;
}

输出结果：

1 2	1.328000 s 0.843000 s 0.569000 s 0.552000 s 0.529000 s 0.582000 s 0.630000 s 0.671000 s 0.661000 s 0.584000 s average time is : 0.6949

方法3：考虑将第三个版本的矩阵拉直成一个向量，即1002*1002的向量，这样在第二个版本的基础上减少一次for循环（不考虑矩阵拉直和还原的时间），代码如下：

#define MATRIXLENGTH 1002
#define TIME 100
double matrix[MATRIXLENGTH*MATRIXLENGTH] = { 0.0 };
double newMatrix[MATRIXLENGTH*MATRIXLENGTH] = { 0.0 };
void calculate(double matrix[MATRIXLENGTH*MATRIXLENGTH], double newMatrix[MATRIXLENGTH*MATRIXLENGTH], int ompNum) {
    omp_set_num_threads(ompNum);
#pragma omp parallel
    {
#pragma omp for
        for (int i = MATRIXLENGTH; i < MATRIXLENGTH*(MATRIXLENGTH - 1); i++) {
            newMatrix[i] = (matrix[i - MATRIXLENGTH] + matrix[i + MATRIXLENGTH] + matrix[i - 1] + matrix[i + 1]) / 4;
        }
    }
};

注意

在并行里未注明则是共享变量。如果有需要可以设置成私有变量（如优化方法1所示）。

1 2	int num = 0; #pragma omp for private(num)

OpenMP对于嵌套循环应该添加多少个parallel for

一个原则是：应该尽量少的使用parallelfor, 因为parallel for也需要时间开销。即:

（1）如果外层循环次数远远小于内层循环次数，内层循环较多时，将parallel for加在内层循环。
（2）否则将parallel for 加在最外层循环，一般情况都是这样。二者在实际情况可对比性能进行选择。