当数组大小均匀时,为什么这段代码会变慢?

警告:实际上它不是由于2的幂而是奇偶校验.请参阅编辑部分.

我找到了一个显示相当奇怪行为的代码.

代码使用2D数组(大小x大小).当size为2的幂时,代码的速度在10%到40%之间,最慢的是size = 32.

我用英特尔编译器获得了这些结果.如果我使用gcc 5.4编译,2问题的能力消失但代码慢3倍.在不同的cpu上测试它,我认为它应该足够可重复.

码:

#define N 10000000

unsigned int tab1[TS][TS];

void simul1() {

  for(int i=0; i<TS; i++)
  for(int j=0; j<TS; j++) {

    if(i > 0)
      tab1[i][j] += tab1[i-1][j];

    if(j > 0)
      tab1[i][j] += tab1[i][j-1];
  }
}

int main() {

  for(int i=0; i<TS; i++)
  for(int j=0; j<TS; j++)
    tab1[j][i] = 0;


  for(int i=0; i<N; i++) {
    tab1[0][0] = 1;
    simul1();
  }

  return tab1[10][10];
}

汇编:

icc:
        icc main.c -O3 -std=c99 -Wall -DTS=31 -o i31
        icc main.c -O3 -std=c99 -Wall -DTS=32 -o i32
        icc main.c -O3 -std=c99 -Wall -DTS=33 -o i33

gcc:
        gcc main.c -O3 -std=c99 -Wall -DTS=31 -o g31
        gcc main.c -O3 -std=c99 -Wall -DTS=32 -o g32
        gcc main.c -O3 -std=c99 -Wall -DTS=33 -o g33

Xeon E5的结果:

time ./icc31
4.549s

time ./icc32
6.557s

time ./icc33
5.188s


time ./gcc31
13.335s

time ./gcc32
13.669s

time ./gcc33
14.399s

我的问题:

>当阵列大小正好是32时,为什么会出现性能损失?
>(可选:为什么icc比gcc快3倍?)

编辑:实际上这是由于奇偶校验,并且仅适用于大小> = 32.偶数和奇数之间的性能差异是一致的,并且当大小变大时减小.这是一个更详细的基准:

(y轴是以秒为单位的每秒元素数,以TS²* N / size / 1000000获得)

我的cpu一个32KB的L1缓存和256 KB的L2

解决方法

Why is icc 3x faster than gcc here?

GCC无法对内循环进行矢量化,因为它报告说,数据引用之间存在依赖关系.英特尔的编译器非常智能,可以将内部循环分成两个独立的部分:

for (int j = 1; j < TS; j++)
    tab1[i][j] += tab1[i-1][j];  // this one is vectorized

for (int j = 1; j < TS; j++)
    tab1[i][j] += tab1[i][j-1];

尝试1:循环重组

通过将simul1重写为以下内容,您可能会在GCC中获得更好的性能

void simul1(void)
{
    for (int j = 1; j < TS; j++)
        tab1[0][j] += tab1[0][j-1];

    for (int i = 1; i < TS; i++) {
        for (int j = 0; j < TS; j++)
            tab1[i][j] += tab1[i-1][j];
        for (int j = 1; j < TS; j++)
            tab1[i][j] += tab1[i][j-1];
    }
}

我在GCC 6.3.0下的结果是-O3 -march-native,TS = 32由Intel Core i5 5200U供电:

原始版本:

real    0m21.110s
user    0m21.004s
sys 0m0.000s

改性:

real    0m8.588s
user    0m8.536s
sys 0m0.000s

尝试2:前缀和的矢量化

在一些reaseach之后,我发现有可能通过向量加法和移位来向量化第二个内环.算法呈现here.

版本A:128位宽矢量

#include "emmintrin.h"

void simul1(void)
{
    for (int j = 1; j < TS; j++)
        tab1[0][j] += tab1[0][j-1];

    for (int i = 1; i < TS; i++) {
        for (int j = 0; j < TS; j++)
            tab1[i][j] += tab1[i-1][j];

        for (int stride = 0; stride < TS; stride += 4) {
            __m128i v;
            v = _mm_loadu_si128((__m128i*) (tab1[i] + stride));
            v = _mm_add_epi32(v,_mm_slli_si128(v,sizeof(int)));
            v = _mm_add_epi32(v,2*sizeof(int)));
            _mm_storeu_si128((__m128i*) (tab1[i] + stride),v);
        }

        for (int stride = 4; stride < TS; stride += 4)
            for (int j = 0; j < 4; j++)
                tab1[i][stride+j] += tab1[i][stride-1];
    }
}

结果:

real    0m7.541s
user    0m7.496s
sys 0m0.004s

版本B:256位宽矢量

这个更复杂.考虑整数的八元素向量:

V = (a,b,c,d,e,f,g,h)

我们可以将它视为两个打包向量:

(a,d),(e,h)

首先,算法执行两个独立的求和:

(a,h)
+
(0,a,c),(0,g)
=
(a,a+b,b+c,c+d),e+f,f+g,g+h)
+
(0,a+b),e+f)
=
(a,a+b+c,a+b+c+d),e+f+g,e+f+g+h)

然后它将第一个向量的最后一个元素传播到第二个向量的每个元素中,因此它最终产生:

(a,(a+b+c+d+e,a+b+c+d+e+f,a+b+c+d+e+f+g,a+b+c+d+e+f+g+h)

我怀疑,这些内在函数可以写得更好,因此有可能有所改进.

#include "immintrin.h"

void simul1(void)
{
    for (int j = 1; j < TS; j++)
        tab1[0][j] += tab1[0][j-1];

    for (int i = 1; i < TS; i++) {
        for (int j = 0; j < TS; j++)
            tab1[i][j] += tab1[i-1][j];

        for (int stride = 0; stride < TS; stride += 8) {
            __m256i v;
            v = _mm256_loadu_si256((__m256i*) (tab1[i] + stride));
            v = _mm256_add_epi32(v,_mm256_slli_si256(v,sizeof(int)));
            v = _mm256_add_epi32(v,2*sizeof(int)));

            __m256i t = _mm256_setzero_si256();
            t = _mm256_insertf128_si256(t,_mm_shuffle_epi32(_mm256_castsi256_si128(v),0xFF),1);

            v = _mm256_add_epi32(v,t);
            _mm256_storeu_si256((__m256i*) (tab1[i] + stride),v);
        }

        for (int stride = 8; stride < TS; stride += 8)
            for (int j = 0; j < 8; j++)
                tab1[i][stride+j] += tab1[i][stride-1];
    }
}

结果(Clang 3.8):

real    0m5.644s
user    0m5.364s
sys 0m0.004s

相关文章

本程序的编译和运行环境如下(如果有运行方面的问题欢迎在评...
水了一学期的院选修,万万没想到期末考试还有比较硬核的编程...
补充一下,先前文章末尾给出的下载链接的完整代码含有部分C&...
思路如标题所说采用模N取余法,难点是这个除法过程如何实现。...
本篇博客有更新!!!更新后效果图如下: 文章末尾的完整代码...
刚开始学习模块化程序设计时,估计大家都被形参和实参搞迷糊...