矢量化模运算

Question

我正在尝试编写一些合理快速的组件向量加法代码.我正在使用(签名,我相信)64位整数.

功能是

void addRq (int64_t* a, const int64_t* b, const int32_t dim, const int64_t q) {
    for(int i = 0; i < dim; i++) {
        a[i] = (a[i]+b[i])%q; // LINE1
    }
}

我在icc -std=gnu99 -O3IvyBridge(SSE4.2和AVX,但不是AVX2)上编译(icc以便我以后可以使用SVML).

我的基线是%q从LINE1中删除.100(迭代)函数调用dim=11221184需要1.6秒.ICC自动矢量化SSE代码; 大.

我真的想做模块化的补充.使用%q,ICC不会自动向量化代码,它在11.8秒(!)内运行.即使忽略了之前尝试的自动矢量化,这似乎仍然过分.

由于我没有AVX2,SSE的矢量化需要SVML,这也许就是为什么ICC没有自动矢量化的原因.无论如何,这是我尝试对内循环进行矢量化:

__m128i qs = _mm_set1_epi64x(q);
for(int i = 0; i < dim; i+=2) {
    __m128i xs = _mm_load_si128((const __m128i*)(a+i));
    __m128i ys = _mm_load_si128((const __m128i*)(b+i));
    __m128i zs = _mm_add_epi64(xs,ys);
    zs = _mm_rem_epi64(zs,qs);
    _mm_store_si128((__m128i*)(a+i),zs);
}

主循环的汇编是:

..B3.4:                         # Preds ..B3.2 ..B3.12
    movdqa    (%r12,%r15,8), %xmm0                          #59.22
    movdqa    %xmm8, %xmm1                                  #60.14
    paddq     (%r14,%r15,8), %xmm0                          #59.22
    call      __svml_i64rem2                                #61.9
    movdqa    %xmm0, (%r12,%r15,8)                          #61.36
    addq      $2, %r15                                      #56.30
    cmpq      %r13, %r15                                    #56.24
    jl        ..B3.4        # Prob 82%                      #56.24

因此代码按预期进行矢量化.我知道由于SVML,我可能不会获得2倍的加速,但是代码在12.5秒内运行,比没有矢量化的速度慢!这真的是最好的吗？

Answer 1

SSE2和AVX2都没有整数除法指令.英特尔对于调用SVML函数内在函数是不诚实的,因为它们中的许多都是复杂的函数,它们映射到几个指令而不仅仅是几个指令.

有一种方法可以使用SSE2或AVX2进行更快的分割(和模数).参见本文通过不变整数改进除法.基本上你预先计算一个除数然后进行乘法运算.预计算除数需要时间,但是对于dim代码中的某些值,它应该胜出.我在这里更详细地描述了这个方法的SSE整数除法？ 我还在素数查找器中成功实现了此方法.使用SIMD-SSE/AVX查找素数列表

Agner Fog 使用该论文中描述的方法在Vector类中实现32位(但不是64位)除法.如果您想要一些代码,那么这将是一个很好的起点,但您必须将其扩展到64位.

编辑:根据Mysticial的评论并假设输入已经减少,我为SSE制作了一个版本. 如果这是在MSVC中编译的,那么它需要处于64位模式,因为32位模式不支持_mm_set1_epi64x.这可以修复为32位模式模式,但我不想这样做.

#ifdef _MSC_VER 
#include <intrin.h>
#endif
#include <nmmintrin.h>                 // SSE4.2
#include <stdint.h>
#include <stdio.h>

void addRq_SSE(int64_t* a, const int64_t* b, const int32_t dim, const int64_t q) {
    __m128i q2 = _mm_set1_epi64x(q);
    __m128i t2 = _mm_sub_epi64(q2,_mm_set1_epi64x(1));
    for(int i = 0; i < dim; i+=2) {
        __m128i a2 = _mm_loadu_si128((__m128i*)&a[i]);
        __m128i b2 = _mm_loadu_si128((__m128i*)&b[i]);
        __m128i c2 = _mm_add_epi64(a2,b2);
        __m128i cmp = _mm_cmpgt_epi64(c2, t2);
        c2 = _mm_sub_epi64(c2, _mm_and_si128(q2,cmp));
        _mm_storeu_si128((__m128i*)&a[i], c2);
    }
}

int main() {
    const int64_t dim = 20;
    int64_t a[dim];
    int64_t b[dim];
    int64_t q = 10;

    for(int i=0; i<dim; i++) {
        a[i] = i%q; b[i] = i%q;
    }
    addRq_SSE(a, b, dim, q);
    for(int i=0; i<dim; i++) {
        printf("%d\n", a[i]);
    }   
}