mpi+openmp矩阵乘法代码

MPI (Message Passing Interface) 和 OpenMP 是两个并行计算模型，常用于大规模分布式内存环境下的并行程序设计。在矩阵乘法中，可以结合它们来加速计算。下面是一个简单的示例，说明如何利用 MPI 进行进程间的通信，而 OpenMP 则负责线程内的并行化：

#include <mpi.h>
#include <omp.h>
#define ROWS_A MPI_SIZE_T(100)
#define COLS_A MPI_SIZE_T(100)
#define ROWS_B MPI_SIZE_T(100)
#define COLS_B MPI_SIZE_T(100)

double** mat_a;
double** mat_b;
double** result;

void init_matrices() {
  // 初始化矩阵
}

void matrix_multiply(int my_rank, int size) {
  double* row_a = new double[COLS_B];
  double* row_b = new double[ROWS_A];

  #pragma omp parallel for private(row_a, row_b)
  for (int i = 0; i < ROWS_A; ++i) {
    #pragma omp task firstprivate(row_b)
    {
      for (int j = 0; j < COLS_B; ++j) {
        row_b[j] = mat_b[i][j]; // 从B矩阵取值
      }
      #pragma omp critical
      {
        if (my_rank == 0) { // 使用OpenMP的critical区段同步，避免数据竞争
          for (int k = 0; k < COLS_A; ++k) {
            result[i][k] += mat_a[i][k] * row_b[k];
          }
        }
      }
    }
  }

  delete[] row_a;
  delete[] row_b;
}

int main(int argc, char* argv[]) {
  MPI_Init(&argc, &argv);
  MPI_Comm_size(MPI_COMM_WORLD, &size);
  MPI_Comm_rank(MPI_COMM_WORLD, &my_rank);

  // 分配矩阵内存，初始化等操作...

  matrix_multiply(my_rank, size); // 根据rank分配任务

  // MPI finalize 和其他收发过程...

  return 0;
}

在这个例子中，每个进程处理一部分矩阵A的行，并通过OpenMP对B矩阵的列进行并行计算。`#pragma omp task`指令将每个B矩阵列的计算作为一个独立的任务，然后在`#pragma omp critical`区段内更新结果矩阵。