使用数组结构数组 (AoSoA) 并一次处理八个点。在下面的代码中point8
是包含八个点的数组结构。功能rotate_point8
旋转八个点并与函数具有相同的代数结构rotate_point
旋转单个点。功能rotate_all8
使用 AoSoA 旋转 32 点point8*
.
单点旋转代码执行 4 次乘法、1 次加法和 1 次减法。
如果我们看一下大会为rotate_point8 https://godbolt.org/g/EYXNN5我们看到 GCC 展开循环并在每次展开时执行 4 次 SIMD 乘法、1 次 SIMD 加法和 1 次 SIMD 减法。这是你能做的最好的事情:用一件的价格买八件。
#include <x86intrin.h>
#include <stdio.h>
#include <math.h>
struct point8 {
__m256 x;
__m256 y;
};
struct point {
float x;
float y;
};
static point rotate_point(point p, float a, float b) {
point r;
r.x = p.x*a - p.y*b;
r.y = p.x*b + p.y*a;
return r;
}
static point8 rotate_point8(point8 p, float a, float b) {
__m256 va = _mm256_set1_ps(a), vb = _mm256_set1_ps(b);
point8 r;
r.x = _mm256_sub_ps(_mm256_mul_ps(p.x,va), _mm256_mul_ps(p.y,vb));
r.y = _mm256_add_ps(_mm256_mul_ps(p.x,vb), _mm256_mul_ps(p.y,va));
return r;
}
void rotate_all(point* points, point* r, float angle) {
float a = cos(angle), b = sin(angle);
for(int i=0; i<32; i++) r[i] = rotate_point(points[i], a, b);
}
void rotate_all8(point8* points, point8* r8, float angle) {
float a = cos(angle), b = sin(angle);
for(int i=0; i<4; i++) r8[i] = rotate_point8(points[i], a, b);
}
int main(void) {
float x[32], y[32];
point p[32], r[32];
point8 p8[4], r8[4];
float angle = 3.14159f/4;
for(int i=0; i<32; i++) y[i] = 1.0*i/31, x[i] = sqrt(1-y[i]*y[i]);
for(int i=0; i<32; i++) p[i].x = x[i], p[i].y = y[i];
for(int i=0; i<4; i++) p8[i].x = _mm256_load_ps(&x[8*i]), p8[i].y = _mm256_load_ps(&y[8*i]);
for(int i=0; i<32; i++) printf("%f %f\n", p[i].x, p[i].y); puts("");
rotate_all(p, r, angle);
for(int i=0; i<32; i++) printf("%f %f\n", r[i].x, r[i].y); puts("");
rotate_all8(p8, r8, angle);
for(int i=0; i<4; i++) {
_mm256_storeu_ps(x, r8[i].x), _mm256_storeu_ps(y, r8[i].y);
for(int j=0; j<8; j++) printf("%f %f\n", x[j], y[j]);
}
}