预运算是指将一个多次执行过程中需要反复用到的量储存在一个变量中,来避免多次的计算和调用。例如:
void set_row(double *a, double *b, long i, long n){
    long j;
    for(j = 0; j < n; j++)
        a[n*i+j] = b[j];
}
n*i这个值没有发生变化,但是在循环过程中发生了多次重复计算,所以我们改进为:
void set_row(double *a, double *b, long i, long n){
    long j;
    int ni = n * i;
    for(j = 0; j < n; j++)
        a[ni+j] = b[j];
}
甚至可以改进为指针风格的代码:
    long j;
    int ni = n * i;
    double *rowp = a+ni;
    for(j = 0; j < n; j++)
        *rowp++ = b[j];
我们可以将一些有预见性的运算进行优化:
for(i = 0; i < n; i++){
    int ni = n * i;
    for(j = 0; j < n; j++){
        a[ni + j] = b[j];
    }
}
显然,int ni = n * i;是一个没必要的乘法运算,因为可预见的,只需要每次加n即可。我们可以将其优化为加法运算:
int ni = 0;
for(i = 0; i < n; i++){
    for(j = 0; j < n; j++)
        a[ni + j] = b[j];
    ni += n;
}
ni += n;与原代码的乘法运算等效,但是使用加法运算还是可预测的优化了代码。
对于通用的表达式,我们也可以像上面一样,储存下来用于优化代码。以访问二维数组的上下左右为例:
up    = val[(i-1)*n + j  ];
down  = val[(i+1)*n + j  ];
left  = val[ i   *n + j-1];
right = val[ i   *n + j+1];
sum = up + down + left + right
显然,i*n+j的运算反复出现,我们可以优化为:
long inj = i*n + j;
up    = val[inj - n];
down  = val[inj + n];
left  = val[inj - 1];
right = val[inj + 1];
example:
void lower(char *s){
    size_t i;
    for(i = 0; i < strlen(s); i++)
        if(s[i] >= 'A' && s[i] <= 'Z')
            s[i] -= ('A' - 'a');
}
这段代码在字符串较大时,时间会以平方形式增加。原因在于:strlen函数是通过测试字符串是否到达末尾来计算字符串长度的,也就是自身带有O(n)的复杂度。
我们可以考虑,用len = strlen(s)提前储存字符串长度(因为我们只改变字符串内容而不改变长度),用来优化这个特定的程序。
但是对于编译器来说,一方面由于我们修改了字符串。编译器不得不小心谨慎地重新计算strlen(s);另一方面,编译器并不知道调用哪个版本的strlen函数。计算机必须假设strlen函数是一个黑箱,不能对其副作用做任何假定。所以编译器不会对此做任何优化。
另一个例子是将一维数组b的第i个数设置为二维数组a第i行的和:
void sum_row(double *a, double *b, long n){
    long i, j;
    for(i = 0; i < n; i++){
        b[i] = 0;
        for(j = 0; j < n; j++)
            b[i] += a[i*n + j];
    }
}
除了将i*n提前作为优化,观察编译代码,我们发现对于b[i]来说,每次都要将其从寄存器%xmm0中提取、加上一个数、再存回寄存器。为什么要在内存和寄存器之间反复传递数据?
原因在于我们无法确定在C中有没有内存别名使用。假设我们有如下代码:
double A[9] = 
    {  0,  1,  2,
       4,  8, 16,
      32, 64,128 };
double B[3] = A+3; // B = [4, 8, 16]
// 应该是 double *B = A + 3;
sum_row(A, B, 3);
此时,我们会发现,A的第二行和B使用了同一块内存,于是有了如下结果:
| 状态 | val of B | 描述 | 
|---|---|---|
| init | [4, 8, 16] | A的第二行和B使用了同一块内存 | 
| i = 0 | [3, 8, 16] | 正常求和,但此时 A[1][0]也被修改为3 | 
| i = 1 | [3, 22, 16] | A[1][1]=B[1]先被初始化为0,在j循环中,先是加上了已经被修改的A[1][0]=3,此时A[1][1]=B[1]=3,而B[1]又加上了A[1][1]等于6,最后6+16=22。 | 
| i = 2 | [3, 22, 224] | 正常求和,但是 A[1]也发生了改变,B的结果也不准确。 | 
编译器必须假设两个存储位置可能相互重叠,所以只能小心地一遍又一遍读取、修改、写入。比较好的方法是引入局部变量:
void sum_row2(double *a, double *b, long n){
    long i, j;
    for(i = 0; i < n; i++){
        double val = 0;
        for(j = 0; j < n; j++)
            val += a[i*n + j];
        b[i] = val;
    }
}
这样就避免了对内存的读写限制程序的性能。
需要对现代处理器设计有基本的了解:
硬件可以并行执行多条指令。
性能受到数据依赖性的限制:
指令之间的依赖关系可能会阻碍并行执行,从而影响性能。
简单的转换可显著提升性能:
编译器通常无法进行这些转换。 浮点运算中存在的结合律和分配律的缺失限制了某些优化。
CPE每个元素的周期数,它描述了处理器在执行特定计算任务时,每个元素所需的平均时钟周期数。CPE值越小,表示处理器执行相同计算任务时的效率越高。
在图像上一般反映为斜率。当然,有时候也会有不可避免的Overhead。
对于代码:
void combine1(vec_ptr v, data_t *dest){
    long int i;
    *dest = IDENT;
    for (i = 0; i < vec_length(v); i++) {
        data_t val;
        get_vec_element(v, i, &val);
        *dest = *dest OP val;
    }
}
只需要将编译器优化级别设置为1,就可以花费不优化时一半的时间。(当然数据类型和运算对时间几乎没有影响)
当然还可以进行更多的优化,例如:将 vec_length 移出循环、避免在每次循环中进行边界检查、在临时变量中累积等:
void combine4(vec_ptr v, data_t *dest){
    long i;
    long length = vec_length(v);
    data_t *d = get_vec_start(v);
    data_t t = IDENT;
    for (i = 0; i < length; i++)
        t = t OP d[i];
    *dest = t;
}
这样的优化结果可以将整数加法压缩到1个周期、乘法和浮点加法3个周期、浮点乘法5个周期(仅使用编译器优化时均为10个周期左右)。
至于为什么会有1、3、5出现,需要了解一些硬件结构。计算机会将指令按顺序执行,而且会提前读取一部分指令。当计算机发现后面的一部分指令与前文并没有依赖关系,就可以同时处理这些指令,即指令级并行性。
计算机的硬件会将指令拆分为低级操作,并确定相互之间的依赖关系。而一系列功能单元(in CPU)可以执行这些操作。
也就是说,当你的资源足够时,计算机可以并行这些代码,但这需要以某种方式构建程序,以使用这些资源。
这种CPU称为Superscalar Processor超标量处理器,可以在一个时钟周期内发射并执行多条指令。这些指令从一个顺序指令流中提取,通常以动态方式进行调度。
现代大多数CPU都是超标量处理器,在不需要额外编程工作的情况下,超标量处理器可以利用大多数程序中存在的指令级并行性。
在这些功能单元中,采用的是流水线处理技术,其基本思想是将计算分解为不同的阶段,当一个操作从一个阶段进行到下一个阶段时,上一个阶段就空出来,可以填入新的数据。
在单核中的并行性是低于多核并行性的,除了低功率嵌入式处理器,其他处理器均提供这种并行性。通过流水线实现的并行性是在资源有限下对时间的高效利用,不同于指令级并行性。
Haswell CPUP有多种功能不同的功能单元,可以同时处理:2次读内存、一次写内存、4次整数运算、2次浮点乘法、1次浮点加法和1次浮点除法。
在机器延迟的限制下,才会有1、3、5出现。这种限制是基于一个操作从开始到结束所需的时间。
循环展开也是一种优化方法,即增加每一次循环的步长,在循环体内补上这种操作,使其在前后功能完善不变的情况下减少循环次数,进而减少循环开销,借助CPU内部流水线来实现优化。
这种方法会导致代码量增加,但速度不一定快,反而可能因为icache存不下而导致指令访问变慢。
例如上述的函数中的循环体,我们可以做2x1循环展开:
for(i = 0; i < limit; i += 2){
    x = (x OP d[i]) OP d[i+1];
}
这种方法仅对整数加法优化了一点。但是用一种简单的方法就可以极大的优化这个程序 (Loop Unrolling with Reassoiation):
for(i = 0; i < limit; i += 2){
    x = x OP (d[i] OP d[i+1]);
}
这下运行时间从1、3、3、5改进到了1.01、1.5、1.5、2.5,是极大的优化。因为在这种优化下,d[i] OP d[i+1]和x OP ...是可以并行的,正好是之前运算结构的一半,CPE = D/2。
延迟界限是指在一系列操作严格顺序执行时,一条指令花费的 全部时间。而一个更加基本的界限,即吞吐量界限。这个限制基于硬件的数量和性能,基于功能单元的原始计算能力。例如:我可以同时进行4个整数乘法,但是可以读取内存的功能单元只有两个,此时性能受此限制。
吞吐量界限是CPE的最小界限。
除此之外,另一种循环展开(Loop Unrolling with Accumulators),将奇数位和偶数位运算相互独立:
for (i = 0; i < limit; i+= 2){
    x0 = x0 OP d[i];
    x1 = x1 OP d[i+1];
}
...
return x0 OP x1; 
这种优化下,运行时间达到了0.8、1.5、1.5、2.5,对整数加法表现出来优化。
最后,我们可以合并两种操作:
理念(Idea)
可以按照任意展开度(degree)L 展开循环。可以并行累积 K 个结果。L 必须是 K 的倍数。
限制(Limitations)
1.收益递减(Diminishing returns):
2.对于短长度会导致较大的开销(Large overhead for short lengths):
在浮点运算中常使用的寄存器%ymm是%xmm大小的两倍,而AVX512版本寄存器为512位,即64字节。其中有一些称为矢量加法的指令,例如:vaddsd %ymm0,%ymm1,%ymm1具有执行8次单精度浮点数加法(or 4次双精度浮点数)的效果。