c - 使用 SSE 内在函数矢量化 2D 模板

Question

我正在尝试仅使用对齐、加载和存储来矢量化 2D 模板。为此，我基本上想使用 _mm_load_ps 和 _mm_shuffle_ps 来获取所需的地址。

我的代码的标量版本是:

    void FDTD_base (float *V, float *U, int dx, int dy, float c0, float c1, float c2, float c3, float c4)
    {
    int i, j, k;

            for (j = 4; j < dy-4; j++)
            {
                    for (i = 4; i < dx-4; i++)
                    {

                            U[j*dx+i] = (c0 * (V[j*dx+i]) //center
                                    + c1 * (V[j*dx+(i-1)] + V[(j-1)*dx+i] + V[j*dx+(i+1)] + V[(j+1)*dx+i] )
                                    + c2 * (V[j*dx+(i-2)] + V[(j-2)*dx+i] + V[j*dx+(i+2)] + V[(j+2)*dx+i] )
                                    + c3 * (V[j*dx+(i-3)] + V[(j-3)*dx+i] + V[j*dx+(i+3)] + V[(j+3)*dx+i] )
                                    + c4 * (V[j*dx+(i-4)] + V[(j-4)*dx+i] + V[j*dx+(i+4)] + V[(j+4)*dx+i] ));

                    }
            }

      }

我的 vector 到目前为止看到的代码版本:

     for (j = 4; j < dy-4; j++)
    {
            for (i = 4; i < dx-4; i+=4)
            {
                    __m128 b = _mm_load_ps(&V[j*dx+i]);
                    center = _mm_mul_ps(b,c0_i);
                    a = _mm_load_ps(&V[j*dx+(i-4)]);
                    c = _mm_load_ps(&V[j*dx+(i+4)]);

                    d = _mm_load_ps(&V[(j-4)*dx+i]);
                    e = _mm_load_ps(&V[(j+4)*dx+i]);

                    u_i2 = _mm_shuffle_ps(a,b,_MM_SHUFFLE(1,0,3,2));//i-2
                    u_i6 = _mm_shuffle_ps(b,c,_MM_SHUFFLE(1,0,3,2));//i+2

                    u_i1 = _mm_shuffle_ps(u_i2,b,_MM_SHUFFLE(2,1,2,1));//i-1
                    u_i5 = _mm_shuffle_ps(b,u_i6,_MM_SHUFFLE(2,1,2,1));//i+1

                    u_i3 = _mm_shuffle_ps(a,u_i2,_MM_SHUFFLE(2,1,2,1));//i-3
                    u_i7 = _mm_shuffle_ps(u_i6,c,_MM_SHUFFLE(2,1,2,1));//i+3

                    u_i4 = a; //i-4
                    u_i8 = c; //i+4

谁能帮助我获得 j-1,j+1 .....j-4,j+4 的位置。

这不起作用:

                    u_j2 = _mm_shuffle_ps(d,b,_MM_SHUFFLE(1,0,3,2));//j-2 (this is incorrect)
                    u_j6 = _mm_shuffle_ps(b,e,_MM_SHUFFLE(1,0,3,2));//j+2

                    u_j1 = _mm_shuffle_ps(u_j2,b,_MM_SHUFFLE(2,1,2,1));//j-1
                    u_j5 = _mm_shuffle_ps(b,u_j6,_MM_SHUFFLE(2,1,2,1));//j+1

                    u_j3 = _mm_shuffle_ps(d,u_j2,_MM_SHUFFLE(2,1,2,1));//j-3
                    u_j7 = _mm_shuffle_ps(u_j6,e,_MM_SHUFFLE(2,1,2,1));//j+3

                    u_j4 = d; //j-4 (this is fine)
                    u_j8 = e; //j+4

我只需要帮助确定如何获取(j-1)*dx+i,(j+1)*dx+1 ..... (j-4)*dx+i 和 (j+4)*dx+i 不使用未对齐的负载。

作为一个潜在的解决方案，我想到了向存储在 d 中的地址添加位移 3*dx 以获得 (j-1)*dx+i。并将e中存储的地址减去3*dx的位移得到(j+1)*dx+i。同理将2*dx添加到d的地址得到j-2等。但我不知道如何使用 SSE 内在函数来实现此策略。

请帮忙。我正在使用英特尔 icc 编译器。

Answer 1

“谁能帮我获得 j-1,j+1 .....j-4,j+4 的位置。” - 这些不需要洗牌；它们已经与您的 SIMD channel 对齐。

u_j2 = _mm_load_ps(&V[(j-2)*dx+i]); 
u_j6 = _mm_load_ps(&V[(j+2)*dx+i]); 
u_j1 = _mm_load_ps(&V[(j-1)*dx+i]); 
u_j5 = _mm_load_ps(&V[(j+1)*dx+i]); 
// and so forth

您绝对不能通过任何可能的重新排列从您标记为 d 和 e 的变量中获取这些，因为 d 中的值(例如)是V[j-4, i], V[j-4, i+1], V[j-4, i+2], V[j-4, i +3]，你无法从中得到 V[j-2, i]。

提示:考虑 SIMD channel ；这清楚地表明您需要水平而不是垂直重新排列。

提示:考虑当内部循环计数器递增 (i+=4) 时会发生什么。上一个循环中的 u_i5 (V[j, i+1..i+5]) 现在是 u_i3 (V[j, i-3..i+1]) 在当前循环中。您正在计算行中数据的每个偏移版本至少两次。您或许可以将循环展开几次并避免执行所有额外的工作。

提示:为什么不使用 AVX？使用 _mm256_permute_ps(如果需要，使用 _mm256_permute2f128_ps)进行随机播放，以及相应的加载指令。它几乎快两倍，因为您拥有两倍宽的 SIMD 寄存器，并且大多数 AVX 指令在现代 CPU 上仍然只需要一个周期，与 SSE 指令相同。