如何让GCC完全展开这个循环(即剥离这个循环)?

有没有办法指示GCC(版本我使用4.8.4)完全展开底部函数中的while循环,即剥离此循环? 循环的迭代次数在编译时是已知的:58。

我先解释一下我的尝试。

通过检查GAS输出:

gcc -fpic -O2 -S GEPDOT.c 

使用12个寄存器XMM0 – XMM11。 如果我将标志-funroll-loops传递给gcc:

 gcc -fpic -O2 -funroll-loops -S GEPDOT.c 

循环只展开两次。 我检查了GCC优化选项。 GCC表示-funroll-loops也会打开-frename-registers ,所以当GCC展开一个循环时,它先前选择的寄存器分配是使用“遗留”寄存器。 但是XMM12只剩下4个 – XMM15,所以GCC最多只能展开2次。 如果有48个而不是16个XMM寄存器可用,GCC将毫无困难地展开while循环4次。

然而我做了另一个实验。 我首先手动两次展开while循环,获得一个函数GEPDOT_2。 然后两者之间没有任何区别

 gcc -fpic -O2 -S GEPDOT_2.c 

 gcc -fpic -O2 -funroll-loops -S GEPDOT_2.c 

由于GEPDOT_2已用完所有寄存器,因此不执行展开。

GCC确实注册了重命名,以避免引入潜在的错误依赖。 但我确信在我的GEPDOT中没有这样的潜力; 即使有,也不重要。 我尝试自己展开循环,展开4次比展开2次更快,比没有展开更快。 当然我可以手动展开更多次,但这很乏味。 GCC可以帮我吗? 谢谢。

 // C file "GEPDOT.c" #include  void GEPDOT (double *A, double *B, double *C) { __m128d A1_vec = _mm_load_pd(A); A += 2; __m128d B_vec = _mm_load1_pd(B); B++; __m128d C1_vec = A1_vec * B_vec; __m128d A2_vec = _mm_load_pd(A); A += 2; __m128d C2_vec = A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; __m128d C3_vec = A1_vec * B_vec; __m128d C4_vec = A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; __m128d C5_vec = A1_vec * B_vec; __m128d C6_vec = A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; __m128d C7_vec = A1_vec * B_vec; A1_vec = _mm_load_pd(A); A += 2; __m128d C8_vec = A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; int k = 58; /* can compiler unroll the loop completely (ie, peel this loop)? */ while (k--) { C1_vec += A1_vec * B_vec; A2_vec = _mm_load_pd(A); A += 2; C2_vec += A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; C3_vec += A1_vec * B_vec; C4_vec += A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; C5_vec += A1_vec * B_vec; C6_vec += A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; C7_vec += A1_vec * B_vec; A1_vec = _mm_load_pd(A); A += 2; C8_vec += A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; } C1_vec += A1_vec * B_vec; A2_vec = _mm_load_pd(A); C2_vec += A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; C3_vec += A1_vec * B_vec; C4_vec += A2_vec * B_vec; B_vec = _mm_load1_pd(B); B++; C5_vec += A1_vec * B_vec; C6_vec += A2_vec * B_vec; B_vec = _mm_load1_pd(B); C7_vec += A1_vec * B_vec; C8_vec += A2_vec * B_vec; /* [write-back] */ A1_vec = _mm_load_pd(C); C1_vec = A1_vec - C1_vec; A2_vec = _mm_load_pd(C + 2); C2_vec = A2_vec - C2_vec; A1_vec = _mm_load_pd(C + 4); C3_vec = A1_vec - C3_vec; A2_vec = _mm_load_pd(C + 6); C4_vec = A2_vec - C4_vec; A1_vec = _mm_load_pd(C + 8); C5_vec = A1_vec - C5_vec; A2_vec = _mm_load_pd(C + 10); C6_vec = A2_vec - C6_vec; A1_vec = _mm_load_pd(C + 12); C7_vec = A1_vec - C7_vec; A2_vec = _mm_load_pd(C + 14); C8_vec = A2_vec - C8_vec; _mm_store_pd(C,C1_vec); _mm_store_pd(C + 2,C2_vec); _mm_store_pd(C + 4,C3_vec); _mm_store_pd(C + 6,C4_vec); _mm_store_pd(C + 8,C5_vec); _mm_store_pd(C + 10,C6_vec); _mm_store_pd(C + 12,C7_vec); _mm_store_pd(C + 14,C8_vec); } 

更新1

感谢@ user3386109的评论,我想稍微扩展这个问题。 @ user3386109提出了一个非常好的问题。 实际上,当有很多并行指令需要调度时,我确实对编译器的最佳寄存器分配能力有所怀疑。

我个人认为一种可靠的方法是首先在asm内联汇编中编码循环体(这是HPC的关键),然后根据需要多次复制它。 今年早些时候,我有一个不受欢迎的post: 内联汇编 。 代码有点不同,因为循环迭代次数j是一个函数参数,因此在编译时是未知的。 在这种情况下,我无法完全展开循环,所以我只复制了汇编代码两次,并将循环转换为标签并跳转。 事实certificate,我的书面汇编的结果性能比编译器生成的汇编高约5%,这可能表明编译器无法以我们预期的最​​佳方式分配寄存器。

我(并且仍然)是汇编编码的宝贝,所以这为我提供了一个很好的案例研究,让我对x86汇编有所了解。 但从长远来看,我并不倾向于对GEPDOT进行大规模的assembly编码。 主要有三个原因:

  1. asm内联汇编因不可移植而受到批评。 虽然我不明白为什么。 也许是因为不同的机器有不同的寄存器?
  2. 编译器也越来越好。 所以我仍然更喜欢算法优化和更好的C编码习惯,以帮助编译器产生良好的输出;
  3. 最后一个原因更为重要。 迭代次数可能并不总是58.我正在开发一个高性能矩阵分解子程序。 对于高速缓存阻塞因子nb ,迭代次数将是(nb-2) 。 我不打算将nb作为函数参数,就像我在之前的post中所做的那样。 这是一个特定于机器的参数,将其定义为宏。 因此,迭代次数在编译时已知,但可能因机器而异。 猜猜在手动循环展开各种nb时我需要做多少繁琐的工作。 因此,如果有一种方法可以简单地指示编译器剥离循环,那就太好了。

如果您还可以分享一些生产高性能,便携式库的经验,我将不胜感激。

尝试调整优化器参数:

 gcc -O3 -funroll-loops --param max-completely-peeled-insns=1000 --param max-completely-peel-times=100 

这应该可以解决问题。

这不是一个答案,但可能对其他试图用GCC进行矩阵乘法矢量化的人感兴趣。

下面,我假设c是行主要顺序的4×4矩阵, a是列主要顺序的4行n列矩阵(转置), b是行中的4列, n行矩阵- 主要顺序,计算的操作是c = a × b + c ,其中×表示矩阵乘法。

实现这一目标的天真function是

 void slow_4(double *c, const double *a, const double *b, size_t n) { size_t row, col, i; for (row = 0; row < 4; row++) for (col = 0; col < 4; col++) for (i = 0; i < n; i++) c[4*row+col] += a[4*i+row] * b[4*i+col]; } 

GCC使用SSE2 / SSE3生成相当好的代码

 #if defined(__SSE2__) || defined(__SSE3__) typedef double vec2d __attribute__((vector_size (2 * sizeof (double)))); void fast_4(vec2d *c, const vec2d *a, const vec2d *b, size_t n) { const vec2d *const b_end = b + 2L * n; vec2d s00 = c[0]; vec2d s02 = c[1]; vec2d s10 = c[2]; vec2d s12 = c[3]; vec2d s20 = c[4]; vec2d s22 = c[5]; vec2d s30 = c[6]; vec2d s32 = c[7]; while (b < b_end) { const vec2d b0 = b[0]; const vec2d b2 = b[1]; const vec2d a0 = { a[0][0], a[0][0] }; const vec2d a1 = { a[0][1], a[0][1] }; const vec2d a2 = { a[1][0], a[1][0] }; const vec2d a3 = { a[1][1], a[1][1] }; s00 += a0 * b0; s10 += a1 * b0; s20 += a2 * b0; s30 += a3 * b0; s02 += a0 * b2; s12 += a1 * b2; s22 += a2 * b2; s32 += a3 * b2; b += 2; a += 2; } c[0] = s00; c[1] = s02; c[2] = s10; c[3] = s12; c[4] = s20; c[5] = s22; c[6] = s30; c[7] = s32; } #endif 

对于AVX,GCC可以做得更好

 #if defined(__AVX__) || defined(__AVX2__) typedef double vec4d __attribute__((vector_size (4 * sizeof (double)))); void fast_4(vec4d *c, const vec4d *a, const vec4d *b, size_t n) { const vec4d *const b_end = b + n; vec4d s0 = c[0]; vec4d s1 = c[1]; vec4d s2 = c[2]; vec4d s3 = c[3]; while (b < b_end) { const vec4d bc = *(b++); const vec4d ac = *(a++); const vec4d a0 = { ac[0], ac[0], ac[0], ac[0] }; const vec4d a1 = { ac[1], ac[1], ac[1], ac[1] }; const vec4d a2 = { ac[2], ac[2], ac[2], ac[2] }; const vec4d a3 = { ac[3], ac[3], ac[3], ac[3] }; s0 += a0 * bc; s1 += a1 * bc; s2 += a2 * bc; s3 += a3 * bc; } c[0] = s0; c[1] = s1; c[2] = s2; c[3] = s3; } #endif 

使用gcc-4.8.4( -O2 -march=x86-64 -mtune=generic -msse3 )生成的程序集的SSE3版本本质上是

 fast_4: salq $5, %rcx movapd (%rdi), %xmm13 addq %rdx, %rcx cmpq %rcx, %rdx movapd 16(%rdi), %xmm12 movapd 32(%rdi), %xmm11 movapd 48(%rdi), %xmm10 movapd 64(%rdi), %xmm9 movapd 80(%rdi), %xmm8 movapd 96(%rdi), %xmm7 movapd 112(%rdi), %xmm6 jnb .L2 .L3: movddup (%rsi), %xmm5 addq $32, %rdx movapd -32(%rdx), %xmm1 addq $32, %rsi movddup -24(%rsi), %xmm4 movapd %xmm5, %xmm14 movddup -16(%rsi), %xmm3 movddup -8(%rsi), %xmm2 mulpd %xmm1, %xmm14 movapd -16(%rdx), %xmm0 cmpq %rdx, %rcx mulpd %xmm0, %xmm5 addpd %xmm14, %xmm13 movapd %xmm4, %xmm14 mulpd %xmm0, %xmm4 addpd %xmm5, %xmm12 mulpd %xmm1, %xmm14 addpd %xmm4, %xmm10 addpd %xmm14, %xmm11 movapd %xmm3, %xmm14 mulpd %xmm0, %xmm3 mulpd %xmm1, %xmm14 mulpd %xmm2, %xmm0 addpd %xmm3, %xmm8 mulpd %xmm2, %xmm1 addpd %xmm14, %xmm9 addpd %xmm0, %xmm6 addpd %xmm1, %xmm7 ja .L3 .L2: movapd %xmm13, (%rdi) movapd %xmm12, 16(%rdi) movapd %xmm11, 32(%rdi) movapd %xmm10, 48(%rdi) movapd %xmm9, 64(%rdi) movapd %xmm8, 80(%rdi) movapd %xmm7, 96(%rdi) movapd %xmm6, 112(%rdi) ret 

生成的程序集的AVX版本( -O2 -march=x86-64 -mtune=generic -mavx )基本上是

 fast_4: salq $5, %rcx vmovapd (%rdi), %ymm5 addq %rdx, %rcx vmovapd 32(%rdi), %ymm4 cmpq %rcx, %rdx vmovapd 64(%rdi), %ymm3 vmovapd 96(%rdi), %ymm2 jnb .L2 .L3: addq $32, %rsi vmovapd -32(%rsi), %ymm1 addq $32, %rdx vmovapd -32(%rdx), %ymm0 cmpq %rdx, %rcx vpermilpd $0, %ymm1, %ymm6 vperm2f128 $0, %ymm6, %ymm6, %ymm6 vmulpd %ymm0, %ymm6, %ymm6 vaddpd %ymm6, %ymm5, %ymm5 vpermilpd $15, %ymm1, %ymm6 vperm2f128 $0, %ymm6, %ymm6, %ymm6 vmulpd %ymm0, %ymm6, %ymm6 vaddpd %ymm6, %ymm4, %ymm4 vpermilpd $0, %ymm1, %ymm6 vpermilpd $15, %ymm1, %ymm1 vperm2f128 $17, %ymm6, %ymm6, %ymm6 vperm2f128 $17, %ymm1, %ymm1, %ymm1 vmulpd %ymm0, %ymm6, %ymm6 vmulpd %ymm0, %ymm1, %ymm0 vaddpd %ymm6, %ymm3, %ymm3 vaddpd %ymm0, %ymm2, %ymm2 ja .L3 .L2: vmovapd %ymm5, (%rdi) vmovapd %ymm4, 32(%rdi) vmovapd %ymm3, 64(%rdi) vmovapd %ymm2, 96(%rdi) vzeroupper ret 

我猜,寄存器调度不是最佳的,但它看起来也不是很糟糕。 我个人对上面的内容感到满意,而不是在此时尝试手动优化它。

在Core i5-4200U处理器(支持AVX2)上,上述function的快速版本在SSE3的1843个CPU周期(中位数)和AVX2的1248个周期内计算两个4×256矩阵的乘积。 每个矩阵条目下降到1.8和1.22个周期。 为了进行比较,未经过滤的慢速版本每个矩阵输入大约需要11个周期。

(循环计数是中值,即一半的测试速度更快。我只进行了大约100k重复的粗略基准测试,所以请将这些数字用一粒盐。)

在这个CPU上,缓存效果是这样的,4×512矩阵大小的AVX2仍然是每个条目1.2个周期,但是在4×1024时,它下降到1.4,在4×4096到1.5,在4×8192到1.8,并且每个条目4×65536到2.2个周期。 SSE3版本每个入口保持1.8个周期,最高可达4×3072,此时它开始减速; 在4×65536时,每个条目也约为2.2个周期。 我相信这个(笔记本电脑!)CPU此时的缓存带宽是有限的。