可以矢量化myNum + = a [b [i]]*c [i]; 在x86_64上？

Question

我将使用什么内在函数来对x86_64上的以下内容进行矢量化(如果它甚至可以进行矢量化)？

double myNum = 0;
for(int i=0;i<n;i++){
    myNum += a[b[i]] * c[i]; //b[i] = int, a[b[i]] = double, c[i] = double
}

Answer 1

这是我的全力以赴,经过全面优化和测试:

#include <emmintrin.h>

__m128d sum = _mm_setzero_pd();
for(int i=0; i<n; i+=2) {
    sum = _mm_add_pd(sum, _mm_mul_pd(
        _mm_loadu_pd(c + i),
        _mm_setr_pd(a[b[i]], a[b[i+1]])
    ));
}

if(n & 1) {
    sum = _mm_add_pd(sum, _mm_set_sd(a[b[n-1]] * c[n-1]));
}

double finalSum = _mm_cvtsd_f64(_mm_add_pd(
    sum, _mm_shuffle_pd(sum, sum, _MM_SHUFFLE2(0, 1))
));

这使用gcc -O2 -msse2(4.4.1)生成非常漂亮的汇编代码.

正如你所知道的那样,拥有一个偶数n会使这个循环变得更快以及一个对齐c.如果你可以对齐c,更改_mm_loadu_pd到_mm_load_pd一个甚至更快的执行时间.

Answer 2

我将从展开循环开始。就像是

double myNum1 = 0, myNum2=0;
for(int i=0;i<n;i+=2)
{
    myNum1 += a[b[ i ]] * c[ i ];
    myNum2 += a[b[i+1]] * c[i+1];
}
// ...extra code to handle the remainder when n isn't a multiple of 2...
double myNum = myNum1 + myNum2;

希望这允许编译器将负载与算术交错；配置文件并查看装配体，看看是否有改进。理想情况下，编译器会生成 SSE 指令，但如果实际发生这种情况，我不会生成 SSE 指令。

进一步展开可能会让你这样做：

__m128d sum0, sum1;
// ...initialize to zero...
for(int i=0;i<n;i+=4)
{
    double temp0 = a[b[ i ]] * c[ i ];
    double temp1 = a[b[i+1]] * c[i+1];
    double temp2 = a[b[i+2]] * c[i+2];
    double temp3 = a[b[i+3]] * c[i+3];
    __m128d pair0 = _mm_set_pd(temp0, temp1);
    __m128d pair1 = _mm_set_pd(temp2, temp3);
    sum0 = _mm_add_pd(sum0, pair0);
    sum1 = _mm_add_pd(sum1, pair1);
}
// ...extra code to handle the remainder when n isn't a multiple of 4...
// ...add sum0 and sum1, then add the result's components...

（对开头和结尾的伪代码表示歉意，我认为重要的部分是循环）。我不确定这是否会更快；这取决于各种延迟以及编译器重新安排一切的能力。确保在之前和之后进行分析，看看是否有实际的改进。

希望有帮助。