太长不看版
template <auto Start, auto End>
constexpr void constexpr_for(auto &&f) {
if constexpr (Start < End) {
f.template operator()<Start>();
constexpr_for<Start + 1, End>(std::forward<decltype(f)>(f));
}
}
结论先行：
- 使用
constexpr_for即可做到显式 ILP。 - 使用 SIMD intrinsic 编程，即使是 -O3 优化也不会做显式 ILP。
- 在编译器优化质量较低时，显式 ILP 可以掩盖性能问题。
- 注意 std::ranges 的抽象税。
下面会用 SIMD 求和的示例程序来说明。限于篇幅，本文不会解释 ILP 的基本概念。
SIMD sum
#include <x86intrin.h>
#include <algorithm>
#include <ranges>
// 使用 C++23 的 range_adaptor_closure 进行快速定制适配器
// 从而方便构造 SIMD 的批量视图（非完整块会丢弃）
template <size_t Lane>
struct simdify_t: std::ranges::range_adaptor_closure<simdify_t<Lane>> {
constexpr auto operator()(auto &&view) const noexcept {
auto size = std::ranges::size(view);
return std::forward<decltype(view)>(view)
| std::views::stride(Lane) // 每 Lane 个元素取一个
| std::views::take(size / Lane); // 总共只取 size / Lane 个元素
}
};
template <auto i>
constexpr simdify_t<i> simdify;
int sum_avx2(std::ranges::range auto &&rng) {
constexpr auto lane = sizeof(__m256i) / sizeof(int);
__m256i partial_sum {};
auto simd_view = rng
| simdify<lane>;
for(auto &&simd_v : simd_view) {
auto addr = &simd_v;
auto loadu = _mm256_loadu_si256((__m256i*)addr);
partial_sum = _mm256_add_epi32(partial_sum, loadu);
}
auto scalar_view = rng
| std::views::drop(lane * std::ranges::size(simd_view));
int sum = std::ranges::fold_left(scalar_view, 0, std::plus());
int temp[lane];
// Clang 编译器会优化为寄存器内归约操作
_mm256_storeu_si256((__m256i*)std::ranges::data(temp), partial_sum);
sum = std::ranges::fold_left(temp, sum, std::plus());
return sum;
}
我们看一个 AVX2 版本的 SIMD 求和实现。非常直接，示例程序将输入的整个 range 拆分为 SIMD 视图和 scalar 视图并分别计算，计算完成后合并求和结果。
SIMD+ILP sum
#include <x86intrin.h>
#include <algorithm>
#include <ranges>
// 显式 ILP 的工具
template <auto Start, auto End>
constexpr void constexpr_for(auto &&f) {
if constexpr (Start < End) {
f.template operator()<Start>();
constexpr_for<Start + 1, End>(std::forward<decltype(f)>(f));
}
}
template <size_t Lane>
struct simdify_t: std::ranges::range_adaptor_closure<simdify_t<Lane>> {
constexpr auto operator()(auto &&view) const noexcept {
auto size = std::ranges::size(view);
return std::forward<decltype(view)>(view)
| std::views::stride(Lane)
| std::views::take(size / Lane);
}
};
template <auto i>
constexpr simdify_t<i> simdify;
template <size_t ILP = 4>
int sum_avx2_ilp(std::ranges::range auto &&rng) {
constexpr auto lane = sizeof(__m256i) / sizeof(int);
constexpr auto bulk = lane * ILP;
__m256i partial_sum[ILP] {};
auto process_simd = [&]<size_t Width>(auto simd_view) {
for(auto &&simd_v : simd_view) {
// ILP 展开，最高一次并行处理 lane * ILP * 4 字节的数据
// 并且同时支持基本的 SIMD 操作（Width=1）
constexpr_for<0, Width>([&, addr = &simd_v]<size_t Index> {
auto &partial = partial_sum[Index];
auto loadu = _mm256_loadu_si256(
(__m256i*)(addr + Index * lane));
partial = _mm256_add_epi32(partial, loadu);
});
}
};
// SIMD+ILP
auto bulk_simd_view = rng
| simdify<lane>
| simdify<ILP>; // 这里就知道定制 simdify 的妙处了
process_simd.template operator()<ILP>(bulk_simd_view);
// SIMD
auto single_simd_view = rng
| std::views::drop(bulk * std::ranges::size(bulk_simd_view))
| simdify<lane>;
process_simd.template operator()<1>(single_simd_view);
auto scalar_view = rng
| std::views::drop(bulk * std::ranges::size(bulk_simd_view))
| std::views::drop(lane * std::ranges::size(single_simd_view));
int sum = std::ranges::fold_left(scalar_view, 0, std::plus());
constexpr_for<1, ILP>([&]<size_t Index> {
partial_sum[0] = _mm256_add_epi32(partial_sum[0], partial_sum[Index]);
});
int temp[lane];
_mm256_storeu_si256((__m256i*)std::ranges::data(temp), partial_sum[0]);
sum = std::ranges::fold_left(temp, sum, std::plus());
return sum;
}
显式 ILP 的实现如上。在前面例子的基础上扩展，将输入 range 拆分为三个视图：bulk SIMD / SIMD / scalar，同样是只有无法完整覆盖的块才会被下一级的视图处理。
也就是说 ILP 版本的思路是扩展到流水线的宽度，而不是单个向量寄存器的宽度：处理器端口这么多，大于一总比等于一好。
动机和跑分
可是为什么需要特意写 constexpr_for，编译器不会优化吗？
与其理论，不如直接跑分。我们使用 Clang-20 编译器和 Google Benchmark，开启 -O3 优化和 -march=znver3 目标架构（面向我的 Zen 3 主机）。
Run on (16 X 3193.93 MHz CPU s)
CPU Caches:
L1 Data 32 KiB (x8)
L1 Instruction 32 KiB (x8)
L2 Unified 512 KiB (x8)
L3 Unified 16384 KiB (x1)
Load Average: 0.47, 0.39, 0.25
------------------------------------------------------------------------------------
Benchmark Time CPU Iterations UserCounters...
------------------------------------------------------------------------------------
BM_sum_avx2/35 3.43 ns 3.36 ns 209909024 bytes_per_second=38.7807Gi/s
BM_sum_avx2/350 53.9 ns 52.8 ns 12902065 bytes_per_second=24.6894Gi/s
BM_sum_avx2/3502 623 ns 610 ns 1142980 bytes_per_second=21.3926Gi/s
BM_sum_avx2/35023 6298 ns 6165 ns 111805 bytes_per_second=21.1622Gi/s
BM_sum_avx2/350234 62364 ns 61045 ns 11117 bytes_per_second=21.3731Gi/s
BM_sum_avx2_ilp<4>/35 1.33 ns 1.30 ns 532863017 bytes_per_second=100.414Gi/s
BM_sum_avx2_ilp<4>/350 23.3 ns 22.8 ns 30071771 bytes_per_second=57.157Gi/s
BM_sum_avx2_ilp<4>/3502 167 ns 164 ns 4274452 bytes_per_second=79.6004Gi/s
BM_sum_avx2_ilp<4>/35023 1606 ns 1572 ns 445802 bytes_per_second=83.0013Gi/s
BM_sum_avx2_ilp<4>/350234 18399 ns 18002 ns 42099 bytes_per_second=72.4749Gi/s
总之，我们确实得到了数倍的并行性能提升，也验证了 ILP 策略的有效性。
汇编结果
avx::test_medium(std::array<int, 3502ul>):
// ...
vpaddd (%rsi), %ymm0, %ymm0
// ...
avx_ilp::test_medium(std::array<int, 3502ul>):
// ...
vpaddd (%rcx,%rdx), %ymm1, %ymm1
vpaddd 32(%rcx,%rdx), %ymm2, %ymm2
vpaddd 64(%rcx,%rdx), %ymm0, %ymm0
vpaddd 96(%rcx,%rdx), %ymm3, %ymm3
// ...
我们可以在 godbolt 快速查验编译器的生成质量。一个很明确的结论是：SIMD intrinsic 确实不会在指令并行这方面有优化，至少现在（2025 年）不行。
NOTES:
- intrinsic 有一定程度的优化，但是下面会指出它相对常规指令的优化会怂很多。
- 能不能用
for替代constexpr_for？结论是能但是不可靠，Clang 在 -O2 级别就没法优化到位。
也就是说，程序中的 __m256i partial_sum 是只能作为一个 ymm 寄存器去使用。而 constexpr_for 显式地以 ILP 的形式去展开，所以 ILP=4 时会明确使用 4 路并行的 ymm 寄存器。
基线差距
// 在 Google benchmark 注册一个新的测试
constexpr auto call_std_fold_left = [](auto&& rng) {
return std::ranges::fold_left<4>(std::forward<decltype(rng)>(rng), 0, std::plus());
};
我们再以 std::ranges::fold_left 作为基线对比，它的实现非常朴素：libstdc++/ranges_algo.h。
--------------------------------------------------------------------------------------
Benchmark Time CPU Iterations UserCounters...
--------------------------------------------------------------------------------------
BM_sum_avx2/35 3.24 ns 3.21 ns 223461745 bytes_per_second=40.6818Gi/s
BM_sum_avx2/350 54.0 ns 53.4 ns 12524425 bytes_per_second=24.4263Gi/s
BM_sum_avx2/3502 612 ns 612 ns 1136500 bytes_per_second=21.3104Gi/s
BM_sum_avx2/35023 6172 ns 6164 ns 112560 bytes_per_second=21.1683Gi/s
BM_sum_avx2/350234 61234 ns 61237 ns 11033 bytes_per_second=21.306Gi/s
BM_sum_avx2_ilp<4>/35 0.959 ns 0.961 ns 735636420 bytes_per_second=135.736Gi/s
BM_sum_avx2_ilp<4>/350 15.4 ns 15.4 ns 44097072 bytes_per_second=84.4887Gi/s
BM_sum_avx2_ilp<4>/3502 142 ns 143 ns 4898694 bytes_per_second=91.4652Gi/s
BM_sum_avx2_ilp<4>/35023 1622 ns 1620 ns 443559 bytes_per_second=80.522Gi/s
BM_sum_avx2_ilp<4>/350234 20027 ns 20051 ns 40449 bytes_per_second=65.0717Gi/s
BM_sum_std_fold_left/35 1.00 ns 1.00 ns 677987410 bytes_per_second=130.061Gi/s
BM_sum_std_fold_left/350 8.67 ns 8.68 ns 78784950 bytes_per_second=150.151Gi/s
BM_sum_std_fold_left/3502 91.2 ns 91.3 ns 7553542 bytes_per_second=142.956Gi/s
BM_sum_std_fold_left/35023 1034 ns 1035 ns 657369 bytes_per_second=126.026Gi/s
BM_sum_std_fold_left/350234 12931 ns 13049 ns 54177 bytes_per_second=99.9844Gi/s
可以看出，fold_left 优化明显更强，我们有必要再看一眼汇编去确认它是怎么做的。
rng::test_medium(std::array<int, 3502ul>):
// ...
.LBB10_1:
vpaddd -992(%rcx,%rax,4), %ymm0, %ymm0
vpaddd -960(%rcx,%rax,4), %ymm1, %ymm1
vpaddd -928(%rcx,%rax,4), %ymm2, %ymm2
vpaddd -896(%rcx,%rax,4), %ymm3, %ymm3
vpaddd -864(%rcx,%rax,4), %ymm0, %ymm0
vpaddd -832(%rcx,%rax,4), %ymm1, %ymm1
vpaddd -800(%rcx,%rax,4), %ymm2, %ymm2
vpaddd -768(%rcx,%rax,4), %ymm3, %ymm3
vpaddd -736(%rcx,%rax,4), %ymm0, %ymm0
vpaddd -704(%rcx,%rax,4), %ymm1, %ymm1
vpaddd -672(%rcx,%rax,4), %ymm2, %ymm2
vpaddd -640(%rcx,%rax,4), %ymm3, %ymm3
vpaddd -608(%rcx,%rax,4), %ymm0, %ymm0
vpaddd -544(%rcx,%rax,4), %ymm2, %ymm4
vpaddd -576(%rcx,%rax,4), %ymm1, %ymm1
vpaddd -512(%rcx,%rax,4), %ymm3, %ymm5
vpaddd -480(%rcx,%rax,4), %ymm0, %ymm2
vpaddd -448(%rcx,%rax,4), %ymm1, %ymm3
vpaddd -416(%rcx,%rax,4), %ymm4, %ymm0
vpaddd -384(%rcx,%rax,4), %ymm5, %ymm1
cmpq $3577, %rax
je .LBB10_3
vpaddd -352(%rcx,%rax,4), %ymm2, %ymm2
vpaddd -320(%rcx,%rax,4), %ymm3, %ymm3
vpaddd -288(%rcx,%rax,4), %ymm0, %ymm0
vpaddd -256(%rcx,%rax,4), %ymm1, %ymm1
vpaddd -224(%rcx,%rax,4), %ymm2, %ymm2
vpaddd -160(%rcx,%rax,4), %ymm0, %ymm4
vpaddd -192(%rcx,%rax,4), %ymm3, %ymm3
vpaddd -128(%rcx,%rax,4), %ymm1, %ymm5
vpaddd -96(%rcx,%rax,4), %ymm2, %ymm0
vpaddd -64(%rcx,%rax,4), %ymm3, %ymm1
vpaddd -32(%rcx,%rax,4), %ymm4, %ymm2
vpaddd (%rcx,%rax,4), %ymm5, %ymm3
addq $256, %rax
jmp .LBB10_1
// ...
其实 fold_left 不仅会并行使用多路寄存器，而且还进一步使用（非常暴力的）循环展开。
进一步优化
现在以 fold_left 为目标改进原有的程序。
auto process_simd = [&]<size_t Width>(auto simd_view) {
#pragma clang loop unroll_count(4)
for(auto &&simd_v : simd_view) {
// ...
}
}
首先能想到的是编译器特定的 unroll 选项，结论是对 ranges 起负作用。感兴趣可以自行看下汇编，编译器会生成一堆非常混乱的指令。
// ILP 版本类似，但是问题更加严重
avx::test_medium(std::array<int, 3502ul>):
// ...
.LBB1_1:
vpaddd (%rsi), %ymm0, %ymm0
movq %rax, %rdi
subq %rsi, %rdi
leaq 32(%rsi), %r8
sarq $2, %rdi
addq $-9, %rdi
cmpq $-8, %rdi
cmovaeq %rax, %r8
cmpq %rsi, %rax
cmovneq %r8, %rsi
cmpq %rdx, %rsi
jne .LBB1_1
进一步调查 AVX/ILP 版本的汇编，发现另一个问题是：它总是在计算莫名其妙的边界，一个合理推测是 ranges 本身的编译器优化不到位。也就是说哪怕是 ILP 版本也是有大量无用指令的。一个简单的解决方案是使用最朴素的 for 循环，不交抽象税。
template <size_t Lane>
struct simdify_t: std::ranges::range_adaptor_closure<simdify_t<Lane>> {
constexpr auto operator()(auto &&r) const noexcept {
auto v = std::forward<decltype(r)>(r) | std::views::all;
auto new_size = std::ranges::size(v) / Lane;
return std::views::transform(
std::views::iota(size_t{0}, new_size),
[v](size_t i) -> decltype(auto) {
return v[i * Lane];
}
);
}
};
template <auto i>
constexpr simdify_t<i> simdify;
当然个人原因我不打算放弃 ranges，重写一份让编译器能更加优化友好的版本（DoD 风格？）。这样我们程序的核心算法完全不需要修改。
avx::test_medium(std::array<int, 3502ul>):
movl $704, %eax
leaq 8(%rsp), %rcx
vpxor %xmm0, %xmm0, %xmm0
.LBB1_1:
vpaddd -704(%rcx,%rax), %ymm0, %ymm0
vpaddd -672(%rcx,%rax), %ymm0, %ymm0
vpaddd -640(%rcx,%rax), %ymm0, %ymm0
vpaddd -608(%rcx,%rax), %ymm0, %ymm0
vpaddd -576(%rcx,%rax), %ymm0, %ymm0
vpaddd -544(%rcx,%rax), %ymm0, %ymm0
vpaddd -512(%rcx,%rax), %ymm0, %ymm0
vpaddd -480(%rcx,%rax), %ymm0, %ymm0
vpaddd -448(%rcx,%rax), %ymm0, %ymm0
vpaddd -416(%rcx,%rax), %ymm0, %ymm0
vpaddd -384(%rcx,%rax), %ymm0, %ymm0
vpaddd -352(%rcx,%rax), %ymm0, %ymm0
vpaddd -320(%rcx,%rax), %ymm0, %ymm0
vpaddd -288(%rcx,%rax), %ymm0, %ymm0
vpaddd -256(%rcx,%rax), %ymm0, %ymm0
vpaddd -224(%rcx,%rax), %ymm0, %ymm0
vpaddd -192(%rcx,%rax), %ymm0, %ymm0
vpaddd -160(%rcx,%rax), %ymm0, %ymm0
vpaddd -128(%rcx,%rax), %ymm0, %ymm0
vpaddd -96(%rcx,%rax), %ymm0, %ymm0
vpaddd -64(%rcx,%rax), %ymm0, %ymm0
vpaddd -32(%rcx,%rax), %ymm0, %ymm0
vpaddd (%rcx,%rax), %ymm0, %ymm0
addq $736, %rax
cmpq $14688, %rax
jne .LBB1_1
// ...
avx_ilp::test_medium(std::array<int, 3502ul>):
movl $992, %eax
leaq 8(%rsp), %rcx
vpxor %xmm0, %xmm0, %xmm0
vpxor %xmm1, %xmm1, %xmm1
vpxor %xmm2, %xmm2, %xmm2
vpxor %xmm3, %xmm3, %xmm3
.LBB5_1:
vpaddd -992(%rcx,%rax), %ymm2, %ymm2
vpaddd -960(%rcx,%rax), %ymm1, %ymm1
vpaddd -928(%rcx,%rax), %ymm0, %ymm0
vpaddd -896(%rcx,%rax), %ymm3, %ymm3
vpaddd -864(%rcx,%rax), %ymm2, %ymm2
vpaddd -832(%rcx,%rax), %ymm1, %ymm1
vpaddd -800(%rcx,%rax), %ymm0, %ymm0
vpaddd -768(%rcx,%rax), %ymm3, %ymm3
vpaddd -736(%rcx,%rax), %ymm2, %ymm2
vpaddd -704(%rcx,%rax), %ymm1, %ymm1
vpaddd -672(%rcx,%rax), %ymm0, %ymm0
vpaddd -640(%rcx,%rax), %ymm3, %ymm3
vpaddd -608(%rcx,%rax), %ymm2, %ymm2
vpaddd -544(%rcx,%rax), %ymm0, %ymm0
vpaddd -576(%rcx,%rax), %ymm1, %ymm1
vpaddd -512(%rcx,%rax), %ymm3, %ymm4
vpaddd -448(%rcx,%rax), %ymm1, %ymm3
vpaddd -480(%rcx,%rax), %ymm2, %ymm2
vpaddd -416(%rcx,%rax), %ymm0, %ymm0
vpaddd -384(%rcx,%rax), %ymm4, %ymm1
cmpq $14304, %rax
je .LBB5_3
vpaddd -352(%rcx,%rax), %ymm2, %ymm2
vpaddd -320(%rcx,%rax), %ymm3, %ymm3
vpaddd -288(%rcx,%rax), %ymm0, %ymm0
vpaddd -256(%rcx,%rax), %ymm1, %ymm1
vpaddd -224(%rcx,%rax), %ymm2, %ymm2
vpaddd -160(%rcx,%rax), %ymm0, %ymm0
vpaddd -192(%rcx,%rax), %ymm3, %ymm3
vpaddd -128(%rcx,%rax), %ymm1, %ymm4
vpaddd -64(%rcx,%rax), %ymm3, %ymm1
vpaddd -96(%rcx,%rax), %ymm2, %ymm2
vpaddd -32(%rcx,%rax), %ymm0, %ymm0
vpaddd (%rcx,%rax), %ymm4, %ymm3
addq $1024, %rax
jmp .LBB5_1
// ...
修改后可以确认循环展开和多路寄存器都是有的，而且不再有诡异的边界计算。一次性解决了两个问题。简单对比，可以认为 ILP 版本和 fold_left 基本一致，AVX 版本去掉了多路 YMM 并行（只依靠硬件调度）。
--------------------------------------------------------------------------------------
Benchmark Time CPU Iterations UserCounters...
--------------------------------------------------------------------------------------
BM_sum_avx2/35 0.996 ns 0.997 ns 708214430 bytes_per_second=130.729Gi/s
BM_sum_avx2/350 9.88 ns 9.91 ns 69465475 bytes_per_second=131.578Gi/s
BM_sum_avx2/3502 96.7 ns 96.8 ns 7250274 bytes_per_second=134.724Gi/s
BM_sum_avx2/35023 1048 ns 1046 ns 646094 bytes_per_second=124.79Gi/s
BM_sum_avx2/350234 15356 ns 15354 ns 47235 bytes_per_second=84.9777Gi/s
BM_sum_avx2_ilp<4>/35 0.959 ns 0.957 ns 696597539 bytes_per_second=136.235Gi/s
BM_sum_avx2_ilp<4>/350 9.75 ns 9.74 ns 70929427 bytes_per_second=133.853Gi/s
BM_sum_avx2_ilp<4>/3502 53.3 ns 53.2 ns 12808066 bytes_per_second=245.455Gi/s
BM_sum_avx2_ilp<4>/35023 1038 ns 1035 ns 641606 bytes_per_second=126.009Gi/s
BM_sum_avx2_ilp<4>/350234 15468 ns 15420 ns 45387 bytes_per_second=84.6099Gi/s
BM_sum_std_fold_left/35 0.989 ns 0.986 ns 705266891 bytes_per_second=132.25Gi/s
BM_sum_std_fold_left/350 8.69 ns 8.65 ns 78872493 bytes_per_second=150.665Gi/s
BM_sum_std_fold_left/3502 92.9 ns 92.1 ns 7300433 bytes_per_second=141.638Gi/s
BM_sum_std_fold_left/35023 1041 ns 1030 ns 672576 bytes_per_second=126.681Gi/s
BM_sum_std_fold_left/350234 15312 ns 15218 ns 54346 bytes_per_second=85.7346Gi/s
NOTE: 我不是很确定 ILP 版本的 245 Gi/s 是怎么做到的，手写 benchmark 的话是和 fold_left 基本没有差异。
另外可以看出，此时 AVX 版本已经和 ILP 版本没有差异。也即是说，即使只用单个 YMM 寄存器，Zen 3 处理器的硬件调度仍然是非常可靠的，在编译器优化到位的前提下，不需要显式手写 ILP 设计。
总结
基本上没别的事情了，因为向量求和算法很固定，所以自动向量化本来挺好的。
但是如果考虑手写/混用 SIMD intrinsic 的话，开发者还是要有手段去控制编译器的优化行为。constexpr_for 是一个极好的控制 ILP 的工具，它是现代且可移植的写法。
处理器的硬件调度能力也是个考虑点，Zen 3 的表现就是优化（暴力展开）到位了，显式 ILP 无所谓。但是前面在你不知道性能瓶颈的时候，简单的显式 ILP 设计就能把性能提升数倍，还不需要做剖析（本文有基线，省了 perf 环节）。考虑到其他更加复杂的算法和延迟更高的向量指令，显式 ILP 的做法仍然是非常的甜点。
std::ranges 的组合玩法很好，但是不看汇编的话你不知道有什么坑。本文中的 ILP 程序性能瓶颈就是 ranges 导致的。