是否可以说服 clang 在不使用内在函数的情况下自动矢量化此代码?
Is it possible to convince clang to auto-vectorize this code without using intrinsics?
想象一下,我有这个天真的函数来检测球体重叠。这个问题的重点并不是真正讨论在球体上进行命中测试的最佳方法,所以这只是为了说明。
inline bool sphere_hit(float x1, float y1, float z1, float r1,
float x2, float y2, float z2, float r2) {
float xd = (x1 - x2);
float yd = (y1 - y2);
float zd = (z1 - z2);
float max_dist = (r1 + r2);
return xd * xd + yd * yd + zd * zd < max_dist * max_dist;
}
并且我在嵌套循环中调用它,如下:
std::vector<float> xs, ys, zs, rs;
int n_spheres;
// <snip>
int n_hits = 0;
for (int i = 0; i < n_spheres; ++i) {
for (int j = i + 1; j < n_spheres; ++j) {
if (sphere_hit(xs[i], ys[i], zs[i], rs[i],
xs[j], ys[j], zs[j], rs[j])) {
++n_hits;
}
}
}
std::printf("total hits: %d\n", n_hits);
现在,clang(-O3 -march=native
)足够聪明,可以弄清楚如何将这个循环矢量化(并展开)为 256 位 avx2 指令。太棒了!
但是,如果我做任何比增加命中次数更复杂的事情,例如调用某个任意函数 handle_hit(i, j)
,clang 会发出一个朴素的标量版本。
命中应该非常罕见,所以我认为应该检查每个矢量化循环迭代是否 any 通道的值为真,然后跳转到某个标量慢路径如果是这样。这应该可以通过 vcmpltps
后跟 vmovmskps
来实现。但是,即使我用 __builtin_expect(..., 0)
.
包围对 sphere_hit
的调用,我也无法让 clang 发出此代码
确实可以说服 clang 向量化此代码。使用编译器选项
-Rpass-analysis=loop-vectorize -Rpass=loop-vectorize -Rpass-missed=loop-vectorize
,clang声称浮点运算是矢量化的,Godbolt output证实了这一点。 (红色下划线 for
不是错误,而是矢量化报告)。
通过将 sphere_hit
的结果作为字符存储到临时数组 hitx8
中,可以实现向量化。
之后,通过从内存中读取 8 个字符作为一个 uint64_t a
,每次迭代测试 8 sphere_hit
个结果。这应该是非常有效的,因为条件 a!=0
(见下面的代码)仍然很少见,因为球体命中非常罕见。此外,数组 hitx8
大部分时间可能在 L1 或 L2 缓存中。
我没有测试代码的正确性,但至少 auto-vectorization 这个想法应该可行。
/* clang -Ofast -Wall -march=broadwell -Rpass-analysis=loop-vectorize -Rpass=loop-vectorize -Rpass-missed=loop-vectorize */
#include<string.h>
char sphere_hit(float x1, float y1, float z1, float r1,
float x2, float y2, float z2, float r2);
void handle_hit(int i, int j);
void vectorized_code(float* __restrict xs, float* __restrict ys, float* __restrict zs, float* __restrict rs, char* __restrict hitx8, int n_spheres){
unsigned long long int a;
for (int i = 0; i < n_spheres; ++i) {
for (int j = i + 1; j < n_spheres; ++j){
/* Store the boolean results temporarily in char array hitx8. */
/* The indices of hitx8 are shifted by i+1, so the loop */
/* starts with hitx8[0] */
/* char array hitx8 should have n_spheres + 8 elements */
hitx8[j-i-1] = sphere_hit(xs[i], ys[i], zs[i], rs[i],
xs[j], ys[j], zs[j], rs[j]);
}
for (int j = n_spheres; j < n_spheres+8; ++j){
/* Add 8 extra zeros at the end of hitx8. */
hitx8[j-i-1] = 0; /* hitx8 is 8 elements longer than xs */
}
for (int j = i + 1; j < n_spheres; j=j+8){
memcpy(&a,&hitx8[j-i-1],8);
/* Check 8 sphere hits in parallel: */
/* one `unsigned long long int a` contains 8 boolean values here */
/* The condition a!=0 is still rare since sphere hits are very rare. */
if (a!=0ull){
if (hitx8[j-i-1+0] != 0) handle_hit(i,j+0);
if (hitx8[j-i-1+1] != 0) handle_hit(i,j+1);
if (hitx8[j-i-1+2] != 0) handle_hit(i,j+2);
if (hitx8[j-i-1+3] != 0) handle_hit(i,j+3);
if (hitx8[j-i-1+4] != 0) handle_hit(i,j+4);
if (hitx8[j-i-1+5] != 0) handle_hit(i,j+5);
if (hitx8[j-i-1+6] != 0) handle_hit(i,j+6);
if (hitx8[j-i-1+7] != 0) handle_hit(i,j+7);
}
}
}
}
inline char sphere_hit(float x1, float y1, float z1, float r1,
float x2, float y2, float z2, float r2) {
float xd = (x1 - x2);
float yd = (y1 - y2);
float zd = (z1 - z2);
float max_dist = (r1 + r2);
return xd * xd + yd * yd + zd * zd < max_dist * max_dist;
}
想象一下,我有这个天真的函数来检测球体重叠。这个问题的重点并不是真正讨论在球体上进行命中测试的最佳方法,所以这只是为了说明。
inline bool sphere_hit(float x1, float y1, float z1, float r1,
float x2, float y2, float z2, float r2) {
float xd = (x1 - x2);
float yd = (y1 - y2);
float zd = (z1 - z2);
float max_dist = (r1 + r2);
return xd * xd + yd * yd + zd * zd < max_dist * max_dist;
}
并且我在嵌套循环中调用它,如下:
std::vector<float> xs, ys, zs, rs;
int n_spheres;
// <snip>
int n_hits = 0;
for (int i = 0; i < n_spheres; ++i) {
for (int j = i + 1; j < n_spheres; ++j) {
if (sphere_hit(xs[i], ys[i], zs[i], rs[i],
xs[j], ys[j], zs[j], rs[j])) {
++n_hits;
}
}
}
std::printf("total hits: %d\n", n_hits);
现在,clang(-O3 -march=native
)足够聪明,可以弄清楚如何将这个循环矢量化(并展开)为 256 位 avx2 指令。太棒了!
但是,如果我做任何比增加命中次数更复杂的事情,例如调用某个任意函数 handle_hit(i, j)
,clang 会发出一个朴素的标量版本。
命中应该非常罕见,所以我认为应该检查每个矢量化循环迭代是否 any 通道的值为真,然后跳转到某个标量慢路径如果是这样。这应该可以通过 vcmpltps
后跟 vmovmskps
来实现。但是,即使我用 __builtin_expect(..., 0)
.
sphere_hit
的调用,我也无法让 clang 发出此代码
确实可以说服 clang 向量化此代码。使用编译器选项
-Rpass-analysis=loop-vectorize -Rpass=loop-vectorize -Rpass-missed=loop-vectorize
,clang声称浮点运算是矢量化的,Godbolt output证实了这一点。 (红色下划线 for
不是错误,而是矢量化报告)。
通过将 sphere_hit
的结果作为字符存储到临时数组 hitx8
中,可以实现向量化。
之后,通过从内存中读取 8 个字符作为一个 uint64_t a
,每次迭代测试 8 sphere_hit
个结果。这应该是非常有效的,因为条件 a!=0
(见下面的代码)仍然很少见,因为球体命中非常罕见。此外,数组 hitx8
大部分时间可能在 L1 或 L2 缓存中。
我没有测试代码的正确性,但至少 auto-vectorization 这个想法应该可行。
/* clang -Ofast -Wall -march=broadwell -Rpass-analysis=loop-vectorize -Rpass=loop-vectorize -Rpass-missed=loop-vectorize */
#include<string.h>
char sphere_hit(float x1, float y1, float z1, float r1,
float x2, float y2, float z2, float r2);
void handle_hit(int i, int j);
void vectorized_code(float* __restrict xs, float* __restrict ys, float* __restrict zs, float* __restrict rs, char* __restrict hitx8, int n_spheres){
unsigned long long int a;
for (int i = 0; i < n_spheres; ++i) {
for (int j = i + 1; j < n_spheres; ++j){
/* Store the boolean results temporarily in char array hitx8. */
/* The indices of hitx8 are shifted by i+1, so the loop */
/* starts with hitx8[0] */
/* char array hitx8 should have n_spheres + 8 elements */
hitx8[j-i-1] = sphere_hit(xs[i], ys[i], zs[i], rs[i],
xs[j], ys[j], zs[j], rs[j]);
}
for (int j = n_spheres; j < n_spheres+8; ++j){
/* Add 8 extra zeros at the end of hitx8. */
hitx8[j-i-1] = 0; /* hitx8 is 8 elements longer than xs */
}
for (int j = i + 1; j < n_spheres; j=j+8){
memcpy(&a,&hitx8[j-i-1],8);
/* Check 8 sphere hits in parallel: */
/* one `unsigned long long int a` contains 8 boolean values here */
/* The condition a!=0 is still rare since sphere hits are very rare. */
if (a!=0ull){
if (hitx8[j-i-1+0] != 0) handle_hit(i,j+0);
if (hitx8[j-i-1+1] != 0) handle_hit(i,j+1);
if (hitx8[j-i-1+2] != 0) handle_hit(i,j+2);
if (hitx8[j-i-1+3] != 0) handle_hit(i,j+3);
if (hitx8[j-i-1+4] != 0) handle_hit(i,j+4);
if (hitx8[j-i-1+5] != 0) handle_hit(i,j+5);
if (hitx8[j-i-1+6] != 0) handle_hit(i,j+6);
if (hitx8[j-i-1+7] != 0) handle_hit(i,j+7);
}
}
}
}
inline char sphere_hit(float x1, float y1, float z1, float r1,
float x2, float y2, float z2, float r2) {
float xd = (x1 - x2);
float yd = (y1 - y2);
float zd = (z1 - z2);
float max_dist = (r1 + r2);
return xd * xd + yd * yd + zd * zd < max_dist * max_dist;
}