缓存友好的离线随机读取
Cache friendly offline random read
在 C++ 中考虑这个函数：
void foo(uint32_t *a1, uint32_t *a2, uint32_t *b1, uint32_t *b2, uint32_t *o) {
while (b1 != b2) {
// assert(0 <= *b1 && *b1 < a2 - a1)
*o++ = a1[*b1++];
}
}
目的应该很明确了吧。不幸的是，b1 包含 随机 数据并丢弃缓存，使 foo 成为我程序的瓶颈。无论如何我可以优化它吗？
这是一个 SSCCE，应该类似于我的实际代码：
#include <iostream>
#include <chrono>
#include <algorithm>
#include <numeric>
namespace {
void foo(uint32_t *a1, uint32_t *a2, uint32_t *b1, uint32_t *b2, uint32_t *o) {
while (b1 != b2) {
// assert(0 <= *b1 && *b1 < a2 - a1)
*o++ = a1[*b1++];
}
}
constexpr unsigned max_n = 1 << 24, max_q = 1 << 24;
uint32_t data[max_n], index[max_q], result[max_q];
}
int main() {
uint32_t seed = 0;
auto rng = [&seed]() { return seed = seed * 9301 + 49297; };
std::generate_n(data, max_n, rng);
std::generate_n(index, max_q, [rng]() { return rng() % max_n; });
auto t1 = std::chrono::high_resolution_clock::now();
foo(data, data + max_n, index, index + max_q, result);
auto t2 = std::chrono::high_resolution_clock::now();
std::cout << std::chrono::duration<double>(t2 - t1).count() << std::endl;
uint32_t hash = 0;
for (unsigned i = 0; i < max_q; i++)
hash += result[i] ^ (i << 8) ^ i;
std::cout << hash << std::endl;
}
这不是 ，它询问随机写入并假设 b 是一个排列。
您可以将索引划分到索引的较高位相同的桶中。请注意，如果索引不是随机的，则桶会溢出。
#include <iostream>
#include <chrono>
#include <cassert>
#include <algorithm>
#include <numeric>
#include <vector>
namespace {
constexpr unsigned max_n = 1 << 24, max_q = 1 << 24;
void foo(uint32_t *a1, uint32_t *a2, uint32_t *b1, uint32_t *b2, uint32_t *o) {
while (b1 != b2) {
// assert(0 <= *b1 && *b1 < a2 - a1)
*o++ = a1[*b1++];
}
}
uint32_t* foo_fx(uint32_t *a1, uint32_t *a2, uint32_t *b1, uint32_t *b2, const uint32_t b_offset, uint32_t *o) {
while (b1 != b2) {
// assert(0 <= *b1 && *b1 < a2 - a1)
*o++ = a1[b_offset+(*b1++)];
}
return o;
}
uint32_t data[max_n], index[max_q], result[max_q];
std::pair<uint32_t, uint32_t[max_q / 8]>index_fx[16];
}
int main() {
uint32_t seed = 0;
auto rng = [&seed]() { return seed = seed * 9301 + 49297; };
std::generate_n(data, max_n, rng);
//std::generate_n(index, max_q, [rng]() { return rng() % max_n; });
for (size_t i = 0; i < max_q;++i) {
const uint32_t idx = rng() % max_n;
const uint32_t bucket = idx >> 20;
assert(bucket < 16);
index_fx[bucket].second[index_fx[bucket].first] = idx % (1 << 20);
index_fx[bucket].first++;
assert((1 << 20)*bucket + index_fx[bucket].second[index_fx[bucket].first - 1] == idx);
}
auto t1 = std::chrono::high_resolution_clock::now();
//foo(data, data + max_n, index, index + max_q, result);
uint32_t* result_begin = result;
for (int i = 0; i < 16; ++i) {
result_begin = foo_fx(data, data + max_n, index_fx[i].second, index_fx[i].second + index_fx[i].first, (1<<20)*i, result_begin);
}
auto t2 = std::chrono::high_resolution_clock::now();
std::cout << std::chrono::duration<double>(t2 - t1).count() << std::endl;
std::cout << std::accumulate(result, result + max_q, 0ull) << std::endl;
}
首先我们来看一下上面代码的实际表现：
$ sudo perf stat ./offline-read
0.123023
1451229184
Performance counter stats for './offline-read':
184.661547 task-clock (msec) # 0.997 CPUs utilized
3 context-switches # 0.016 K/sec
0 cpu-migrations # 0.000 K/sec
717 page-faults # 0.004 M/sec
623,638,834 cycles # 3.377 GHz
419,309,952 instructions # 0.67 insn per cycle
70,803,672 branches # 383.424 M/sec
16,895 branch-misses # 0.02% of all branches
0.185129552 seconds time elapsed
我们的 IPC 很低，为 0.67，这可能几乎完全是由于 DRAM 加载未命中造成的5。让我们确认一下：
sudo ../pmu-tools/ocperf.py stat -e cycles,LLC-load-misses,cycle_activity.stalls_l3_miss ./offline-read
perf stat -e cycles,LLC-load-misses,cpu/event=0xa3,umask=0x6,cmask=6,name=cycle_activity_stalls_l3_miss/ ./offline-read
0.123979
1451229184
Performance counter stats for './offline-read':
622,661,371 cycles
16,114,063 LLC-load-misses
368,395,404 cycle_activity_stalls_l3_miss
0.184045411 seconds time elapsed
所以 620k 中的约 370k 周期直接因未命中而停滞。事实上，foo() 中以这种方式停滞的周期部分要高得多，接近 90%，因为 perf 也在测量 init 和 accumulate 代码，它们占用了大约三分之一的运行时间（但没有明显的 L3 未命中）。
这并不意外，因为我们知道随机读取模式 a1[*b1++] 基本上将具有零局部性。事实上，LLC-load-misses的数量是1600万1，几乎正好对应于a1的1600万随机读取。2
如果我们假设 foo() 的 100% 花费在等待内存访问上，我们可以了解每次未命中的总成本：0.123 sec / 16,114,063 misses == 7.63 ns/miss。在我的机器上，最佳情况下内存延迟约为 60 ns，因此每次未命中少于 8 ns 意味着我们已经提取了大量内存级并行性 (MLP)：大约 8 次未命中必须重叠并在 -平均飞行以实现这一目标（甚至完全忽略来自 b1 的流式负载和 o 的流式写入的额外流量）。
所以我认为没有多少调整可以应用于简单循环以做得更好。不过，两种可能性是：
- 写入
o 的非临时存储，如果您的平台支持的话。这将删除 RFO 隐含的正常存储的读取。这应该是一场直接的胜利，因为 o 再也不会被阅读（在计时部分内！）。
软件预取。仔细调整 a1 或 b1 的预取可能会对 有点帮助 。然而，影响将相当有限，因为如上所述，我们已经接近 MLP 的极限。此外，我们希望 b1 的线性读取能够被硬件预取器几乎完美地预取。 a1 的随机读取似乎可以进行预取，但实际上循环中的 ILP 通过乱序处理导致足够的 MLP（至少在像最近的 x86 这样的大型 OoO 处理器上）。
在评论中，用户 harold 已经提到他尝试过预取，但效果很小。
因此，由于简单的调整不太可能产生太多成果，因此您只能对循环进行改造。一个 "obvious" 转换是对索引 b1（以及索引元素的原始位置）进行排序，然后按排序顺序从 a1 读取。这将 a1 的读取从完全随机转换为几乎 3 线性，但现在写入都是随机的，这也好不到哪儿去。
排序然后取消排序
关键问题是 b1 控制下的 a1 的读取是随机的，并且 a1 很大，基本上每次读取都会丢失到 DRAM。我们可以通过排序 b1 来解决这个问题，然后读取 a1 以获得排列后的结果。现在需要"un-permute"结果a1得到最终顺序的结果，这只是另一种排序，这次在"output index".
这是一个使用给定输入数组 a、索引数组 b 和输出数组 o 以及 i 的工作示例，它是每个元素：
i = 0 1 2 3
a = [00, 10, 20, 30]
b = [ 3, 1, 0, 1]
o = [30, 10, 00, 10] (desired result)
首先，对数组 b 进行排序，将原始数组位置 i 作为辅助数据（或者您可以将其视为排序元组 (b[0], 0), (b[1], 1), ...），这将为您提供排序后的 b数组b'和排序后的索引列表i'如图：
i' = [ 2, 1, 3, 0]
b' = [ 0, 1, 1, 3]
现在你可以在b'的控制下从a读取排列结果数组o'。此读取按顺序严格递增，并且应该能够以接近 memcpy 的速度运行。事实上，您可以利用广泛连续的 SIMD 读取和一些随机播放来进行多次读取，然后将 4 字节元素移动到正确的位置（复制一些元素并跳过其他元素）：
a = [00, 10, 20, 30]
b' = [ 0, 1, 1, 3]
o' = [00, 10, 10, 30]
最后，您取消排列 o' 以获得 o，从概念上讲，只需对排列后的索引 i':
进行排序 o'
i' = [ 2, 1, 3, 0]
o' = [00, 10, 10, 30]
i = [ 0, 1, 2, 3]
o = [30, 10, 00, 10]
完成！
现在这是该技术最简单的想法，并不是特别适合缓存（从概念上讲，每次遍历都遍历一个或多个 2^26 字节数组），但它至少完全使用了它读取的每个缓存行（与仅从缓存行读取单个元素的原始循环不同，这就是为什么即使数据仅占用 100 万个缓存行，您也会有 1600 万次未命中！）。所有读取或多或少都是线性的，因此硬件预取会有很大帮助。
您获得多少加速可能取决于您将如何实现排序：它们需要快速且缓存敏感。几乎可以肯定，某种类型的缓存感知基数排序最有效。
这里有一些关于进一步改进方法的说明：
优化排序量
您实际上不需要完全排序 b。您只想对其进行排序 "enough" ，以便在 b' 的控制下对 a 的后续读取或多或少是线性的。例如，16 个元素适合一个缓存行，因此您根本不需要根据最后 4 位进行排序：无论如何都会读取相同线性序列的缓存行。您还可以对更少的位进行排序：例如，如果您忽略了 5 个最低有效位，您将以 "almost linear" 方式读取缓存行，有时会从完全线性模式交换两个缓存行，例如：0, 1, 3, 2, 5, 4, 6, 7。在这里，您仍然可以获得 L1 缓存的全部好处（后续读取缓存行将始终命中），我怀疑这样的模式仍然可以很好地预取，如果不能，您总是可以通过软件预取来帮助它。
您可以在您的系统上测试最佳忽略位数是多少。忽略位有两个好处：
- 基数搜索中要做的工作更少，因为需要的遍数更少，或者在一次或多次遍历中需要更少的桶（这有助于缓存）。
- 对"undo"最后一步中的排列所做的工作可能更少：如果通过检查原始索引数组
b撤消，忽略位意味着您在撤消时获得相同的节省搜索。
缓存块工作
上面的描述在几个连续的、不相交的通道中列出了所有内容，每个通道都在整个数据集上工作。实际上，您可能希望交错使用它们以获得更好的缓存行为。例如，假设您使用 MSD radix-256 排序，您可能会执行第一遍，将数据排序到 256 个桶中，每个桶大约有 256K 个元素。
然后，您可能只完成对第一个（或前几个）桶的排序，而不是进行完整的第二遍，然后根据 [=49= 的结果块继续读取 a ].你保证这个块是连续的（即，最终排序序列的后缀）所以你不会放弃读取中的任何位置，并且你的读取通常会被缓存。您也可以进行第一次反排列 o'，因为 o' 的块在缓存中也很热（也许您可以将后两个阶段组合成一个循环）。
智能反排列
优化的一个方面是 o' 的去排列的具体实现方式。在上面的描述中，我们假设一些索引数组 i 最初具有值 [0, 1, 2, ..., max_q]，它与 b 一起排序。这就是 概念上 它的工作原理，但您可能不需要立即具体化 i 并将其分类为辅助数据。例如，在基数排序的第一遍中，i 的值是隐式已知的（因为您正在遍历数据），因此可以免费计算 4 并在第一遍中写出，而没有按排序顺序出现。
可能还有比维护完整索引更有效的方法来执行 "unsort" 操作。例如，原始的未排序 b 数组在概念上具有执行取消排序所需的所有信息，但我很清楚如何使用它来有效地取消排序。
会更快吗？
那么这实际上会比天真的方法更快吗？它在很大程度上取决于实施细节，尤其是实施排序的效率。在我的硬件上，天真的方法是每秒处理大约 1.4 亿个元素。缓存感知基数排序的在线描述似乎从 200 到 6 亿 elements/s 不等，并且由于您需要其中两个，如果您相信这些数字，则大幅加速的机会似乎有限。另一方面，这些数字来自较旧的硬件，并且用于更一般的搜索（例如，对于密钥的所有 32 位，而我们可以使用少至 16 位）。
只有仔细实施才能确定是否可行，可行性还取决于硬件。例如，在无法支持尽可能多的 MLP 的硬件上，排序-取消排序方法变得相对更有利。
最佳方法还取决于 max_n 和 max_q 的相对值。例如，如果 max_n >> max_q，那么即使采用最佳排序，读取也将是 "sparse"，因此朴素的方法会更好。另一方面，如果 max_n << max_q，那么同一个索引通常会被读取多次，因此排序方法将具有良好的读取局部性，排序步骤本身将具有更好的局部性，并且进一步优化显式处理重复读取可能有可能。
多核
从问题中不清楚您是否有兴趣将其并行化。 foo() 的简单解决方案已经承认 "straightforward" 并行化，您只需将 a 和 b 数组划分为大小相等的块，对每个线程，这似乎提供完美的加速。不幸的是，你可能会发现你比线性缩放更糟糕，因为你会 运行 进入内存控制器和相关的 uncore/offcore 资源的资源争用，这些资源在套接字上的所有内核之间共享.因此，尚不清楚当您添加更多内核时，纯并行随机读取内存的吞吐量会增加多少6.
对于基数排序版本，大部分瓶颈（存储吞吐量、总指令吞吐量）都在内核中，因此我希望它可以通过额外的内核进行合理扩展。正如 Peter 在评论中提到的那样，如果您使用超线程，排序可能会在核心本地 L1 和 L2 缓存中具有良好的局部性，从而有效地让每个兄弟线程使用整个缓存，而不是将有效容量减半.当然，这涉及仔细管理您的线程亲和性，以便兄弟线程实际使用附近的数据，而不是让调度程序做任何它做的事情。
1 你可能会问为什么 LLC-load-misses 不是 32 或 4800 万，因为我们还必须读取 [=24] 的所有 1600 万个元素=] 然后 accumulate() 调用读取所有 result。答案是 LLC-load-misses 只计算 L3 中 实际上 未命中的需求未命中数。其他提到的读取模式是完全线性的，因此预取器将始终在需要之前将行带入 L3。根据 perf 使用的定义，这些不算作 "LLC misses"。
2 您可能想知道我是如何知道 加载未命中的全部来自 a1 的读取foo：我只是使用 perf record 和 perf mem 来确认未命中来自预期的汇编指令。
3几乎线性因为b1不是所有索引的排列，所以原则上可以跳过和重复索引。然而，在缓存行级别，很可能每个缓存行都会按顺序读取，因为每个元素都有大约 63% 的机会被包含，并且一个缓存行有 16 个 4 字节元素，所以只有任何给定缓存具有 零 个元素的概率约为千万分之一。因此，在高速缓存行级别工作的预取将正常工作。
4 这里我的意思是这个值的计算是免费的或者几乎是免费的，当然写入还是要花钱的。这仍然比 "up-front materialization" 方法好得多，但是，它首先创建需要 max_q 写入的 i 数组 [0, 1, 2, ...] 然后再次需要另一个 max_q 写入在第一个基数排序通道中对其进行排序。隐式具体化仅引起第二次写入。
5其实实际计时段foo()的IPC要低很多：根据我的大概是0.15计算。整个流程上报的IPC是定时段IPC的平均值，前后的初始化累加代码IPC高很多。
6 值得注意的是，这与依赖负载延迟限制工作流的扩展方式不同：一种负载正在进行随机读取，但只能有一个正在进行的负载，因为每个负载取决于上次扩展到多个内核的结果，因为负载的串行性质不会使用很多下游资源（但从概念上讲，即使在单个内核上也可以通过更改内核循环来处理更多资源来加速此类负载多于一个并联的相关负载流）。
在 C++ 中考虑这个函数：
void foo(uint32_t *a1, uint32_t *a2, uint32_t *b1, uint32_t *b2, uint32_t *o) {
while (b1 != b2) {
// assert(0 <= *b1 && *b1 < a2 - a1)
*o++ = a1[*b1++];
}
}
目的应该很明确了吧。不幸的是，b1 包含 随机 数据并丢弃缓存，使 foo 成为我程序的瓶颈。无论如何我可以优化它吗？
这是一个 SSCCE，应该类似于我的实际代码：
#include <iostream>
#include <chrono>
#include <algorithm>
#include <numeric>
namespace {
void foo(uint32_t *a1, uint32_t *a2, uint32_t *b1, uint32_t *b2, uint32_t *o) {
while (b1 != b2) {
// assert(0 <= *b1 && *b1 < a2 - a1)
*o++ = a1[*b1++];
}
}
constexpr unsigned max_n = 1 << 24, max_q = 1 << 24;
uint32_t data[max_n], index[max_q], result[max_q];
}
int main() {
uint32_t seed = 0;
auto rng = [&seed]() { return seed = seed * 9301 + 49297; };
std::generate_n(data, max_n, rng);
std::generate_n(index, max_q, [rng]() { return rng() % max_n; });
auto t1 = std::chrono::high_resolution_clock::now();
foo(data, data + max_n, index, index + max_q, result);
auto t2 = std::chrono::high_resolution_clock::now();
std::cout << std::chrono::duration<double>(t2 - t1).count() << std::endl;
uint32_t hash = 0;
for (unsigned i = 0; i < max_q; i++)
hash += result[i] ^ (i << 8) ^ i;
std::cout << hash << std::endl;
}
这不是 b 是一个排列。
您可以将索引划分到索引的较高位相同的桶中。请注意，如果索引不是随机的，则桶会溢出。
#include <iostream>
#include <chrono>
#include <cassert>
#include <algorithm>
#include <numeric>
#include <vector>
namespace {
constexpr unsigned max_n = 1 << 24, max_q = 1 << 24;
void foo(uint32_t *a1, uint32_t *a2, uint32_t *b1, uint32_t *b2, uint32_t *o) {
while (b1 != b2) {
// assert(0 <= *b1 && *b1 < a2 - a1)
*o++ = a1[*b1++];
}
}
uint32_t* foo_fx(uint32_t *a1, uint32_t *a2, uint32_t *b1, uint32_t *b2, const uint32_t b_offset, uint32_t *o) {
while (b1 != b2) {
// assert(0 <= *b1 && *b1 < a2 - a1)
*o++ = a1[b_offset+(*b1++)];
}
return o;
}
uint32_t data[max_n], index[max_q], result[max_q];
std::pair<uint32_t, uint32_t[max_q / 8]>index_fx[16];
}
int main() {
uint32_t seed = 0;
auto rng = [&seed]() { return seed = seed * 9301 + 49297; };
std::generate_n(data, max_n, rng);
//std::generate_n(index, max_q, [rng]() { return rng() % max_n; });
for (size_t i = 0; i < max_q;++i) {
const uint32_t idx = rng() % max_n;
const uint32_t bucket = idx >> 20;
assert(bucket < 16);
index_fx[bucket].second[index_fx[bucket].first] = idx % (1 << 20);
index_fx[bucket].first++;
assert((1 << 20)*bucket + index_fx[bucket].second[index_fx[bucket].first - 1] == idx);
}
auto t1 = std::chrono::high_resolution_clock::now();
//foo(data, data + max_n, index, index + max_q, result);
uint32_t* result_begin = result;
for (int i = 0; i < 16; ++i) {
result_begin = foo_fx(data, data + max_n, index_fx[i].second, index_fx[i].second + index_fx[i].first, (1<<20)*i, result_begin);
}
auto t2 = std::chrono::high_resolution_clock::now();
std::cout << std::chrono::duration<double>(t2 - t1).count() << std::endl;
std::cout << std::accumulate(result, result + max_q, 0ull) << std::endl;
}
首先我们来看一下上面代码的实际表现：
$ sudo perf stat ./offline-read
0.123023
1451229184
Performance counter stats for './offline-read':
184.661547 task-clock (msec) # 0.997 CPUs utilized
3 context-switches # 0.016 K/sec
0 cpu-migrations # 0.000 K/sec
717 page-faults # 0.004 M/sec
623,638,834 cycles # 3.377 GHz
419,309,952 instructions # 0.67 insn per cycle
70,803,672 branches # 383.424 M/sec
16,895 branch-misses # 0.02% of all branches
0.185129552 seconds time elapsed
我们的 IPC 很低，为 0.67，这可能几乎完全是由于 DRAM 加载未命中造成的5。让我们确认一下：
sudo ../pmu-tools/ocperf.py stat -e cycles,LLC-load-misses,cycle_activity.stalls_l3_miss ./offline-read
perf stat -e cycles,LLC-load-misses,cpu/event=0xa3,umask=0x6,cmask=6,name=cycle_activity_stalls_l3_miss/ ./offline-read
0.123979
1451229184
Performance counter stats for './offline-read':
622,661,371 cycles
16,114,063 LLC-load-misses
368,395,404 cycle_activity_stalls_l3_miss
0.184045411 seconds time elapsed
所以 620k 中的约 370k 周期直接因未命中而停滞。事实上，foo() 中以这种方式停滞的周期部分要高得多，接近 90%，因为 perf 也在测量 init 和 accumulate 代码，它们占用了大约三分之一的运行时间（但没有明显的 L3 未命中）。
这并不意外，因为我们知道随机读取模式 a1[*b1++] 基本上将具有零局部性。事实上，LLC-load-misses的数量是1600万1，几乎正好对应于a1的1600万随机读取。2
如果我们假设 foo() 的 100% 花费在等待内存访问上，我们可以了解每次未命中的总成本：0.123 sec / 16,114,063 misses == 7.63 ns/miss。在我的机器上，最佳情况下内存延迟约为 60 ns，因此每次未命中少于 8 ns 意味着我们已经提取了大量内存级并行性 (MLP)：大约 8 次未命中必须重叠并在 -平均飞行以实现这一目标（甚至完全忽略来自 b1 的流式负载和 o 的流式写入的额外流量）。
所以我认为没有多少调整可以应用于简单循环以做得更好。不过，两种可能性是：
- 写入
o的非临时存储，如果您的平台支持的话。这将删除 RFO 隐含的正常存储的读取。这应该是一场直接的胜利，因为o再也不会被阅读（在计时部分内！）。 软件预取。仔细调整
a1或b1的预取可能会对 有点帮助 。然而，影响将相当有限，因为如上所述，我们已经接近 MLP 的极限。此外，我们希望b1的线性读取能够被硬件预取器几乎完美地预取。a1的随机读取似乎可以进行预取，但实际上循环中的 ILP 通过乱序处理导致足够的 MLP（至少在像最近的 x86 这样的大型 OoO 处理器上）。在评论中，用户 harold 已经提到他尝试过预取，但效果很小。
因此，由于简单的调整不太可能产生太多成果，因此您只能对循环进行改造。一个 "obvious" 转换是对索引 b1（以及索引元素的原始位置）进行排序，然后按排序顺序从 a1 读取。这将 a1 的读取从完全随机转换为几乎 3 线性，但现在写入都是随机的，这也好不到哪儿去。
排序然后取消排序
关键问题是 b1 控制下的 a1 的读取是随机的，并且 a1 很大，基本上每次读取都会丢失到 DRAM。我们可以通过排序 b1 来解决这个问题，然后读取 a1 以获得排列后的结果。现在需要"un-permute"结果a1得到最终顺序的结果，这只是另一种排序，这次在"output index".
这是一个使用给定输入数组 a、索引数组 b 和输出数组 o 以及 i 的工作示例，它是每个元素：
i = 0 1 2 3
a = [00, 10, 20, 30]
b = [ 3, 1, 0, 1]
o = [30, 10, 00, 10] (desired result)
首先，对数组 b 进行排序，将原始数组位置 i 作为辅助数据（或者您可以将其视为排序元组 (b[0], 0), (b[1], 1), ...），这将为您提供排序后的 b数组b'和排序后的索引列表i'如图：
i' = [ 2, 1, 3, 0]
b' = [ 0, 1, 1, 3]
现在你可以在b'的控制下从a读取排列结果数组o'。此读取按顺序严格递增，并且应该能够以接近 memcpy 的速度运行。事实上，您可以利用广泛连续的 SIMD 读取和一些随机播放来进行多次读取，然后将 4 字节元素移动到正确的位置（复制一些元素并跳过其他元素）：
a = [00, 10, 20, 30]
b' = [ 0, 1, 1, 3]
o' = [00, 10, 10, 30]
最后，您取消排列 o' 以获得 o，从概念上讲，只需对排列后的索引 i':
o'
i' = [ 2, 1, 3, 0]
o' = [00, 10, 10, 30]
i = [ 0, 1, 2, 3]
o = [30, 10, 00, 10]
完成！
现在这是该技术最简单的想法，并不是特别适合缓存（从概念上讲，每次遍历都遍历一个或多个 2^26 字节数组），但它至少完全使用了它读取的每个缓存行（与仅从缓存行读取单个元素的原始循环不同，这就是为什么即使数据仅占用 100 万个缓存行，您也会有 1600 万次未命中！）。所有读取或多或少都是线性的，因此硬件预取会有很大帮助。
您获得多少加速可能取决于您将如何实现排序：它们需要快速且缓存敏感。几乎可以肯定，某种类型的缓存感知基数排序最有效。
这里有一些关于进一步改进方法的说明：
优化排序量
您实际上不需要完全排序 b。您只想对其进行排序 "enough" ，以便在 b' 的控制下对 a 的后续读取或多或少是线性的。例如，16 个元素适合一个缓存行，因此您根本不需要根据最后 4 位进行排序：无论如何都会读取相同线性序列的缓存行。您还可以对更少的位进行排序：例如，如果您忽略了 5 个最低有效位，您将以 "almost linear" 方式读取缓存行，有时会从完全线性模式交换两个缓存行，例如：0, 1, 3, 2, 5, 4, 6, 7。在这里，您仍然可以获得 L1 缓存的全部好处（后续读取缓存行将始终命中），我怀疑这样的模式仍然可以很好地预取，如果不能，您总是可以通过软件预取来帮助它。
您可以在您的系统上测试最佳忽略位数是多少。忽略位有两个好处：
- 基数搜索中要做的工作更少，因为需要的遍数更少，或者在一次或多次遍历中需要更少的桶（这有助于缓存）。
- 对"undo"最后一步中的排列所做的工作可能更少：如果通过检查原始索引数组
b撤消，忽略位意味着您在撤消时获得相同的节省搜索。
缓存块工作
上面的描述在几个连续的、不相交的通道中列出了所有内容，每个通道都在整个数据集上工作。实际上，您可能希望交错使用它们以获得更好的缓存行为。例如，假设您使用 MSD radix-256 排序，您可能会执行第一遍，将数据排序到 256 个桶中，每个桶大约有 256K 个元素。
然后，您可能只完成对第一个（或前几个）桶的排序，而不是进行完整的第二遍，然后根据 [=49= 的结果块继续读取 a ].你保证这个块是连续的（即，最终排序序列的后缀）所以你不会放弃读取中的任何位置，并且你的读取通常会被缓存。您也可以进行第一次反排列 o'，因为 o' 的块在缓存中也很热（也许您可以将后两个阶段组合成一个循环）。
智能反排列
优化的一个方面是 o' 的去排列的具体实现方式。在上面的描述中，我们假设一些索引数组 i 最初具有值 [0, 1, 2, ..., max_q]，它与 b 一起排序。这就是 概念上 它的工作原理，但您可能不需要立即具体化 i 并将其分类为辅助数据。例如，在基数排序的第一遍中，i 的值是隐式已知的（因为您正在遍历数据），因此可以免费计算 4 并在第一遍中写出，而没有按排序顺序出现。
可能还有比维护完整索引更有效的方法来执行 "unsort" 操作。例如，原始的未排序 b 数组在概念上具有执行取消排序所需的所有信息，但我很清楚如何使用它来有效地取消排序。
会更快吗？
那么这实际上会比天真的方法更快吗？它在很大程度上取决于实施细节，尤其是实施排序的效率。在我的硬件上，天真的方法是每秒处理大约 1.4 亿个元素。缓存感知基数排序的在线描述似乎从 200 到 6 亿 elements/s 不等，并且由于您需要其中两个，如果您相信这些数字，则大幅加速的机会似乎有限。另一方面，这些数字来自较旧的硬件，并且用于更一般的搜索（例如，对于密钥的所有 32 位，而我们可以使用少至 16 位）。
只有仔细实施才能确定是否可行，可行性还取决于硬件。例如，在无法支持尽可能多的 MLP 的硬件上，排序-取消排序方法变得相对更有利。
最佳方法还取决于 max_n 和 max_q 的相对值。例如，如果 max_n >> max_q，那么即使采用最佳排序，读取也将是 "sparse"，因此朴素的方法会更好。另一方面，如果 max_n << max_q，那么同一个索引通常会被读取多次，因此排序方法将具有良好的读取局部性，排序步骤本身将具有更好的局部性，并且进一步优化显式处理重复读取可能有可能。
多核
从问题中不清楚您是否有兴趣将其并行化。 foo() 的简单解决方案已经承认 "straightforward" 并行化，您只需将 a 和 b 数组划分为大小相等的块，对每个线程，这似乎提供完美的加速。不幸的是，你可能会发现你比线性缩放更糟糕，因为你会 运行 进入内存控制器和相关的 uncore/offcore 资源的资源争用，这些资源在套接字上的所有内核之间共享.因此，尚不清楚当您添加更多内核时，纯并行随机读取内存的吞吐量会增加多少6.
对于基数排序版本，大部分瓶颈（存储吞吐量、总指令吞吐量）都在内核中，因此我希望它可以通过额外的内核进行合理扩展。正如 Peter 在评论中提到的那样，如果您使用超线程，排序可能会在核心本地 L1 和 L2 缓存中具有良好的局部性，从而有效地让每个兄弟线程使用整个缓存，而不是将有效容量减半.当然，这涉及仔细管理您的线程亲和性，以便兄弟线程实际使用附近的数据，而不是让调度程序做任何它做的事情。
1 你可能会问为什么 LLC-load-misses 不是 32 或 4800 万，因为我们还必须读取 [=24] 的所有 1600 万个元素=] 然后 accumulate() 调用读取所有 result。答案是 LLC-load-misses 只计算 L3 中 实际上 未命中的需求未命中数。其他提到的读取模式是完全线性的，因此预取器将始终在需要之前将行带入 L3。根据 perf 使用的定义，这些不算作 "LLC misses"。
2 您可能想知道我是如何知道 加载未命中的全部来自 a1 的读取foo：我只是使用 perf record 和 perf mem 来确认未命中来自预期的汇编指令。
3几乎线性因为b1不是所有索引的排列，所以原则上可以跳过和重复索引。然而，在缓存行级别，很可能每个缓存行都会按顺序读取，因为每个元素都有大约 63% 的机会被包含，并且一个缓存行有 16 个 4 字节元素，所以只有任何给定缓存具有 零 个元素的概率约为千万分之一。因此，在高速缓存行级别工作的预取将正常工作。
4 这里我的意思是这个值的计算是免费的或者几乎是免费的，当然写入还是要花钱的。这仍然比 "up-front materialization" 方法好得多，但是，它首先创建需要 max_q 写入的 i 数组 [0, 1, 2, ...] 然后再次需要另一个 max_q 写入在第一个基数排序通道中对其进行排序。隐式具体化仅引起第二次写入。
5其实实际计时段foo()的IPC要低很多：根据我的大概是0.15计算。整个流程上报的IPC是定时段IPC的平均值，前后的初始化累加代码IPC高很多。
6 值得注意的是，这与依赖负载延迟限制工作流的扩展方式不同：一种负载正在进行随机读取，但只能有一个正在进行的负载，因为每个负载取决于上次扩展到多个内核的结果，因为负载的串行性质不会使用很多下游资源（但从概念上讲，即使在单个内核上也可以通过更改内核循环来处理更多资源来加速此类负载多于一个并联的相关负载流）。