使用 OpenMP 优化矩阵转置函数答案

【问题标题】：Optimizing a matrix transpose function with OpenMP使用 OpenMP 优化矩阵转置函数
【发布时间】：2022-01-25 01:38:21
【问题描述】：

我有这段代码使用循环平铺策略转置矩阵。

void transposer(int n, int m, double *dst, const double *src) {
   int blocksize;
   for (int i = 0; i < n; i += blocksize) {
      for (int j = 0; j < m; j += blocksize) {
          // transpose the block beginning at [i,j]
          for (int k = i; k < i + blocksize; ++k) {
              for (int l = j; l < j + blocksize; ++l) {
                  dst[k + l*n] = src[l + k*m];
               }
          }
      }
   }
}

我想通过使用 OpenMP 的多线程来优化这一点，但是当有这么多嵌套的 for 循环时，我不确定该怎么做。我想过只添加#pragma omp parallel for，但这不只是并行化外循环吗？

【问题讨论】：

这能回答你的问题吗？ How does OpenMP handle nested loops?
@JanezKuhar 有点，但是我更想知道是否可以进行任何进一步的优化。在您链接的线程中，他们只真正谈论崩溃等

标签： c multithreading performance parallel-processing openmp

【解决方案1】：

当您尝试并行化循环嵌套时，您应该问自己有多少级别是无冲突的。如：每次迭代都写入不同的位置。如果两次迭代（可能）写入同一个位置，您需要 1. 使用归约 2. 使用临界区或其他同步 3. 确定此循环不值得并行化，或 4. 重写您的算法。

在您的情况下，写入位置取决于k,l。由于k<n 和l*n，没有写入同一位置的k.l / k',l' 对。此外，没有两个内部迭代具有相同的k 或l 值。所以四个循环都是并行的，而且它们是完美嵌套的，所以你可以使用collapse(4)。

您也可以通过抽象地考虑算法得出这个结论：在矩阵转置中，每个目标位置都只写入一次，因此无论您如何遍历目标数据结构，它都是完全并行的。

【讨论】：

@MarcusF“还有别的吗”：一定要设置OMP_PROC_BIND=true。这在双插槽设计中当然很重要，但我会默认设置它。如果您不这样做，您的操作系统可能会迁移线程并且您会丢失缓存位置。
嗯，我正在尝试完全理解collapse，并且我读过循环不能相互依赖。但是最里面的两个 for 循环有 k = i 和 l = j，这不意味着它们是依赖的吗？你能澄清一个数据依赖的例子吗？

【解决方案2】：

您可以使用折叠说明符在两个循环上并行化。

#   pragma omp parallel for collapse(2) 
    for (int i = 0; i < n; i += blocksize) {
        for (int j = 0; j < m; j += blocksize) {
            // transpose the block beginning at [i,j]
            for (int k = i; k < i + blocksize; ++k) {
                for (int l = j; l < j + blocksize; ++l) {
                    dst[k + l*n] = src[l + k*m];
                }
            }
        }
    }

作为旁注，我认为您应该交换最里面的两个循环。通常，当您在顺序写入和顺序读取之间进行选择时，写入对性能更重要。

【讨论】：

【解决方案3】：

我想过只是添加#pragma omp parallel for 但不是这样只是并行化外循环？

是的。要并行化多个连续循环，可以使用 OpenMP 的 collapse 子句。但是请记住：

（正如Victor Eijkhout 指出的那样）。即使这并不直接适用于您的代码 sn-p，通常情况下，对于每个要并行化的新循环，都应该考虑潜在的较新的 race-conditions 例如，并行化可能已添加。例如，不同的线程同时写入同一个 dst 位置。
在某些情况下，并行嵌套循环可能会导致比并行单个循环更慢的执行时间。因为，collapse 子句的具体实现使用更复杂的启发式（比简单的循环并行化）在线程之间划分循环的迭代，这可能导致开销高于它的收益提供。

您应该尝试使用单个并行循环进行基准测试，然后使用两个，依此类推，并相应地比较结果。

void transposer(int n, int m, double *dst, const double *src) {
    int blocksize;
    #pragma omp parallel for collapse(...)
    for (int i = 0; i < n; i += blocksize)
       for (int j = 0; j < m; j += blocksize)
           for (int k = i; k < i + blocksize; ++k
               for (int l = j; l < j + blocksize; ++l)
                  dst[k + l*n] = src[l + k*m];
}

根据线程数、内核数、矩阵大小以及其他因素，顺序运行实际上可能比并行版本更快。在您的代码中尤其如此

【讨论】：