c - AVX2稀疏矩阵乘法

我正在尝试利用新的AVX2 GATHER指令来加速稀疏矩阵- vector 乘法。矩阵为CSR(或Yale)格式，带有指向行索引数组的行指针，该数组又保存着列。
这样的mat-vec mul的C代码看起来像这样:

for (int row = 0; row < n_rows - 1; row++) {
    double rowsum = 0;
    for (int col = row_ptr[row]; col < row_ptr[row + 1]; col++) {
        rowsum += values[col] * x[col_indices[col]];
    }
    result[row] = rowsum;
}

现在，我的目标是使用AVX2内部函数来加速此过程。以下代码基于https://blog.fox-toolkit.org/?p=174与最新的Intel或GCC配合使用。我删除了其余部分，因为无论如何我的行都对齐4个 double (列％4 == 0)(幸运的是，我)。如果有人感兴趣，我也有处理剩余部分的代码，但是重点是，该代码实际上稍微慢一些。我检查了反汇编，对于以上版本，仅生成FP指令，并且对于我的AVX2代码，所有AVX2操作均按预期方式显示。即使使用适合高速缓存的小型矩阵，AVX2版本也不是一件好事。我在这里感到困惑...

double* value_base = &values[0];
double* x_base = &x[0];
int*    index_base = &col_indices[0];


for (int row = 0; row < n_rows - 1; row++) {
    int col_length   = row_ptr[row + 1] - row_ptr[row];

    __m256d rowsum = _mm256_set1_pd(0.);
    for (int col4 = 0; col4 < col_length; col4 += 4) {
        // Load indices for x vector(const __m128i*)
        __m128i idxreg     = _mm_load_si128((const __m128i*)index_base);
        // Load 4 doubles from x indexed by idxreg (AVX2)
        __m256d x_     = _mm256_i32gather_pd(x_base, idxreg, 8);
        // Load 4 doubles linear from memory (value array)
        __m256d v_     = _mm256_load_pd(value_base);
        // FMA: rowsum += x_ * v_
        rowsum = _mm256_fmadd_pd(x_, v_, rowsum);

        index_base += 4;
        value_base += 4;
    }
    __m256d s = _mm256_hadd_pd(rowsum, rowsum);
    result[row] = ((double*)&s)[0] + ((double*)&s)[2];
    // Alternative (not faster):
    // Now we split the upper and lower AVX register, and do a number of horizontal adds
    //__m256d hsum = _mm256_add_pd(rowsum, _mm256_permute2f128_pd(rowsum, rowsum, 0x1));
    //_mm_store_sd(&result[row], _mm_hadd_pd( _mm256_castpd256_pd128(hsum), _mm256_castpd256_pd128(hsum) ) );
}

欢迎任何建议。

非常感谢，
克里斯

最佳答案

聚集在Haswell上的速度很慢。我以几种不同的方式实现了对16位值的8位索引LUT查找(对于GF16乘以par2)，以找出最快的方法。在Haswell上，VPGATHERDD版本的时间是movd / pinsrw版本的1.7倍。 (只需要几个VPUNPCK/shift指令即可。)code here, if anyone wants to run the benchmark。

第一次引入指令时，通常不会使用大量的硅来使其变得超快。在那里只是为了获得硬件支持，因此可以编写代码来使用它。为了在所有CPU上获得理想的性能，则需要执行x264对pshufb的操作:为Core2之类的CPU提供SLOW_SHUFFLE标志，并将其作为最佳例程查找功能指针设置的因素，而不是仅对CPU起作用。支持。

对于不太热衷于为可以运行的每个CPU调整asm版本的项目，引入无加速版本的指令将使人们更快地使用它，因此当下一个设计出现时，其更快，更快的代码速度得到了提高。发布像Haswell这样的设计，其中聚集实际上是一个缓慢的过程，这有点麻烦。也许他们想看看人们会如何使用它？它确实增加了代码密度，这在聚集不紧密的情况下有帮助。

Broadwell应该有一个更快的聚集实现，但是我无法使用它。列出延迟/吞吐量指令的英特尔手册说，Broadwell的收集速度大约快1.6倍，因此，它仍然比手工制作的循环慢一些，该循环将GP reg中的索引移位/解包，并将它们用于PINSRW转换为 vector 。

如果gather可以利用多个元素具有相同索引，甚至指向同一32B fetch块的索引的情况，则根据输入数据的不同，速度可能会大大提高。

希望Skylake能够进一步改善。我以为我读过一些东西，说可以，但是经检查，我什么都没找到。

RE:稀疏矩阵:没有可以复制数据的格式，因此您可以对行或列进行连续读取吗？这不是我必须编写的代码，但是我认为我已经在一些答案中提到它了。

关于c - AVX2稀疏矩阵乘法，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/31430198/