展开循环并使用矢量化进行独立求和
对于以下循环,如果我告诉它使用关联数学,例如使用-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()
替换大部分内部函数。