C++ 编译器优化与SIMD指令集

C++ 编译器优化与SIMD指令集

Reference:

  1. 【公开课】编译器优化与SIMD指令集(#4)
  2. 10分钟速览 C++20 新增特性
  3. 编译器实验用:https://godbolt.org/

相关文章:

  1. C++ 字节对齐
  2. 向量化运算 和 EIGEN_MAKE_ALIGNED_OPERATOR_NEW

本次的主题是编译器优化。编译器就是从源代码生成汇编语言。所以这里要从汇编语言的角度来理解编译器是怎么优化的,以及我们如何用好它。

1. 汇编语言

首先讲解下什么是汇编语言

下图为某 x86 的 64 64 64 位架构,它有以下寄存器:
在这里插入图片描述
寄存器在 CPU 内,它的速度比内存快。如果读到寄存器里,然后再计算就比较高效。所以上图的 x64 就提供了这么多的寄存器,其中白色的这些是 32 32 32 位模式下就有的;而 64 64 64 位以后,不仅把原来 32 32 32 位的这几个扩充到 64 64 64 位,它还额外追加了 8 8 8 个寄存器,这样我们用起来就更方便了。

图中 RIP 是它当前执行的代码的地址,也是扩充到了 64 64 64 位。还有 MMX、TMM、XMM,都是用于存储浮点数的寄存器。一个 XMM 有 128 128 128 位宽(YMM 256 256 256 位,后面还有 ZMM 512 512 512 位,适用于 AVX-512 寄存器),而 float 有 32 32 32 位宽,所以这里可以塞下四个 float 或者塞两个 double。这样它们在运算加法的时候就可以两个 double 一起来算,这样它效率就更高了。

1.1 通用寄存器:32位时代

开始说的 x86 架构中 32 32 32 位,只有以下 8 8 8 个寄存器:

  • eax, ecx, edx, ebx, esi, edi, esp, ebp

其中:

  • esp 是堆栈指针寄存器,和函数的调用与返回相关;
  • eax 是用于保存返回值的寄存器。

1.2 通用寄存器:64位时代

到了 64 64 64 位以后,又新增了 8 8 8 个寄存器。 64 64 64 位 x86 架构中的通用寄存器有:

  • rax, rcx, rdx, rbx, rsi, rdi, rsp, rbp, r8, r9, r10, r11, …, r15

其中:

  • r8 到 r15 是 64 64 64 位 x86 新增的寄存器,给了汇编程序员更大的空间,降低了编译器处理寄存器翻车(register spill)的压力;

寄存器多有什么好处呢?比如局部变量有 16 16 16 个,那它们就都能够存进寄存器里,而不需要存到内存上,这样它读写就更快了。因此 64 64 64 位比 32 32 32 位机器,除了内存突破 4 4 4GB 的限制外,也有一定性能优势。

1.3 8位,16位,32位,64位版本

从前两小节可以看到, 32 32 32 位的叫 eax,而到了 64 64 64 位就变成了 rax,它们之间有什么关系呢?它们是共用前 32 32 32 位的关系,就是 eax 这个地方的值和 rax 的低 32 32 32 位是共用的。同理还有 16 16 16 位的 ax 和 eax 共用的低 16 16 16 位,然后 ah 是 ax 的高 8 8 8 位、al 是低 8 8 8 位。

当然以 r 开头的这些数字编号寄存器也是一样的,它们通过 b、w、d、都没有 来区分的。
在这里插入图片描述

1.4 AT&T 汇编语言

汇编语言大致分为两种:

  1. 一种是 Intel 模式的汇编,它写寄存器就直接写 eax;写操作数,比如 eax 赋值给 edx,就把 edx 写在前面,eax 写在后面,然后用 mov 指令把 eax 赋给 edx。
  2. 另一种为 AT&T 汇编。它恰恰相反,它读的数字写在前面,写的数字写在后面,也就是说下图内的 movl %eax, %edx 是往右移动的,Intel 是往左移动的。

在这里插入图片描述

  • 操作数长度:Intel 就直接 mov,而 AT&T 后面要加一个 movl,代表是 l(long),也就是 32 32 32 位、b(byte) 代表的是 8 8 8 位、w(word) 代表的是 16 16 16 位。
  • 访问地址:Intel 是更直观的一个 [],然后后面接偏移量 lea ax, [dx + 0x10],而 AT&T 把偏移量写在 () 前面,后面来一个 (%dx) 代表它偏移的寄存器。

总之 AT&T 非常复杂,但是 GCC 用的就是这一种。

1.5 返回值:通过 eax 传出

比如说有一个函数 fun(),右图中 ret 是函数的返回指令,movl 是赋值指令。所以可以看到函数被编译成了把 42 42 42 复制给 eax 然后返回,就可以看出这个返回值是通过 eax 传出去的。

int func() {
    return 42;
}

在这里插入图片描述

这里是使用下面指令来编译的一个 cpp 文件,变成一个 .S 文件,-S 代表生成的不是 .o 文件,而是 .S 的汇编语言文件。

-fomit-frame-pointer 可以让生成代码更简洁。(作用为:如果函数不需要 frame pointer,就不要将 frame pointer 保留在寄存器中。当打开优化选项:-O, -O2, -O3, -Os 时或者对某些平台不打开任何优化选项时,-fomit-frame-pointer 会被默认打开

-fverbose-asm 让上图的每一条汇编前面有一个注释,来提示这一条指令代表的是源文件当中的第几行。比如 return 42 就代表了 movl $42 %eax 这一行。

gcc -fomit-frame-pointer -fverbose-asm -S main.cpp -o /tmp/main.S

编译器会自动设置寄存器优化,不需要手动指定。

1.6 前6个参数:分别通过 edi,esi,edx,ecx,r8d,r9d 传入

下面可以看到有一个有很多参数的函数,它有 6 6 6 个参数,它们通过寄存器传入的:
在这里插入图片描述

可以看到这里显示把所有传入的寄存器先给存到一个堆栈上面,rsp 就代表堆栈,而这个 - 就代表是堆栈上的某一个地址。它的返回值不是返回了 a 嘛,而 edi 它不是存到了 a 变量吗?它现在又取出来设为这个返回值,调用它的人一看 eax 就知道它返回了多少。

movl %edi, -4(%rsp) 相当于: *(rsp - 4) = edi

1.7 开启优化:-O3

对先前的指令稍加改进,通过加一个 -O3 选项,加了这个 flag 以后,它就比刚才更优化

gcc -fomit-frame-pointer -fverbose-asm -S main.cpp -o /tmp/main.S

刚才还把 b、c、d、e、f 都传到栈里了,现在直接简单了,直接把 a 传过去,就给送到 eax。它只需要一行指令就够了,这就是编译器的自动优化。它发现这几个存到栈,后来没使用过哇,所以就直接把这些指令全砍了。然后它又一想:欸?我这写入了一遍栈,又读了一遍栈,没意思啊,直接把 edi 复制给 eax 不就好了嘛,所以它就优化成了这样,movl %edi, %eax 相当于: eax = edi

int func(int a, int b, int c, int d, int e, int f) {
    return a;
}

在这里插入图片描述

所以开启这个选项之后,它就能够生成更精炼更高效的代码。

1.8 32位乘法运算:imull

除了可以 movl 来复制之外,汇编语言还支持乘法、加法之类的指令。比如乘法就是 imul,和 movl 是一样的,代表它的操作数是 32 32 32 位的。看看下面的这个例子做了些什么:

int func(int a, int b, int c, int d, int e, int f) {
   return a;
}

在这里插入图片描述
edi 我们知道代表的第一个参数 a,esi 代表的第二个参数 b。x86 汇编就是这样,它就是就地乘,没有写入到另一个寄存器,它是一种 CISC 指令的设计。这里首先把 edi 赋值给 eax,所以它就相当于,先把 eax 变成了 a,再把 b 乘到 a 里面,就变成所需的返回值了。

1.9 64位乘法运算:imulq

然后就是 64 64 64 位的乘法。之前提到了 l 代表 32 32 32 位,而 q 代表 64 64 64 位,所以这里的 long long 它是 64 64 64 位的。 64 64 64 位寄存器也从 e 寄存器变成 r 系列寄存器,就更宽了,能存储 64 64 64 位。

int func(int a, int b, int c, int d, int e, int f) {
    return a;
}

在这里插入图片描述

1.10 整数加法:被优化成 leal 了

这里的整数加法我们以为它会生成一个 addl,也就是 add,没想到它却变成了一个 lea 指令。

lea 是什么呢?可以看到 (%rdi,%rsi) 是一个地址的操作数,也就是这个地址其实相当于 *(rdi+rsi),而 lea 是加载表达式的地址,相当于 &,不是把表达式的值给取出来。它取的是这个地址,所以这个就相当于 eax = &*(rdi + rsi),也就是这两个会抵消,也就相当于 eax = rdi + rsi。所以说这是在妙用加载地址指令,把它当作 add 来操作,就可以更简练。
在这里插入图片描述

1.11 整数加常数乘整数:都可以被优化成 leal

leal 除了可以直接执行加法,它还可以执行任何一次函数。比如下面的 (%rdi,%rsi,8),就代表把第二个数乘以了 8 8 8。因为这种线性变换经常用于地址访问中,所以 x86 把它做进了地址访问的操作符内。同上一节,lea 却不读这个地址的值,而是将这个地址读出来。所以说这就变成了,我们可以把任何一个一次函数写在一条指令里。所以说一次函数在 x86 上都是很高效的。

在这里插入图片描述

1.12 指针访问对象:线性访问地址

那么线性变换为什么要做成指令呢?就是为了针对下面这种情况,比如 func(int* a, int*b),a 是一个指针,这里要返回 a[b],你可能会想 a 是怎么知道 a[b] 的地址,然后读取它的:
在这里插入图片描述

这里不能简单的 a+b,因为 a 是一个 int 类型的指针,如果直接在汇编里面 a+b,它会变成一个 char 类型的加。所以说这里要 char 类型的加,加上 sizeof(int) 乘以 b。因为 int 类型的大小是 4 4 4,所以这个偏移量也要乘以 4 4 4 才能访问到正确的地址。所以这里就用了 x86 提供的这个很方便的线性函数来访问。

这里的 moveslq 语句又是什么意思呢?因为我们刚刚使用的 int 是 32 32 32 位的,而指针和地址都是 64 64 64 位的,所以说要用这个 64 64 64 位的 a 去加上一个 32 32 32 位的 b 之前,要把 b 转换成 64 64 64 位的。可以看到这里是把 esi 这个 32 32 32 位变成 rsi 这个 64 64 64 位的了。这样以后才能进行线性访问。

1.13 指针的索引:尽量用 size_t

如果要避免上一节的转换呢?就可以使用 <sdtdint> 里提供的 std::size_t。这个 size_t 能够保证在 64 64 64 位系统上就是 int64;而 32 32 32 位系统上就是 int32,所以就不需要再进行刚才先从 32 32 32 扩展到 64 64 64,它直接就是 64 64 64 位的。所以说它有可能更搞笑,而且它还有一个好处:就比如数组的大小超过 INT_MAX,也就是数据大到超过 2 31 2^{31} 231,这时候 int 就不仅是效率问题,它是会出错的。所以我们就要用 64 64 64 位的 size_t 来表示超过了那个大小的索引。

在这里插入图片描述

所以建议如果用于下标的话,最好用 size_t,或者说循环体的 int i=0; i<xxx,也推荐使用 size_t。

1.14 浮点作为参数和返回:xmm 系列寄存器

高性能编程中经常会用到浮点数,所以说 CPU 也对它们做了专门的指令。比如 add 指令,刚才提到 addl 代表两个 32 32 32 位整数相加。这里的 a d d s s addss addsss 代表 single,也就是单精度浮点数,所以 a d d s s addss addss 就代表两个数相加。它传参数也不是 edi、esi 了,它直接用 xmm 系列,就是 xmm0 传 a;xmm1 传 b,然后我们把 xmm1 加到 xmm0,因为正好 xmm0 是用于返回值的,所以我们这里直接返回了,就得到 a+b 了。

在这里插入图片描述

1.15 什么是 xmm 系列寄存器?

刚才提到了 xmm 这个系列的寄存器,它们都有 128 128 128 位宽,可以容纳 4 4 4 个 float 或者 2 2 2 个 double。
在这里插入图片描述

刚才的例子中,因为只有一个 float 存在一个 128 128 128 位的寄存器内,所以只用到了它最低的 32 32 32 位。但是这样也没问题,因为我们刚才说的是,addss 它只会加最低位。这就要说到下一节将提到的 addss 了。

1.16 addss 是什么意思?

addss 可以拆成三个部分:add、s、s

  1. add 表示执行加法操作。
  2. 第一个 s 表示标量(scalar),只对 xmm 的最低位进行运算;也可以是 p 表示矢量(packed),一次对 xmm 中所有位进行运算
  3. 第二个 s 表示单精度浮点数(single),即 float 类型;也可以是 d 表示双精度浮点数(double),即 double 类型。

举个例子:

  • addss:一个 float 加法。
  • addsd:一个 double 加法。
  • addps:四个 float 加法。
  • addpd:两个 double 加法。

比如下图的 xmm0 和 xmm1:
在这里插入图片描述
addps 像左边这样同时计算四个的加;但如果是 addss 的话,它只会对最低位进行一个加法,而其他位都保留 xmm0 原来的值。

如果你看到编译器生成的汇编里,有大量 ss 结尾的指令则说明矢量化失败;如果看到大多数都是 ps 结尾则说明矢量化成功,这样生成代码就是比较高效的。

1.17 为什么需要 SIMD?单个指令处理四个数据

这种单个指令处理多个数据的技术称为 SIMD(single-instruction multiple-data),它可以大大增加计算密集型程序的吞吐量。

为什么我们需要这些指令,将四个打包到一起?这样速度更快。

因为 SIMD 把 4 4 4 个 float 打包到一个 xmm 寄存器里同时运算,很像数学中矢量的逐元素加法。因此 SIMD 又被称为矢量,而原始的一次只能处理 1 1 1 个 float 的方式,则称为标量

在一定条件下,编译器能够把一个处理标量 float 的代码,转换成一个利用 SIMD 指令的,处理矢量 float 的代码,从而增强你程序的吞吐能力!

通常认为利用同时处理 4 4 4 个 float 的 SIMD 指令可以加速 4 4 4 倍。但是如果你的算法不适合 SIMD,则可能加速达不到 4 4 4 倍;也有因为 SIMD 让访问内存更有规律,节约了指令解码和指令缓存的压力等原因,出现加速超过 4 4 4 倍的情况。

2. 化简

然后来看一看编译器还有哪些能做的优化。

2.1 编译器优化:代数化简

首先一件事是,要使它能够执行一些基本代数化。比如看下图左边的例子,这里 c c c a + b a+b a+b d d d a − b a-b ab,那么它们相加起来是不是 b b b 就被抵消了?所以编译器很聪明,它直接把 a a a 作为返回值。
在这里插入图片描述

2.2 编译器优化:常量折叠

如果我看到 a a a b b b 都是常量,然后它们还相加起来。我不是已经知道两个常量了,把它们拿过来,然后相加变成 42 42 42,它就直接优化成 return 42 了。

在这里插入图片描述

2.3 编译器优化:举个例子

它更疯狂的是,甚至可以弄一个 for 循环,它看到初始值是一个常数,而 i 每次的值也都是常数,那么它可以把 1 1 1 加到 100 100 100 给优化成直接返回 5050 5050 5050 都没有问题。
在这里插入图片描述

2.4 编译器优化:我毕竟不是万能的

当然编译器也没那么聪明,就是如果把代码写的很复杂,比如用了一些 STL 容器库。当它在看到 arr.push_back() 的时候,就停止思考了。这时它就会放弃优化,干脆去一个个调用 vector。所以说像 vector 这种会在堆上分配内存的容器,编译器通常是不会去进行优化的。所以说尽量把代码写的简单一点,编译器能够看得懂,它才能帮你优化;但如果用的很复杂的话,编译器就认为优化很复杂的抽象语法树需要花很长时间,编译器就觉得不划算,于是放弃优化-----不要把答案藏的太深,它的搜索范围是有限度的,它不会无限制的搜索。所以要把代码简化,从而它的搜索时间就能变短,就能找到正确的优化。

在这里插入图片描述
结论:尽量避免代码复杂化,把代码写的简单一点,避免使用会造成 new/delete 的容器。简单的代码,比什么优化手段都强。

2.5 造成 new/delete 的容器:即内存分配在堆上的容器

下表可以看到哪些内存分配在堆上,哪些又在栈上:

存储在堆上(妨碍优化)存储在栈上(利于优化)
vector, map, set, string, function, anyarray, bitset, glm::vec, string_view
unique_ptr, shared_ptr, weak_ptrpair, tuple, optional, variant

比如 vectormapset 这种,一般来说这些老的容器都是在堆上的。

为什么要区分堆和栈呢,全部存在栈上不好嘛?事实是有区别的:比如说 array,它是固定大小的,这是它的一个缺点,但是它的优点是它可以存在栈上。那为什么 vector 不能被存在栈上?因为它可以动态的 push_back(),也就是它的大小是变化的,而栈上只能一次性开辟一定的大小,而不能动态扩充大小。所以这时候像 mapsetstring 这种可以动态扩展的,得存储在堆上。像 unique_ptr 这种,它肯定是去调用了 new 和 delete 的,所以是在堆上。 (vector 这种肯定也是存在堆上的)

2.6 那改用 array 试试?

既然 array 存储在堆上,那么将之前的示例改写成 array 的,是不是就能优化成功了呢?
在这里插入图片描述
从图中可以看到,还是优化失败。

2.6.1 那改用手写的 reduce?

这里再改一下,改用手写的 reduce,它还是优化失败:
在这里插入图片描述

2.6.2 那改小到 10?成功了!

那么再修改一下,把这里的 100 100 100 改成 10 10 10,它优化成功了:
在这里插入图片描述
这是因为如果代码过于复杂,涉及的语句数量过多时,编译器会放弃优化

还是那句话,如果代码很简单,那么其实不需要什么优化手段,编译器能够理解你在干什么,那么它就能优化了。

2.7 constexpr:强迫编译器在编译期求值

当然上面的例子中,如果的确想要 100 100 100 这么大的能自动在编译器求值,可以考虑使用 constexpr 语法。这时候这个函数强制编译器在运行时求值,这会让编译器花费较长的时间,即编译变慢。所以说如果想迫使编译器优化,即进行常量折叠,就使用 constexpr 关键字。

不过,constexpr 函数中无法使用非 constexpr 的容器,比如前面提到的 vectormapsetstring。因为要让编译器理解全部,它必须全部分配在栈上,而不能分配在堆上。据说 C++20 开始 vector 可以在编译期分配了,newdelete 也可以 constexpr 了 (10分钟速览 C++20 新增特性)。但是 C++17 还是不行的,C++17 的 constexpr 还是只能用 std::array 这些在栈上的容器。
在这里插入图片描述你甚至可以把它作为一个模板,即使这个数再大,它也会去算,但是会让编译变得很慢,因为下面示例中的 50000 50000 50000 次迭代实在编译期进行的。
在这里插入图片描述

3. 内联

3.1 调用外部函数:call 指令

函数分为两种:外部的函数、内部的函数。

什么是外部函数?就是像下面这种,只要一个申明,然后实际的实现在另一个文件的这种就叫外部函数。这种编译期没办法优化,它只能生成一个 call 指令 来调用这个函数:
在这里插入图片描述

@PLTProcedure Linkage Table 的缩写,即函数链接表。链接器会查找其他 .o 文件中是否定义了 _Z5otheri 这个符号,如果定义了则把这个 @PLT 替换为它的地址。可以看到 PLT 在 call _Z5otheri@PLT 出现,代表这个函数是外部函数。它会在链接的时候,如果另一个文件定义了 _Z5otheri 这个符号,它就会把那个符号的地址填充到这里。所以如果别的函数里没有定义 _Z5otheri,它这里就没办法填充,从而链接器就会报错。

gcc -fomit-frame-pointer -fverbose-asm -S main.cpp -o /tmp/main.S

3.1.1 编译器优化:call 变 jmp

上面一小节是没开优化的结果,开优化以后就没有 call 指令了。它直接跳转到 jmp 的 other 这个地址。在前面一小节,调用完 func 一样是要返回的,不如让 other 代为返回,这样也是可以的
在这里插入图片描述

gcc -O3 -fomit-frame-pointer -fverbose-asm -S main.cpp -o /tmp/main.S

3.2 多个函数定义在同一个文件中

刚才提到外部函数会让编译器无法优化,但如果是内部函数呢?内部函数是声明和定义在同一个文件,就是它定义在 func 调用它的相同文件。编译器在编译 func 的同时,它是看得到 other 的定义的。也就是说它知道 other 里面什么都没做,直接返回了 a a a,从而它就可以优化了。下图内上面一个红框为 other 的函数,下面的是 func 的函数。
在这里插入图片描述也就是说,这样的话,它本来不优化的话是会直接去调用的,而开了优化之后,定义在同一个文件里是没有 @PLT 的,这样链接器也不用到时候把它替换了。

gcc -fomit-frame-pointer -fverbose-asm -S main.cpp -o /tmp/main.S

3.2.1 编译器优化:内联化

如果在同一个文件里,而且开启了 -O3,这时候我们发现 func() 根本没有调用 other(),它直接返回了 233,这是为什么呢?
在这里插入图片描述
因为编译器看到 other 就是直接返回 a a a,所以这里就直接替换成了 return 233 了。

这就是内联:当编译器看得到被调用函数(other)实现的时候,会直接把函数实现贴到调用它的函数(func)里。

没使用关键字 inline 而开启了 -O3 也会内联,编译器只要看得到这个函数体就行了,它自动帮你 inline,不需要去指定。

gcc -O3 -fomit-frame-pointer -fverbose-asm -S main.cpp -o /tmp/main.S

注意:

  1. 只有定义在同一个文件的函数可以被内联!否则编译器看不见函数体里的内容是没法内联的。
  2. 为了效率我们可以尽量把常用函数定义在头文件里,然后声明为 static(这里的重点是要声明为 static)。这样调用它们的时候编译器看得到它们的函数体,从而有机会内联。

3.2.2 局部可见函数:static

如果 other 声明为 static,就代表只有这文件可见,编译器就干脆不生成 other 函数了,这样的话还是可以内联。如果不是 static,就像上一小节那样,它也是可以内联的,不过会将额外生成的 other 暴露出来。
在这里插入图片描述
与上一小节相比,可以看见编译器内少了 _Z5otheri 的定义。这是因为因为 static 声明表示不会暴露 other 给其他文件,而且 func 也已经内联了 other,所以编译器干脆不定义 other 了。

gcc -O3 -fomit-frame-pointer -fverbose-asm -S main.cpp -o /tmp/main.S

3.2.3 inline 关键字?不需要!

在现代编译器中,它们都足够只智能,知道哪些函数要内联,不需要去提醒它,下图为加与不加 inline 的代码:
在这里插入图片描述
可以看到提醒前和提醒后的代码是完全一样的,编译器都帮着内联了。

结论:
在现代编译器的高强度优化下,加不加 inline 无所谓。编译器不是傻子,只要它看得见 other 的函数体定义,就会自动内联。
内联与否和 inline 没关系,内联与否只取决于是否在同文件,且函数体够小。要性能的,定义在头文件声明为 static 即可,没必要加 inline。static 纯粹是为了避免多个 .cpp 引用同一个头文件造成冲突,并不是必须 static 才内联。如果你不确定某修改是否能提升性能,那你最好实际测一下,不要脑内模拟。
inline 在现代 C++ 中有其他含义,但和内联没有关系,它是一个迷惑性的名字。

在线做编译器实验推荐这个网站:https://godbolt.org/,比如在下面写了一个平方函数,可以在右侧看到它自动变成了图上所示的这些指令。
在这里插入图片描述当然也可以在 Compiler options 处加上优化开关,可以看到它就优化成下面这个样子了:
在这里插入图片描述

4. 指针

4.1 指针别名现象(pointer aliasing)

我们知道 C 语言的历史遗产就是指针。下面弄一个 func 函数,可以看到:它将 a a a 指向的值取出来,然后写入到 c c c,再把 b b b 指向的值写入到 c c c。那你肯定想:先把 a a a 复制到 c c c,再把 b b b 复制到 c c c,那第一个语句不就浪费了嘛。我直接将 b b b 复制给 c c c 不就行了,因为 a a a 的值被覆盖了。那么为啥编译器明明是开了优化它还是复制了两遍呢?
在这里插入图片描述考虑这种调用方式, b b b c c c 指向同一个变量:
在这里插入图片描述优化前相当于: b = a ; b = b ; b=a; b=b; b=a;b=b; 最后 b b b 变成了 a a a
如果优化了: b = b ; b=b; b=b; 这样 b b b 就没有发生改变,这就是编译器放弃优化的原因

编译器并不知道 a a a b b b 是否指向了同一个地址,所以如果它优化掉的话,结果就不对了。

编译器是保守的,它宁愿变慢也不要出错。

4.1.1 告诉编译器别怕指针别名:__restrict 关键字

如果你想让编译器放心大胆,告诉它 b b b 永远不会等于 c c c 的,可以这样写:
在这里插入图片描述

这时 gcc 特有的 __restrict 关键字,它在 C语言里是标准,而在 C++ 里却没有标准,所以需要加上 __ 前缀代表它是当前编译器特定的。__restrict 是一个提示性的关键字,是程序员向编译器保证:这些指针之间不会发生重叠。所以 __restrict 提示以后它不会再重复两遍了, 它只会重复一遍。也就是说,向编译器保证 b b b c c c 是不同的值,从而它就可以安全的把 *c=*a 给优化掉

4.1.2 __restrict 关键字:只需加在非 const 的即可

那是不是在程序里每一个出现指针的地方都得加个 __restrict 呢?不一定,只需要非 const 的加一个就行了。因为只有写入才会造成指针别名的问题,而如果只读的话就没有任何影响。即:所有非 const 的指针都声明 __restrict
在这里插入图片描述

4.1.3 禁止优化:volatile

除此之外,还有一个让编译器不要优化的关键字,就是 volatile
在这里插入图片描述

可以看到上面示例中,先是把 a a a 赋值为 42 42 42,然后有返回了 *a。也就是说,这时候先把 42 42 42 写入,然后这里直接被编译器优化成 return 42 了。因为它想刚刚我写入的 42 42 42 嘛,但是万一这个 a a a 指针是指向一个多线程同时在写入的地方,也就是比如 *a=42 我赋值 42 42 42 了以后,又有某个线程进来,把 a a a 修改了,这时候这个返回值就可能是不一样的。所以这时候可以告诉编译器不要优化这个变量,就是 int volatile

在这里插入图片描述这时候它就会硬生生去写一遍再读一遍,生怕函数内再插入了什么东西。如果想优化的话,就不要加 volatile 了。如果在使用多线程后想要测试这个内存读写有多快,就可以使用关键字 volatile,它不会优化读写,从而可以测试它的性能。

4.1.4 注意一下区别

  • volatile
    volatile int *a 或 int volatile *a

  • __restrict
    int *__restrict a

  • 语法上区别:volatile* 前面而 __restrict* 后面。

  • 功能上区别:volatile 是禁用优化,__restrict 是帮助优化。

  • 是否属于标准上区别:volatileconst 一样是 C++ 标准的一部分。restrict 是 C99 标准关键字,但不是 C++ 标准的关键字。__restrict 其实是编译器的“私货”,好在大多数主流编译器都支持。

  • volatile 是可以针对不是指针的,而 __restrict 必定是针对指针。

4.2 编译器优化:合并写入

说道指针,还有一点就是,下面示例中有两个写入,第一个是针对 a[0],然后是 a[1],这两个分别是 int,也就是 32 32 32 位的写入。在编译出来之后就变成了一个 q,也就是变成了一个 64 64 64 位的写入。因为这里的两个 32 32 32 位是连续的地址,所以它能够自动优化成一个 64 64 64 位的写入:
在这里插入图片描述

4.2.1 合并写入:不能跳跃

刚才说的写入能够变成单独一个 64 64 64 位的,但如果地址中间跳了,变成了 a[0]a[2],中间 a[1] 空了,那它没办法优化了,只能变成写入两个 32 32 32 位的了。
在这里插入图片描述所以说我们在设计数据结构的时候,尽可能把它设计的紧凑点,不要把很多无关的变量交错排列,这样它就很难优化了。

5. 矢量化

5.1 更宽的合并写入:矢量化指令(SIMD)

前面说道可以将 2 2 2int 变成一个int64,其实还可以把 4 4 4int 优化成一个 __m128,也就是 xmm 寄存器。比如下图中有 4 4 4 个 int,因为 4 4 4 32 32 32 位就是 128 128 128 位,所以就可以用 xmm 来实现一次写入 4 4 4 个 int:
在这里插入图片描述xmm 是 SSE 引入的,它是 128 128 128 位的。它本来是用于存储 4 4 4 个 float 的,但它也可以存 4 4 4 个 int。

这里还有一个 movups:move unaligned packed single,之前都是没有这个 u 的,这个 u 是特有的,它代表上面 %rdi 这个地址不一定对齐到 16 16 16 字节。也可以是 movaps:move aligned packed single,这时候就是向它保证 %rdi 一定是对齐到 16 16 16 字节的,这样可能有微乎其微的提升。

5.1.1 SIMD 指令:敢不敢再宽一点?

刚才说的是 4 4 4 个可以变成 128 128 128,然后 8 8 8 个就可以变成 256 256 256,这个是 ymm 寄存器。但我为什么明明像下面这样写,它却用了两个 xmm 的寄存器,而没有优化成一个 ymm 呢?
在这里插入图片描述
这是因为编译器不敢保证你的电脑支持 AVX。也就是说编译器它只能生成、只能保证因为所有 64 64 64 位的电脑都支持 xmm,所以它只能保证 cmm 是能用的,而 256 256 256 位的 ymm 只有一些比较新的电脑上才有,它就不敢用了。

5.1.2 让编译器自动检测当前硬件支持的指令集在这里插入图片描述

所以说你要让它敢用,可以加以下 flag,跟它说,你来检测一下我的电脑是否支持 AVX,如果我的电脑支持,那你就可以用 ymm 优化了:

gcc -march=native -O3

比如上面的示例中写了 8 8 8 个,编译器也的确用 ymm 读取了一次,然后写入。这样就更加高效了。

-march=native让编译器自动判断当前硬件支持的指令。不过注意这样编译出的程序,可能放到别人不支持 AVX 的电脑上没法运行。

所以如果用于发布到特定平台的话,最好还是不要用这个 flag。

5.2 数组清零:自动调用标准库的 memset

下面为数组清零的代码,可以看到代码内除了清零什么都没做,也没有调任何标准库函数,而在右边汇编内,它缺变成了一个对 memset 的调用,这是为什么呢?
在这里插入图片描述
这是因为编译器是和标准库绑定的。所以说在这里生成一个像是把整个数组清零操作,它会识别到,然后就能够一个对标准库的调用了。这样就不必为了高效,手动改写成对 memcpy/memset 的调用,影响可读性。编译器会自动分析是在做拷贝还是清零,并优化成对标准库 memcpy/memset 的调用,标准库里面的实现通常是更高效的。

5.3 从 0 到 1024 填充:SIMD 加速

现在有一个从 0 0 0 1024 1024 1024,然后填入到数组 a[] 中:

void func(int *a) {
    for (int i = 0; i < 1024; i++) {
        a[i] = i;
    }
}

它的汇编如下图所示:
在这里插入图片描述
其中 paddd 表示 4 4 4 个 int 的加法;movdqa 表示加载 4 4 4 个 int。

是不是挺难看懂的。实际上它就是把代码优化成下面这个方便理解的伪代码:

void func(int *a) {
    __m128i curr = {0, 1, 2, 3};
    __m128i delta = {4, 4, 4, 4};
    for (int i = 0; i < 1024; i += 4) {
        a[i : i + 4] = curr;
        curr += delta;
    }
}

这里 i++ 变成 i+=4 了,就是一次可以写入 4 4 4 个,然后我们这里在 a [ i ] a[i] a[i] 写入 0 , 1 , 2 , 3 {0,1,2,3} 0,1,2,3 这四个数。在下一步的时候把 curr 给全部加上 4 4 4,curr 就变成 4 , 5 , 6 , 7 4,5,6,7 4,5,6,7 了;这时候再写入就是 4 , 5 , 6 , 7 4,5,6,7 4,5,6,7,再加一下,就是 8 , 9 , 10 , 11 8,9,10,11 8,9,10,11。这样和原来的运行结果完全一致,但却能更加高效,因为它是用 SIMD 指令并行去加的。

5.4 如果不是 4 的倍数?边界特判法

上面的示例中有一个缺点,因为这里是 i+=4,如果上面不是 1024 1024 1024 而是 1023 1023 1023 的话,那它这里就会多写入一个 int,而这个地址可能是别人的数据就把它覆盖掉,这该怎么办呢?

void func(int *a, int n) {
    for (int i = 0; i < n; i++) {
        a[i] = i;
    }
}

可以用边界特判法。刚才是一个固定的 1024 1024 1024,编译器看得到,它是一个 4 4 4 的倍数,可以直接像上一小节那样写。但在边界不是 4 4 4 的倍数的时候,它就生成了一个超复杂的代码:
在这里插入图片描述
它做的事是,假如 n 是 1023 1023 1023 的话,它会先把 n 编程一个刚才优化过的 1020 1020 1020 个元素的填入,这样每次可以处理四个。然后还剩下三个是没办法直接一次写入了,它就变成了传统的标量模式,来一个处理一个。所以这种思想就是,如果 n 不是 4 4 4 的倍数,就先取其中和 4 4 4 整除的部分,然后剩下不整除的部分,再用低效的标量代码去处理。这样大部分 1020 1020 1020 个数据都是矢量化的加,而剩下的 3 3 3 个数据才是标量加。所以它既高效又不会犯错,这就是边界特判。

编译器做优化时,它会自己去判断。但如果要自己手动去写 mm 什么的 SIMD 指令的话,需要手动去判断这个边界。

5.4.1 n 总是 4 的倍数?避免边界特判

如果能保证 n n n 总归是 4 4 4 的倍数,不想让它这么复杂的生成两个版本,可以这样写:

void func(int *a, int n) {
    n = n / 4 * 4;
    for (int i = 0; i < n; i++) {
        a[i] = i;
    }
}

用 n 先除以 4 4 4,因为除是整除,就让 n 会落入 4 4 4 的倍数中,这样的话编译器能够发现 n 始终是 4 4 4 的倍数,即 n % 4 = 0 n \% 4=0 n%4=0,从而它就不会生成边界特判的分支了,相当智能:
在这里插入图片描述

5.5 假定指针是 16 字节对齐的:assume_aligned

GCC 内置的函数,都是以两个下划线 __,然后 builtin 开头的。(int *)__builtin_assume_aligned(a, 16) 是在告诉编译器,保证 a a a 16 16 16 组合对齐的。

void func(int *a, int n) {
    n = n / 4 * 4;
    a = (int *)__builtin_assume_aligned(a, 16);
    for (int i = 0; i < n; i++) {
        a[i] = i;
    }
}

来看看这与上面的示例有什么区别:
在这里插入图片描述
先前这里是 movups,这里变成了 movaps。这就让 CPU 知道,保证是 16 16 16 字节对齐的,传 CPU 可能更快一点,也可能不快。但是这样的话就是调用 GCC 才有的函数,如果在微软编译器上,就不能通过。

所以为了通用,C++20 在 <memory> 头文件内引入了这个模板函数 std::assume_aligned,然后通过模板参数来制定对齐到 16 16 16 字节:

#include <memory>
void func(int *a, int n) {
    n = n / 4 * 4;
    a = std::assume_aligned<16>(a);
    for (int i = 0; i < n; i++) {
        a[i] = i;
    }
}

5.9 数组求和:reduction 的优化

float func(float *a) {
    float ret = 0;
    for (int i = 0; i < 1024; i++) {
        ret += a[i];
    }
    return ret;
}

优化成了下面这个样子:
在这里插入图片描述

伪代码如下:

#include <x86intrin.h>

float func(float *a) {
    __m128 ret = _mm_setzero_ps();
    for (int i = 0; i < 1024; i += 4) {
        __m128 a_i = _mm_loadu_ps(&a[i]);
        ret = _mm_add_ps(ret, a_i);
    }
    float r[4];
    _mm_storeu_ps(r, ret);
    return r[0] + r[1] + r[2] + r[3];
}

总之就是非常高效。

6. 循环

我们的程序中有 热 和 冷 这两个术语,这里的热是指经常被访问的代码。比如上面一节中 float ret=0,它只会调用一次,而这个 for 循环里面被调用了 1024 1024 1024 次。这个调用次数多的地方称为热,而调用数少的地方称为冷。一般都是在热的地方优化,因为热的地方通常它执行很多遍,更容易成为瓶颈,而往往循环的地方都能执行很多遍。所以说优化通常都是针对循环的优化,而循环也往往是优化的重点。

6.1 循环中的矢量化:还在伺候指针别名

这里提供了一个用循环的 func 函数,它的工作就是把 b[] 里面的数 + 1 +1 +1,然后复制到 a[]

void func(float *a, float *b) {
    for (int i = 0; i < 1024; i++) {
        a[i] = b[i] + 1;
    }
}

来看看这会发生什么:
在这里插入图片描述
按照我们的理解,这个代码应该是可以矢量化的,但是编译器却生成了两个版本,这是为什么?编译器还在担心前面章节提到的, a a a b b b 只想的数组是否有重合。考虑 func(a,a+1) 的情况,那样会产生数据依赖链,没法 SIMD 化。

为了优化而不失正确性,编译器索性生成两份代码。 一份是 SIMD 的,一份是传统标量的:
它在运行时检测 a a a, b b b 指针的差是否超过 1024 1024 1024 来判断是否有重叠现象。

  1. 如果没有重叠,则跳转到 SIMD 版本高效运行。
  2. 如果重叠,则跳转到标量版本低效运行,但至少不会错。

6.1.1 循环中的矢量化:解决指针别名

所以我们为了让编译器放心,会加上之前提到的关键字 __restrict,告诉编译器这两个数组保证不会重合,打消编译器的顾虑:

void func(float *__restrict a, float *__restrict b) {
    for (int i = 0; i < 1024; i++) {
        a[i] = b[i] + 1;
    }
}

这样它就只生成了一个 SIMD 版,而不需要在运行时判断是不是有重合了:
在这里插入图片描述
这样的话,因为不需要分支它效率就能提升,代码更为合理。

6.2 循环中的矢量化:OpenMP 强制矢量化

#pragma 的作用是设定编译器的状态或者是指示编译器完成一些特定的动作 。#pragma 指令对每个编译器给出了一个方法,在保持与 C 和 C++ 语言完全兼容的情况下,给出主机或操作系统专有的特征。依据定义,编译指示是机器或操作系统专有的,且对于每个编译器都是不同的。

除了可以用 GCC 特有的 __restrict 让编译器放心做 SIMD,也可以使用比较跨平台的 OpenMP:

void func(float *a, float *b) {
#pragma omp simd
    for (int i = 0; i < 1024; i++) {
        a[i] = b[i] + 1;
    }
}

OpenMP 是高性能计算的一个框架,它可以通过 -fopenmp 来打开:

gcc -fopenmp -O3 -fomit-frame-pointer -fverbose-asm -S main.cpp -o /tmp/main.S

使用起来只需要用 #pragma omp simd 语句就行了。它会告诉你:下面代码有可能出现数据依赖,但是不要担心。这样编译器就知道这里都显式的说 SIMD 了,那应该能默认这两个指针不会打架,从而它也是只会生成一份代码:
在这里插入图片描述

6.3 循环中的矢量化:编译器提示忽略指针别名

除了可以用上面的两种方法外,如果只用 GCC,还可以用 ivdep 这个 pragma,它是 ignore vector dependency 的简称,即忽略矢量依赖关系,表示忽略下方 for 循环内可能的指针别名现象。就是说在下面的代码中,它可以忽略 b b b 可能依赖于 a a a 的情况,这个是 GCC 特有的:

void func(float *a, float *b) {
#pragma GCC ivdep
    for (int i = 0; i < 1024; i++) {
        a[i] = b[i] + 1;
    }
}

其他编译器这个 pragma 可能不一样,需要自己去查,这里只是拿 GCC 举例。

6.4 循环中的 if 语句:挪到外面来

还有下面这种循环的优化。这里作者的用意很明显,在 is_mul 为真时执行 a*=b,否则执行 a+=b

void func(float *__restrict a, float *__restrict b, bool is_mul) {
    for (int i = 0; i < 1024; i++) {
        if (is_mul) {
            a[i] = a[i] * b[i];
        } else {
            a[i] = a[i] + b[i];
        }
    }
}

然而有 if 分支的循环体是难以 SIMD 矢量化的。所以说编译器会把这 if 语句挪到外面,具体是过程如下:
在这里插入图片描述
可以看到 corl %eax, %eax 内的 eax,是参数 is_mul 的寄存器,它会先判断这个是 true 还是 false。如果为 true 就跳到 L2,可以看到生成 addps 的一个循环,反之跳到 L3,可以看到 mulps,也就是乘法。

它会去提前判断,也就是它相当于进行了这么一个优化:

void func(float *__restrict a, float *__restrict b, bool is_mul) {
 	if (is_mul) {
    	for (int i = 0; i < 1024; i++) {
            a[i] = a[i] * b[i];
       }
    } else {
  		for (int i = 0; i < 1024; i++) {
            a[i] = a[i] + b[i];
        }
    }
}

也就是编译器看到 is_mul 是一个常量,把 if 挪到了 for 外面,相当于生成了两个版本:一个乘法、一个加法。这样就可以里面没有 if 语句,从而可以矢量化了。这样写起来比较麻烦,所以有时候也像前面那样写。只要编译器知道 is_mul 是循环中不会改变的量,它就能够这样优化。

6.5 循环中的不变量:挪到外面来

除了上面可以挪到外面的 if,还有可以挪到外面的表达式。比如优化之前是长下面这样的:

void func(float *__restrict a, float *__restrict b, float dt) {
    for (int i = 0; i < 1024; i++) {
        a[i]=a[i]+b[i]*(dt*dt);
    }
}

优化之后它发现,dt*dt 在循环中是不会改变的。那如果先把 dt*dt 计算好,之后直接乘就可以了,就不要再重新计算 1024 1024 1024 遍了,也就是这样写:

void func(float *__restrict a, float *__restrict b, float dt) {
    float dt2 = dt * dt;
    for (int i = 0; i < 1024; i++) {
        a[i]=a[i]+b[i]*dt2;
    }
}

所以再汇编中可以看到,它先是进行了一个 dt*dt 的运算 mulss %xmm0, %xmm0,后面就直接应用这个乘好的 %mm 了,就可以节省 1024 1024 1024 次乘法运算:
在这里插入图片描述
这种思想就是,把热的代码尽量移到冷的代码里。

6.5.1 挪到外面来:优化失败

如果把上面代码的括号去掉会怎样?

void func(float *__restrict a, float *__restrict b, float dt) {
    for (int i = 0; i < 1024; i++) {
        a[i]=a[i]+b[i]*dt*dt;
    }
}

只要去掉 (dt*dt) 的括号就会优化失败:因为乘法是左结合的,就相当于 (b[i]*dt)*dt。编译器识别不到不变量,从而优化失败。由于浮点运算不满足结合律,因为浮点不能精确表示实数。
在这里插入图片描述
因此,要么帮编译器打上括号帮助它识别,要么手动提取不变量到循环体外

6.6 调用不在另一个文件的函数:SIMD 优化失败

还有一个会导致 SIMD 优化失败的例子,就是调用外部函数。前面提到过,外部函数是优化失败的罪魁祸首,比如下面代码有一个 for 循环,本来这有这个是可以优化成功的,而在代码内调用了一个在其他文件的函数声明,编译器没办法矢量化,这里并不知道 other 是不是会改变 a a a,因为编译器并不知道 other 里干了啥哪怕 other 在定义它的文件里是个空函数,它也不敢优化。

void other();

float func(float *a) {
    float ret = 0;
    for (int i = 0; i < 1024; i++) {
        ret += a[i];
        other();
    }
    return ret;
}

可以看到汇编内全是 ss,说明它优化失败了,都是标量的。如果全是 ps 的话,才说明优化成功。

在这里插入图片描述

6.6.1 解决方案:放在同一个文件里

所以说,尽可能把 other 放在同一个文件里定义,这样编译器就能够看到。比如下面的代码:

void other() {
}

float func(float *a) {
    float ret = 0;
    for (int i = 0; i < 1024; i++) {
        ret += a[i];
        other();
    }
    return ret;
}

这里讲 other 放到了和 func 同一个 .cpp 文件里,这样编译器看得到 other 的函数体,就可以内联化该代码。(再次提醒,是不需要 inline 的)
在这里插入图片描述
所以说在热的 for 循环里,尽量不要调用外部函数,把它们移到同一个文件里,或者放在同文件声明为 static 函数。

6.7 循环中的下标

6.7.1 随机访问

像下面这种随机访问也会影响 SIMD,使优化失败:

void func(float *a, int *b) {
    for (int i = 0; i < 1024; i++) {
        a[b[i]] += 1;
    }
}

这是因为 b[i] 每次的值是不确定的,它可能是跳跃的在访问,这时候编译器就没有相应的指令能够生成优化。(理论上这种没什么好的优化方式)
在这里插入图片描述

6.7.2 循环中的下标:跳跃访问

下面这种情况也比较恶劣,就是 i*2,也就是说先访问 0 0 0、再访问 2 2 2 4 4 4

void func(float *a) {
    for (int i = 0; i < 1024; i++) {
        a[i * 2] += 1;
    }
}

就如同刚刚说的,这种跳跃的访问是不利于 SIMD 的,但是编译器就非常自然了,它用了一大堆 shufps 指令,还是成功的对 a a a 进行了一定的 SIMD 优化。(矢量化部分成功,但是非常艰难)
在这里插入图片描述

6.7.3 连续访问

所以最好的情况是什么呢?也就是直接 a[i],也就是顺序访问。这种情况下,编译器能够只用一个指令进行读写。

void func(float *a) {
    for (int i = 0; i < 1024; i++) {
        a[i] += 1;
    }
}

也就是说,不管是编译器还是 CPU,都喜欢顺序的连续访问。
在这里插入图片描述

6.12 编译器指令:循环展开

还有一个优化方式是循环展开。像下面这种循环体非常简单:

void func(float *a) {
    for (int i = 0; i < 1024; i++) {
        a[i] = 1;
    }
}

生成的代码一行是 movups,接着的一行是 addq,然后是 cmpq。就是说实际执行计算的指令只有一个,而进行比较和 i++ 的指令却用了两个,还有一个跳转指令。这样就导致大部分时间都花在跳转和 + 上了,而不是在实际的计算上。
在这里插入图片描述
对于 GCC 编译器,可以用 #pragma GCC unroll 4,表示把循环体展开为四个。添加后的代码为:

void func(float *a) {
#pragma GCC unroll 4
    for (int i = 0; i < 1024; i++) {
        a[i] = 1;
    }
}

比如下面是把四个循环体放在一个循环内,然后 i++ 改成 i+=4。注意,这里不建议手动这样写,会妨碍编译器的 SIMD 矢量化:

void func(float *a) {
    for (int i = 0; i < 1024; i += 4) {
        a[i + 0] = 1;
        a[i + 1] = 1;
        a[i + 2] = 1;
        a[i + 3] = 1;
    }
}

这样的话它就有 4 4 4movups 3 3 3 个其他指令,从而实际执行计算的时间就更多了,从而就可以避免反复比较的 overhead。(也就是本来是 1 1 1movups + 3 3 3 个其他指令,优化后变成 4 4 4movups + 3 3 3 个其他指令,其他指令使用的次数变少了)

在这里插入图片描述
注意,指令 #pragma GCC unroll 4 是 GCC 特定的,笔者目前没有找到什么通用的替代。

4 4 4 这个数字可以改,这里代表的是把 4 4 4 个循环体变成了一个,可以根据实际的情况决定要用多少。可以根据实际的情况决定要用多少,就是小的循环器一般 unroll 一下是划算的,但最好不要 unroll 大的循环体,否则会造成指令缓存的压力反而变慢

7. 结构体

这章是重点。

7.1 字节对齐

7.1.1 两个 float:对齐到 8 字节

比如我们有一个 MyVec,它里面有 x , y x, y x,y 两个值。

struct MyVec {
    float x;
    float y;
};

MyVec a[1024];

void func() {
    for (int i = 0; i < 1024; i++) {
        a[i].x *= a[i].y;
    }
}

func 内会使用到 x x x y y y 两个成员,可以看到非常复杂的代码:
在这里插入图片描述
里面有 shufps 什么的,但最终还是用到了 mulps,也就是它的确矢量化成功了。

7.1.2 三个 float:对齐到 12 字节

但是如果我们突然想把这个改成三维的了,就凭空加一个 z z z,然后我们的 func 甚至没有任何变化:

struct MyVec {
    float x;
    float y;
    float z;
};

MyVec a[1024];

void func() {
    for (int i = 0; i < 1024; i++) {
        a[i].x *= a[i].y;
    }
}

汇编出来的就全部变成了 ss,也就是矢量化失败了,这是什么原因呢?
在这里插入图片描述
这是因为刚才之所以能够矢量化,是因为它把两个 MyVec 给读出来,然后读到一个 SIMD 的 128 128 128 位寄存器里了。而在这里有 3 3 3 个 float 的话,它如果要读到一个 SIMD 里, 128 128 128 位就是 xyzx,剩余的一个 y y y 是没有办法读出来的。所以说这时候编译器就尴尬了,它没办法生成一个 SIMD 的读,所以说要怎么办呢?

7.1.3 添加一个辅助对齐的变量:对齐到 16 字节

如果想让它再变回高效,甚至可以像下面这样直接加一个 char padding[4],这是一个没有任何用处的空变量:

struct MyVec {
    float x;
    float y;
    float z;
    char padding[4];
};

MyVec a[1024];

void func() {
    for (int i = 0; i < 1024; i++) {
        a[i].x *= a[i].y;
    }
}

可以看到这时它又能够优化成功了。可以看到这里生成的 ps 的代码,是进行矢量化的:
在这里插入图片描述
可以发现,如果你的 struct 2 2 2 的整数幂大小,这时候是有利于 SIMD 优化的;而如果像刚才那样是三个一组,这时候它就没有办法对齐,从而导致 SIMD 优化失败。

结论:计算机喜欢 2 2 2 的整数幂,2, 4, 8, 16, 32, 64, 128…结构体大小若不是 2 的整数幂,往往会导致 SIMD 优化失败。

7.2 C++11 新语法:alignas

除了上面那种直接加一个从来用不到的变量的方法外,C++11 还新增了 alignas,也就是对齐到 16 16 16 字节。这样 MyVec 的起始地址,它会对齐到 16 16 16,而且它的大小也会变成 16 16 16 的整数倍。(这里的 16 16 16 代表的是字节而不是位)

struct alignas(16) MyVec {
    float x;
    float y;
    float z;
};

MyVec a[1024];

void func() {
    for (int i = 0; i < 1024; i++) {
        a[i].x *= a[i].y;
    }
}

可以看到这时候生成的汇编是一样的,但是又不需要手动加一个 char padding
在这里插入图片描述
那是不是所有结构体打上 alignas(16) 我的程序就会变快?
错了,有可能不仅不变快,反而还变慢!SIMD 和缓存行对齐只是性能优化的一个点,又不是全部。还要考虑结构体变大会导致内存带宽的占用,对缓存的占用等一系列连锁反应,比如导致它更容易变成 memory bound,从而有可能还会变慢,总之,要根据实际情况选择优化方案。实际测了以后如果的确变快了那才去加。

7.5 结构体的内存布局:AOS 与 SOA

  • AOS(Array of Struct)单个对象的属性紧挨着存-----xyzxyzxyzxyz
  • SOA(Struct of Array)属性分离存储在多个数组-----xxxxyyyyzzzz
  • AOS 必须对齐到 2 2 2 的幂才高效,SOA 就不需要;
  • AOS 符合直觉,不一定要存储在数组这种线性结构,而 SOA 可能无法保证多个数组大小一致;
  • SOA 不符合直觉,但通常是更高效的。

刚才说的结构是不是 x、y、z、padding?如下图所示,就是先 xyz,然后图上的黑框处空一个 padding,再 xyz,然后空一个 padding。这种结构类型叫 Array of Struct,它指的是 x、y、z 单个对象的属性紧紧凑在一起。

与之相对的就是 SOA,它会把 1 1 1 个对象拆成 3 3 3 个数组,就是一个对象的 x x x 在一个数组里,一个对象的 y y y 又在另一个数组里,这样的好处就是它不需要去考虑 padding 了。因为这样,如果只访问了 x x x y y y z z z 根本就没有用到,那这时候就相当于只用到了 x x x y y y;而如果用 AOS 的话, z z zpadding 都相当于没有用的那个属性。所以如果刚才只用 x x x y y y,那 CPU 看到就是在不停地跳着两格读,从而效率是最低的;而如果用 SOA,只访问 x x x y y y,没用 z z z,那 CPU 看到的是在访问两个数组,从而它能够用 SIMD 优化。所以说 SOA 虽然它好像不符合我们面向对象的思想,但它其实对 CPU 更友好。
在这里插入图片描述

7.5.1 AOS:紧凑存储多个属性

来看看怎么做的。首先是刚才说的比较烂的 AOS,这符合一般面向对象编程(OOP)的习惯,但常常不利于性能:

struct MyVec {
    float x;
    float y;
    float z;
};

MyVec a[1024];

void func() {
    for (int i = 0; i < 1024; i++) {
        a[i].x *= a[i].y;
    }
}

我们用上面代码生成的就是一个不能 SIMD 的代码:
在这里插入图片描述

7.5.2 SOA:分离存储多个属性

那要怎么变成 SOA 呢?可以将上面代码写成这样:

struct MyVec {
    float x[1024];
    float y[1024];
    float z[1024];
};

MyVec a;

void func() {
    for (int i = 0; i < 1024; i++) {
        a.x[i] *= a.y[i];
    }
}

之前是 a[1024],现在 a[1024] 不要了,把它移动到 MyVec 里,就是 x x x 1024 1024 1024 个, y y y 1024 1024 1024。访问的时候本来是 a[i].x,现在变成 a.x[i],就是括号换了下位置。它们分离开来的存储,可以看到下面就用到了 ps,从而矢量化是成功的:

在这里插入图片描述
这样不符合 面向对象编程(OOP) 的习惯,但常常有利于性能。又称为 面向数据编程(DOP, date-oriented programming),这种在高性能计算中比较常见。

7.5.3 AOSOA:中间方案

但是可以发现,SOA 有时不利于优化,所以说推出了一种中间解决方案,就 4 4 4 个对象打包乘以 SOA,再用一个 n / 4 n/4 n/4 大小的数组打包成一个 AOS,这样既有 AOS 的直观,又有 SOA 的高效:

struct MyVec {
    float x[4];
    float y[4];
    float z[4];
};

MyVec a[1024 / 4];

void func() {
    for (int i = 0; i < 1024 / 4; i++) {
        for (int j = 0; j < 4; j++) {
            a[i].x[j] *= a[i].y[j];
        }
    }
}

这里其实就相当于把 4 4 4 个变成一个 m128 了。但是这种其实也不常用,最常用的还是上面的 SOA。在这里插入图片描述

  • AOSOA 优点:
    SOA 便于 SIMD 优化;AOS 便于存储在传统容器;AOSOA 两者得兼。
  • AOSOA 缺点:
    需要两层 for 循环,不利于随机访问;需要数组大小是 4 4 4 的整数倍,不过可以用边界特判法解决。在代码中可以看到,前面都只需要一个 for 循环,而这里却需要两个,就很麻烦了,而且这对随机访问也非常不友好。

7.9 测试一下加速了多少倍?

这里做个总结测试,比如刚才的程序,有 x x x y y y z z z 3 3 3 个属性:

struct Point {
    float x;
    float y;
    float z;
};

Point ps[N];

void compute() {
    for (int i = 0; i < N; i++) {
        ps[i].x = ps[i].x + ps[i].y + ps[i].z;
    }
}

第一步将它优化成 SOA,也就是把 N 移上去。再加上每个编译器特有的 unroll 语句,如 GCC 是 #progma GCC unroll 32,MSC 也就是微软的比较霸气,直接就是 #progma unroll 32

struct Point {
    float x[N];
    float y[N];
    float z[N];
};

Point ps;

void compute() {
    for (int i = 0; i < N; i++) {
        ps.x[i] = ps.x[i] + ps.y[i] + ps.z[i];
    }
}

下面是测试结果:
在这里插入图片描述

  • 811620 n s 811620ns 811620ns:未经优化的耗时 ;
  • 989485 n s 989485ns 989485ns:用了 16 16 16 字节对齐以后的效率,它耗时反而更长了;
  • 176316 n s 176316ns 176316ns:如果把原来代码一点不变,只是加上一个 parallel for 的效果;
  • 218063 n s 218063ns 218063ns:使用 SOA 优化的耗时;
  • 187422 n s 187422ns 187422ns:SOA+grpgrma omp simd 忽略指针别名后的耗时;
  • 186657 n s 186657ns 186657ns:使用 size_t 作为下标以后的耗时;
  • 158938 n s 158938ns 158938ns:用 size_t 以后再去 unroll,可以看到它甚至比前面直接加一个 parallel 还要快;
  • 59488 n s 59488ns 59488ns:如果优化以后再上一个并行,这样就比谁都快了;
  • 541613 n s 541613ns 541613ns:使用 AOSOA 的优化,它虽然快一点,但还是不如 SOA;

下面是我电脑上的测试效果,差的有点多,首先是对齐后的效率是会更快的,再者 AOSOA 甚至比 SOA 还快。
在这里插入图片描述

图中可以看到,一条绿线是无脑并行的,而粉线是我们深度优化,可以看到再深度优化之后,再低数据量的时候甚至超过了无脑并行,因为并行需要创建很多线程,这也需要时间,所以对于小数据甚至超过了它。
在这里插入图片描述
这里的结论是:SOA 是针对这个案例最高效的数据排布格式。当然不一定所有的案例都可以用 SOA 就能优化,需要看情况。

8. STL 容器

前面都是基于 C 语言的指针,有些老了。C++ 喜欢 RAII 思想,那就得用 STL 的 vector 容器。

8.1 std::vector:也有指针别名问题

首先是一个经典案例,就是第一章说到的指针别名问题,它不能保证下面代码内的 c c c b b b 是同一个 vector,从而它得生成两遍:

#include <vector>

void func(std::vector<int> &a,
          std::vector<int> &b,
          std::vector<int> &c) {
    c[0] = a[0];
    c[0] = b[0];
}

汇编结果如下:
在这里插入图片描述

8.1.1 __restrict:能否用于 std::vector?

那我想把这个引用声明 __restrict 行不行:

#include <vector>

void func(std::vector<int> &__restrict a,
          std::vector<int> &__restrict b,
          std::vector<int> &__restrict c) {
    c[0] = a[0];
    c[0] = b[0];
}

可以从下面汇编看到,没有任何效果。因为这里只是把 vector 本身做了一个 __restrict ,而这个 vector 里的指针,它没有上 __restrict,这时候编译器还是不知道 c 和 a 它是不是同一地址的:
在这里插入图片描述

8.1.2 解决方案:pragma omp simd 或 pragma GCC ivdep

所以结论就是:要么用 #pragma omp simd 让它忽视可能存在的指针别名,或者用 #pragma gcc ivdep 也可以:
#pragma omp simd:

void func(std::vector<int> &a,
          std::vector<int> &b) {
#pragma omp simd
    for (int i = 0; i < 1024; i++) {
        a[i] = b[i] + 1;
    }
}

#pragma gcc ivdep:

void func(std::vector<int> &a,
          std::vector<int> &b) {
#pragma GCC ivdep
    for (int i = 0; i < 1024; i++) {
        a[i] = b[i] + 1;
    }
}

C/C++ 缺点:指针自由度过高,允许多个 immutable reference 指向同一个对象。

所以我们在适当的时候要加上一些提示来帮助编译器优化。

8.4 std::vector:也能实现 SOA!

先前不是说 SOA 比较好嘛,但是如果使用了 vector 该怎么办?可以像下面这样把 vector 移上去:

  • 优化前(AOS):

    struct MyVec {
        float x;
        float y;
        float z;
    };
    std::vector<MyVec> a;
    void func() {
        for (std::size_t i = 0; i < a.size(); i++) {
            a[i].x *= a[i].y;
        }
    }
    
  • 优化后(SOA):

    struct MyVec {
    	std::vector<float> x;
    	std::vector<float> y;
    	std::vector<float> z;
    };
    MyVec a;
    void func() {
    	for (std::size_t i = 0; i < a.x.size(); i++) {
        	a.x[i] *= a.y[i];
    	}
    }
    

    不过需要保证, x x x y y y z z z 当中一个被 push 的时候, 3 3 3 个都被 push。否则的话下面的 for 循环会出错。

9. 数学运算

9.1 数学优化:除法变乘法

下面示例接受 a a a 参数,然后变成 a/2

float func(float a) {
    return a / 2;
}

返回 a/2,生成的汇编如下所示:
在这里插入图片描述
这里 mulss 是一个乘,.long 1056964608 是浮点数 0.5 0.5 0.5 的意思,也就是说它把一个除法优化成了更快的乘法。

9.2 编译器放弃的优化:分离公共除数

虽然乘法更快,但这个示例中它却没有把除法提取成一个倒数出来:

void func(float *a, float b) {
    for (int i = 0; i < 1024; i++) {
        a[i] /= b;
    }
}

为什么它不能优化呢?因为它害怕 b b b 0 0 0 的情况,b 是 0 0 0 的话公式 a[i]/=b 会变成 inf 或者说它会报一个错。
在这里插入图片描述

9.2.1 解决方案1:手动优化

这时候编译器就不敢优化,所以我们可以手动进行优化,就把 1/b 预先计算出来, 然后在前面公式内除法改成乘法:

void func(float *a, float b) {
    float inv_b = 1 / b;
    for (int i = 0; i < 1024; i++) {
        a[i] *= inv_b;
    }
}

这样就能够把 1024 1024 1024 次除法变成 1024 1024 1024 次乘法加一次除法。这样其实是更快的。
在这里插入图片描述
也就是说,乘法比除法更快。提前计算好 b 的倒数避免重复求除法。

9.2.2 解决方案2:-ffast-math

既然前面放弃优化了,这里还有一个解决方案,就是加选项 gcc -ffast-math -O3

void func(float *a, float b) {
    for (int i = 0; i < 1024; i++) {
        a[i] /= b;
    }
}

-ffast-math 选项让 GCC 更大胆地尝试浮点运算的优化,有时能带来 2 2 2 倍左右的提升。作为代价,它对 NaN 和无穷大的处理,可能会和 IEEE 标准规定的不一致。如果能保证程序中永远不会出现 NaN 和无穷大,那么可以放心打开 -ffast-math。(测试发现有时候加了 -ffast-math 是快,但是算出来的值差很多

在这里插入图片描述

9.3 数学函数请加 std:: 前缀!

还有一个问题在于,数学函数 sqrt 只接受 double 的,如果使用 float 去调用 sqrt,它也会返回一个 double,而且计算精度也是基于 double 的,计算速度会慢一些。而如果想要调用 float 的开根号,其实应该用 sqrtf。所以说最好的解决方案应该是 C++ 的 sqrt 函数,它重载了doublefloat 两个版本,就能够自动根据拿什么调用它就自动返回什么类型:

#include <cmath>
#include <type_traits>

void func(float *a) {
    for (int i = 0; i < 1024; i++) {
        auto x = sqrt(a[i]);      // wrong and bad
        static_assert(std::is_same_v<decltype(x), double>);
        auto y = sqrtf(a[i]);     // correct but bad
        static_assert(std::is_same_v<decltype(y), float>);
        auto z = std::sqrt(a[i]); // correct and good
        static_assert(std::is_same_v<decltype(z), float>);
    }
}

也就是推荐不要再用 C 语言全局的 sqrt,多用用 std 里面的。

更坑的是 abs,如果使用全局命名空间的 abs,会发现 abs(1.4f) = 1,这是因为 abs 只接受 int,所以 1.4 1.4 1.4 被隐式的转换成了 int。要用 fabs 才能调用 double 的;fabsf 调用 float 的。而 std::abs 它重载了 intfloatdouble ,就不会出现上面奇怪的现象了。

9.4 -ffast-math 的另一优点

-ffast-math 还有一个优点,比如下面是对 a 的每一个数求根号:

void func(float *a) {
    for (int i = 0; i < 1024; i++) {
        a[i] = std::sqrt(a[i]);
    }
}

9.4.1 开启前:sqrt 矢量化失败

开启前可以看到,这里是没有矢量化的:
在这里插入图片描述
为什么这里没有矢量化呢?因为它生怕 a[i] 是一个负数,那么计算到 a[i] 就会出错,然后停止。但是如果矢量化了,那就变成 a[i] 出错时它已经写入了,这里不能保证写入的顺序 ,所以它就不敢优化。

9.4.2 开启后:sqrt 矢量化成功

在开启 -ffast-math 以后,就告诉它不要怕,我能保证 a 绝对不会出错,它就能优化成下面这个:
在这里插入图片描述

9.5 嵌套循环:直接累加,有指针别名问题

下面的嵌套循环在计算卷积或者矩阵乘法都会用到:

void func(float *a, float *b, float *c) {
    for (int i = 0; i < 1024; i++) {
        for (int j = 0; j < 1024; j++) {
            c[i] += a[i] * b[j];
        }
    }
}

但直接这样写是有问题的,可以看到这里 a a a c c c 是有可能重合的,甚至 b b b c c c 是重合的。编译器担心 c c c a a a 可能会指向同一个地址,而连续判断三个指针是否有重合又过于复杂,所以放弃了矢量化。
在这里插入图片描述
从这里的汇编可以看到,它与前面的章节不同,前面章节还生成了两个版本,而这里只生成了一个版本,即没有生成 SIMD 版本。这就是因为当编译器碰到问题复杂化时,它就直接放弃优化了。所以要么像前面的方式添加 __restrict,要么使用下面的方式:

9.5.1 解决方案1:先读到局部变量,累加完毕后,再写入

可以先把 c c c 读到 tmp,然后对 tmp 进行追加,读完之后再写回到 c c c

void func(float *a, float *b, float *c) {
    for (int i = 0; i < 1024; i++) {
        float tmp = c[i];
        for (int j = 0; j < 1024; j++) {
            tmp += a[i] * b[j];
        }
        c[i] = tmp;
    }
}

这样对寄存器更加友好,这时候编译器认为不会存在指针别名的问题,矢量化成功:
在这里插入图片描述

9.5.2 解决方案2:先累加到初始为 0 的局部变量,再累加到 c

另一种更好的方法是将 tmp 初始化为 0 0 0。这与上面方法是有区别的:上面一种加法是在 c[i] 已经有值的情况下去加,这与如果 c[i] 是一个比较大的值,比如 1000 1000 1000,而 a a a 是一个比较小的值,比如 0.01 0.01 0.01。浮点数的构造有一个指数有一个底数,而指数相差太大的话,它们两个之间相加的精度就会受损。就比如 1000 + 0.000001 1000+0.000001 1000+0.000001,可能还是 1000 1000 1000,就导致加出来的结果可能不对,所以我们最好是把 tmp 初始化为 0 0 0,然后这样再去 c[i] += tmp 的话就会更精确。因为这样的话 tmp 就足够大,从而不会出现指数位差太远,导致精度误差的问题。

void func(float *a, float *b, float *c) {
    for (int i = 0; i < 1024; i++) {
        float tmp = 0;
        for (int j = 0; j < 1024; j++) {
            tmp += a[i] * b[j];
        }
        c[i] += tmp;
    }
}

这里是精度上的提升,而不是性能上的提升。
在这里插入图片描述

10. 优化手法总结

  1. 函数尽量写在同一个文件内;
  2. 避免在 for 循环内调用外部函数;
  3. 非 const 指针加上 __restrict 修饰;
  4. 试着用 SOA 取代 AOS;
  5. 对齐到 16 或 64 字节( 16 16 16 是 SIMD 宽度, 64 64 64 是缓存行宽度);
  6. 简单的代码,不要复杂化(把代码写的简单点让编译器看懂;太复杂的看不懂就没法优化);
  7. 试试看 #pragma omp simd;
  8. 循环中不变的常量挪到外面来;
  9. 对小循环体用 #pragma unroll;
  10. -ffast-math 和 -march=native(对于 GCC 编译器可以加上这两个,让它快 2 2 2 倍甚至 4 4 4 倍都有可能)。

在 CMake 中开启的几个开关

  • CMake 中开启 -O3
    默认的 set(CMAKE_BUILD_TYPE Release)Debug,它会生成一些调试,但是它的优化就会被关掉。而 Release 模式,它就相当于 -O3

  • CMake 中开启 -fopenmp
    因为 fopenmp 是只针对 GCC 的开启方式,微软的可能不一样,所以可以使用这种通用的方式来开启:

    find_package(OpenMP REQUIRED)
    target_link_libraries(testbench PUBLIC OpenMP::OpenMP_CXX)
    

    它是做成了一个库的形式。

  • CMake 中开启 -ffast-math 和 -march=native

    target_compile_options(testbench PUBLIC -ffast-math -march=native)
    

    这个好像没办法变成一个跨平台的,只能针对 GCC 开启这两个。

  • 6
    点赞
  • 24
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

泠山

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值