`#pragma parallel for collapse` 和 `#pragma omp parallel for` 之间的区别

Question

首先，这个问题可能会有点误导，我了解并行区域中的折叠子句与没有区域的区域之间的主要区别。假设我要转置一个矩阵，有以下两种方法，第一种是与SIMD 指令并行用于内部循环，第二种方法使用collapse(2) 子句

#pragma omp parallel for
    for(int i=0; i<rows; i++){
#pragma omp simd
      for(int j=0; j<columns; j++){
         *(output + j * rows + i) = *(input + i * columns + j);
    }
}

#pragma omp parallel for collapse(2)
    for(int i=0; i<rows; i++){
      for(int j=0; j<columns; j++){
         *(output + j * rows + i) = *(input + i * columns + j);
    }

在上述两种方法中，尤其是在缓存方面，哪种方法更有效？
在上述两个中，哪个实现会更高效、更快？有什么方法可以通过查看实现来确定。

鉴于所有循环计数器都是相互独立的，是否可以设置一个关于何时使用的基本准则？
TIA

Answer 1

TL;DR:这两种实现都非常低效。第二个在实践中可能会比第一个慢，尽管理论上它可以更好地扩展。

第一个实现不太可能被向量化，因为访问在内存中不连续。 GCC 10 和 Clang 11 都会生成低效的代码。 关键是 OpenMP 没有提供高级 SIMD 构造 来处理数据转置！因此，如果您想高效地完成此操作，您可能需要亲自动手（或使用为您完成此操作的外部库）。

第二个实现可能比第一个实现慢得多，因为循环迭代器是线性化的，通常会导致在热路径中执行更多指令。一些实现（例如 Clang 11 和 ICC 19 但不是 GCC 10）甚至使用非常慢的模数运算（例如 div 指令）来执行此操作，从而导致更慢的循环。

第二个实现在理论上也应该比第一个实现更好的扩展，因为collapse 子句提供了更多的并行性。实际上，在第一个实现中，n 线程之间只有 rows 行可以共享。因此，如果您在大型并行机器或宽矩形矩阵上工作，n 与rows 相比并不那么小，这可能会导致一些工作不平衡，甚至线程饥饿。

---

为什么两种实现都效率低下

由于内存访问模式，这两种实现效率低下。实际上，在大矩阵上，output 中的写入不是连续的，会导致许多缓存未命中。将写入一个完整的高速缓存行（在大多数常见架构上为 64 字节），而只写入几个字节。如果 columns 是 2 的幂，则会发生缓存抖动并进一步降低性能。

缓解这些问题的一种解决方案是使用平铺。这是一个例子：

// Assume rows and columns are nice for sake of clarity ;)
constexpr int tileSize = 8;
assert(rows % tileSize == 0);
assert(columns % tileSize == 0);

// Note the collapse clause is needed here for scalability and 
// the collapse overhead is mitigated by the inner loop.
#pragma omp parallel for collapse(2)
for(int i=0; i<rows; i+=tileSize)
{
    for(int j=0; j<columns; j+=tileSize)
    {
        for(int ti=i; ti<i+tileSize; ++ti)
        {
            for(int tj=j; tj<j+tileSize; ++tj)
            {
                output[tj * rows + ti] = input[ti * columns + tj];
            }
        }
    }
}

上面的代码应该更快，但不是最优的。成功编写快速转置代码具有挑战性。以下是一些改进代码的建议：

• 使用临时切片缓冲区来改进内存访问模式（因此编译器可以使用快速 SIMD 指令）
• 使用方形切片来提高缓存的使用率
• 使用多级平铺来改进 L2/L3 缓存的使用或使用 Z 平铺方法

或者，您可以简单地使用快速的 BLAS 实现，提供非常优化的矩阵转置函数（并非全部都可以，但 AFAIK OpenBLAS 和 MKL 可以）。

PS：我假设矩阵以行优先顺序存储。