展开循环并使用矢量化进行独立求和

对于以下循环,如果我告诉它使用关联数学,例如使用-Ofast GCC将仅对循环进行矢量化。

 float sumf(float *x) { x = (float*)__builtin_assume_aligned(x, 64); float sum = 0; for(int i=0; i<2048; i++) sum += x[i]; return sum; } 

这是带-Ofast -mavx的程序集

 sumf(float*): vxorps %xmm0, %xmm0, %xmm0 leaq 8192(%rdi), %rax .L2: vaddps (%rdi), %ymm0, %ymm0 addq $32, %rdi cmpq %rdi, %rax jne .L2 vhaddps %ymm0, %ymm0, %ymm0 vhaddps %ymm0, %ymm0, %ymm1 vperm2f128 $1, %ymm1, %ymm1, %ymm0 vaddps %ymm1, %ymm0, %ymm0 vzeroupper ret 

这清楚地表明循环已被矢量化。

但是这个循环也有一个依赖链。 为了克服添加的延迟,我需要在x86_64上展开并执行至少三次部分和(不包括需要展开八次的Skylake并使用需要在Haswell和Broadwell上展开10次的FMA指令进行添加) 。 据我所知,我可以使用-funroll-loops展开-funroll-loops

这是带有-Ofast -mavx -funroll-loops的程序集。

 sumf(float*): vxorps %xmm7, %xmm7, %xmm7 leaq 8192(%rdi), %rax .L2: vaddps (%rdi), %ymm7, %ymm0 addq $256, %rdi vaddps -224(%rdi), %ymm0, %ymm1 vaddps -192(%rdi), %ymm1, %ymm2 vaddps -160(%rdi), %ymm2, %ymm3 vaddps -128(%rdi), %ymm3, %ymm4 vaddps -96(%rdi), %ymm4, %ymm5 vaddps -64(%rdi), %ymm5, %ymm6 vaddps -32(%rdi), %ymm6, %ymm7 cmpq %rdi, %rax jne .L2 vhaddps %ymm7, %ymm7, %ymm8 vhaddps %ymm8, %ymm8, %ymm9 vperm2f128 $1, %ymm9, %ymm9, %ymm10 vaddps %ymm9, %ymm10, %ymm0 vzeroupper ret 

GCC会将循环展开八次。 但是,它不做独立的总和。 它有八个相关的总和。 这是毫无意义的,没有比没有展开更好的。

如何让GCC展开循环并进行独立的部分和?


编辑:

即使没有针对SSE的-funroll-loops ,Clang也会展开四个独立的部分和,但我不确定它的AVX代码效率如何。 无论如何编译器都不需要-funroll-loops-Ofast ,所以很高兴看到Clang至少对SSE这样做了。

Clang 3.5.1 with -Ofast

 sumf(float*): # @sumf(float*) xorps %xmm0, %xmm0 xorl %eax, %eax xorps %xmm1, %xmm1 .LBB0_1: # %vector.body movups (%rdi,%rax,4), %xmm2 movups 16(%rdi,%rax,4), %xmm3 addps %xmm0, %xmm2 addps %xmm1, %xmm3 movups 32(%rdi,%rax,4), %xmm0 movups 48(%rdi,%rax,4), %xmm1 addps %xmm2, %xmm0 addps %xmm3, %xmm1 addq $16, %rax cmpq $2048, %rax # imm = 0x800 jne .LBB0_1 addps %xmm0, %xmm1 movaps %xmm1, %xmm2 movhlps %xmm2, %xmm2 # xmm2 = xmm2[1,1] addps %xmm1, %xmm2 pshufd $1, %xmm2, %xmm0 # xmm0 = xmm2[1,0,0,0] addps %xmm2, %xmm0 retq 

带有-O3 ICC 13.0.1展开为两个独立的部分和。 ICC显然假设只有-O3关联数学。

 .B1.8: vaddps (%rdi,%rdx,4), %ymm1, %ymm1 #5.29 vaddps 32(%rdi,%rdx,4), %ymm0, %ymm0 #5.29 vaddps 64(%rdi,%rdx,4), %ymm1, %ymm1 #5.29 vaddps 96(%rdi,%rdx,4), %ymm0, %ymm0 #5.29 addq $32, %rdx #5.3 cmpq %rax, %rdx #5.3 jb ..B1.8 # Prob 99% #5.3 

一些gcc内在函数和__builtin_产生了这个:

 typedef float v8sf __attribute__((vector_size(32))); typedef uint32_t v8u32 __attribute__((vector_size(32))); static v8sf sumfvhelper1(v8sf arr[4]) { v8sf retval = {0}; for (size_t i = 0; i < 4; i++) retval += arr[i]; return retval; } static float sumfvhelper2(v8sf x) { v8sf t = __builtin_shuffle(x, (v8u32){4,5,6,7,0,1,2,3}); x += t; t = __builtin_shuffle(x, (v8u32){2,3,0,1,6,7,4,5}); x += t; t = __builtin_shuffle(x, (v8u32){1,0,3,2,5,4,7,6}); x += t; return x[0]; } float sumfv(float *x) { //x = __builtin_assume_aligned(x, 64); v8sf *vx = (v8sf*)x; v8sf sumvv[4] = {{0}}; for (size_t i = 0; i < 2048/8; i+=4) { sumvv[0] += vx[i+0]; sumvv[1] += vx[i+1]; sumvv[2] += vx[i+2]; sumvv[3] += vx[i+3]; } v8sf sumv = sumfvhelper1(sumvv); return sumfvhelper2(sumv); } 

哪个gcc 4.8.4 gcc -Wall -Wextra -Wpedantic -std=gnu11 -march=native -O3 -fno-signed-zeros -fno-trapping-math -freciprocal-math -ffinite-math-only -fassociative-math -S变成:

 sumfv: vxorps %xmm2, %xmm2, %xmm2 xorl %eax, %eax vmovaps %ymm2, %ymm3 vmovaps %ymm2, %ymm0 vmovaps %ymm2, %ymm1 .L7: addq $4, %rax vaddps (%rdi), %ymm1, %ymm1 subq $-128, %rdi vaddps -96(%rdi), %ymm0, %ymm0 vaddps -64(%rdi), %ymm3, %ymm3 vaddps -32(%rdi), %ymm2, %ymm2 cmpq $256, %rax jne .L7 vaddps %ymm2, %ymm3, %ymm2 vaddps %ymm0, %ymm1, %ymm0 vaddps %ymm0, %ymm2, %ymm0 vperm2f128 $1, %ymm0, %ymm0, %ymm1 vaddps %ymm0, %ymm1, %ymm0 vpermilps $78, %ymm0, %ymm1 vaddps %ymm0, %ymm1, %ymm0 vpermilps $177, %ymm0, %ymm1 vaddps %ymm0, %ymm1, %ymm0 vzeroupper ret 

第二个辅助函数并不是绝对必要的,但是对向量元素的求和往往会在gcc中产生可怕的代码。 如果您愿意使用依赖于平台的内在函数,则可以使用__builtin_ia32_hadps256()替换大部分内部函数。