1

使用库“immintrin.h”，我可以为简单的循环和操作编写 SSE 指令。但是，如何为所示语句编写 SSE 指令？

for (int i =0; i<n; i++){
   for (int j=0; j<n; j++) {
     x[i] += a[i] + a[j];
}}

x 和 a 是使用 _mm_malloc() 初始化的 float*。内存访问模式可用作 __m128 和 4 字节的展开策略。

对不起，如果我不太清楚，但就像

for (int i = 0; i < vecsize; i+=4) {
        __m128 a  = _mm_load_ps(a+i);
        __m128 x  = _mm_add_ps(x,a);
        _mm_store_ps(x+i, x);
    }

（仅适用于 1 个循环），我想要上面显示的循环类似的东西。

编辑：我（EricPostpischil）正在从评论中注入此文本，因为它对问题陈述很重要。作者 NeilDA 应该对此进行扩展：

......在我的程序中'a'总是在变化，因此我想要'x'随之变化。

我已经成功了！我提交了答案..

4

4 回答 4

3

这只是部分答案，但对于评论来说太长/太详细了。

我怀疑你的问题是否写得正确。如图所示，对于 each x[i]，它相加a[i] n次，每个相加a[j]一次，对于 0≤<em>j< n。所以它相当于：

sum = 0;
for (j = 0; j < n; ++j)
    sum += a[j];
for (i = 0; i < n; ++i)
    x[i] += n*a[i] + sum;

这将使用比其他可能的数组操作更简单的 SSE 代码来实现。当然，像上面那样简单地重写它会产生比原始公式更快的代码。

于 2013-03-11T22:52:39.470 回答
0 
__m128 x_v, a_v, aj_v;
for (int i = 0; i < (vec_size); i+=4) {
            x_v = _mm_load_ps(x + i);
            a_v = _mm_load_ps(a+i);
        for (int j = 0; j < N; j++) {
            aj_v = _mm_set1_ps(a[j]);
            x_v = _mm_add_ps(x_v, _mm_add_ps(aj_v, a_v));
        }
            _mm_store_ps(x + i, x_v);
    }

我不知道它是否可以进一步改进，但这将我的时间从 0.17 秒减少到 0.04 秒：D 任何改进它的评论或方法都会很棒！

于 2013-03-12T02:02:05.527 回答
0

转换：

    for (i = 0; i < n; i++)
        for (j = 0; j < n; j++)
            x[i] += a[i] + a[j];

进入：

    for (i = 0; i < n; i++)
        x[i] += n*a[i] + sum(a));

请参阅f_sse()下面的代码：

#include <stdio.h>
#include <string.h>
#include <immintrin.h>

enum {
    N = 4,
};

float x[N], a[N] = { .1, .2, .3, .4 }, y[N];

void f(float *x, float *a, int n)
{
    int i, j;
    for (i = 0; i < n; i++)
        for (j = 0; j < n; j++)
            x[i] += a[i] + a[j];
}

float array_sum(float *a, int n)
{
    /* could be vectorized as well */
    int i;
    float s;
    for (s = 0, i = 0; i < n; i++)
        s += a[i];
    return s;
}

void f_sse(float *x, float *a, int n)
{
    int i, l;
    float t;
    __m128 sum_a, n_vec;
    t = array_sum(a, n);
    sum_a = _mm_set1_ps(t);
    n_vec = _mm_set1_ps(n);
    l = n / 4;
    for (i = 0; i < l; i += 4) {
        __m128 ai, xi;
        ai = _mm_load_ps(a + i);
        xi = _mm_load_ps(x + i);
        ai = _mm_mul_ps(ai, n_vec);
        ai = _mm_add_ps(ai, sum_a);
        xi = _mm_add_ps(xi, ai);
        _mm_store_ps(x + i, xi);
    }
}

int main()
{
    int i, r;
    f(x, a, N);
    f_sse(y, a, N);
    r = memcmp(x, y, N);
    if (r == 0)
        return 0;
    printf("x: { ");
    for (i = 0; i < N; i++)
        printf("%f, ", x[i]);
    printf("}\n");
    printf("y: { ");
    for (i = 0; i < N; i++)
        printf("%f, ", y[i]);
    printf("}\n");
    return 3;
}

由于您声明a正在同时更新并且您不能进行上面的循环转换。您需要有意识地决定何时从a. 它永远不会匹配原始的非矢量化版本，因为我们一次加载 4 个浮点数。

a[j]在内循环中重新加载：

void f_sse(float *x, float *a, int n)
{
    int i, j, l;
    l = n / 4;
    for (i = 0; i < l; i += 4) {
        __m128 ai, xi;
        ai = _mm_load_ps(a + i);
        xi = _mm_load_ps(x + i);
        for (j = 0; j < n; j++) {
            /* re-load a[j] to get updates */
            xi = _mm_add_ps(xi, _mm_add_ps(ai, _mm_set1_ps(((volatile float *)a)[j])));
        }
        _mm_store_ps(x + i, xi);
    }

}

重新加载a[j]并(a + i)在内循环中：

void f_sse(float *x, float *a, int n)
{
    int i, j, l;
    l = n / 4;
    for (i = 0; i < l; i += 4) {
        __m128 ai, xi;
        xi = _mm_load_ps(x + i);
        for (j = 0; j < n; j++) {
            /* re-load (a + i) to get updates */
            ai = _mm_load_ps(a + i);
            /* re-load a[j] to get updates */
            xi = _mm_add_ps(xi, _mm_add_ps(ai, _mm_set1_ps(((volatile float *)a)[j])));
        }
        _mm_store_ps(x + i, xi);
    }
}
于 2013-03-11T23:44:29.063 回答
0 
__m128 *mx = (__m128*)x;
__m128 *ma = (__m128*)a;
__m128 temp_a;

for (int i = 0; i < (n>>2); ++i) {
    for (int j = 0; j < (n>>2); ++j) {
        temp_a = _mm_add_ps(*(ma+i), *(ma+j));
        *mx = _mm_add_ps(*mx, temp_a);
    }
    mx++;
}

当然，您必须确保它n是 4 的倍数，并确保 x 初始化为 0，以便累积正确。

于 2013-03-11T22:58:19.047 回答