并行 O(√n) 开销 LSD 基数排序

# 并行 O(√n) 开销 LSD 基数排序  
来源：https://arxiv.org/html/2607.05302  
柏林数学应用科学研究所 (Zuse Institute Berlin)，德国  
[email protected]  
NHR-联盟研究生项目，柏林数学应用科学研究所，德国  
[email protected]  
https://orcid.org/0000-0003-4548-788X  
\{CCSXML\}¡ccs2012¿ ¡concept¿ ¡concept_id¿10003752.10003809.10010031.10010033¡/concept_id¿ ¡concept_desc¿计算理论 排序与搜索¡/concept_desc¿ ¡concept_significance¿500¡/concept_significance¿ ¡/concept¿ ¡/ccs2012¿  

本文作者之一曾在卡尔斯鲁厄理工学院 (KIT) 的 Sanders 课题组进行访问研究。作者在此感谢 Marvin Williams 提供的宝贵意见。

## 具有 O(\(\sqrt{n}\)) 开销的并行 LSD 基数排序

###### 摘要

我们介绍了 Radsort，一种 LSD 基数排序的变体，它以 O(\(\sqrt{n}\)) 的额外空间对数据进行排序。Radsort 是稳定的，实现简单且易于并行化。对于大小超过约 2 MiB 的数组，其性能优于传统的非原地 LSD 基数排序。

###### 分类主题

计算理论 排序与搜索

###### 关键词

缓存局部性、基数排序、排序

††runningauthor: R. Clausecker 和 F. Schintke  
††copyright: Robert Clausecker 和 Florian Schintke  
††supplement: https://github.com/clausecker/radsort  
††editors: John Q. Open 和 Joan R. Access  
††event-title: 第42届非常重要的主题会议 (CVIT 2016)  
††event-shorttitle: CVIT 2016  
††event-acronym: CVIT  
††year: 2016  
††event-date: 2016年12月24日至27日  
††event-location: Little Whinging，英国  
††series-volume: 42  
††articleno: 23  

## 1 引言

LSD（最低有效位）基数排序算法用于对一个长度为 \(n\) 的数组 \(A\) 进行排序，其中该数组包含 \(n\) 个键值的 \(n_t\) 元组 \((\textit{key}(A[i],0), \textit{key}(A[i],1), \ldots, \textit{key}(A[i],n_t-1))\)，这些键值来自大小为 \(\sigma\) 的字母表 \(\Sigma\)。它通过依次根据 \(\textit{key}(A[i],n_t-1)\)、然后 \(\textit{key}(A[i],n_t-2)\)，依此类推直到 \(\textit{key}(A[i],0)\) 进行排序，从而将数组排列成字典序。使用稳定排序时，相同的键会根据其后缀排序，从而得到字典序。传统的 LSD 基数排序首先对每个 \(c \in \Sigma\) 在 \(\textit{key}(A[i],t)\) 中出现的次数进行直方图统计 \(H\)。然后对 \(H\) 进行前缀和运算得到桶的起始位置 \(S\)，并据此将数组 \(A\) 排序到一个新数组 \(A'\) 中（算法1 (https://arxiv.org/html/2607.05302#alg1)）。

算法1：非原地 LSD 基数排序一轮
1. 过程 `radixSortStep(A', A, t)`
   ▷ 按键 t 对数组 A 排序，结果存入 A'
2.  \(H \leftarrow 0 \ldots 0\)
3.  对于 \(i \leftarrow 0 \ldots n-1\) 执行
   ▷ 统计 \(\textit{key}(A[i],t)\) 的直方图
4.      \(H[\textit{key}(A[i],t)] \leftarrow H[\textit{key}(A[i],t)] + 1\)
5.  结束循环
6.  \(S[0] \leftarrow 0\)
7.  对于 \(i \leftarrow 1 \ldots \sigma-1\) 执行
8.      \(S[i] \leftarrow S[i-1] + H[i-1]\)
   ▷ 对 H 进行前缀和
9.  结束循环
10. 对于 \(i \leftarrow 0 \ldots n-1\) 执行
11.     \(c \leftarrow \textit{key}(A[i],t)\)
   ▷ 确定桶
12.     \(A'[S[c]] \leftarrow A[i]\)
   ▷ 将 A[i] 排序到 A' 中
13.     \(S[c] \leftarrow S[c] + 1\)
   ▷ 桶 c 的位置前进，为新元素腾位
14. 结束循环
15. 结束过程

当输入数组 \(A\) 被消耗完时，输出被产生在一个单独的数组 \(A'\) 中。虽然在算法的任何时刻，两个数组中总共只占用 \(n\) 个元素，但仍然需要提供 \(2n\) 个元素的存储空间，即除了输入外还需要 \(n\) 个元素的额外开销。我们希望用一种新颖的 LSD 基数排序算法来解决这个缺点，该算法仅需 O(\(\sqrt{n}\)) 的额外空间就能完成工作，且性能相似。我们将此算法命名为 *Radsort*，因其具有“**根本性**开销”。

算法2：Radsort 算法
1. 过程 `radSort(A, n, n_t)`
2.     `setup(A, n)`
   ▷ 参见算法4 (https://arxiv.org/html/2607.05302#alg4)
3.     对于 \(t \leftarrow n_t-1 \ldots 0\) 执行
4.         `sortPhase(t)`
   ▷ 参见算法5 (https://arxiv.org/html/2607.05302#alg5)
5.     `fixupPhase`
   ▷ 参见算法6 (https://arxiv.org/html/2607.05302#alg6)
6.     结束循环
7.     `Finalize`
   ▷ 参见算法7 (https://arxiv.org/html/2607.05302#alg7)
8. 结束过程

## 2 Radsort

Radsort（算法2 (https://arxiv.org/html/2607.05302#alg2)）将输入数组 \(A\) 视为一系列大小为 \(b\) 的*块*（关于 \(b\) 的选择见第3节 (https://arxiv.org/html/2607.05302#S3)）。当一块输入被消耗后，其存储空间会被重复用于后续产生的输出，从而将存储开销减少到固定数量的*暂存块*以及一些簿记信息。

每一轮排序包含两个*阶段*。在*排序阶段*，我们首先为 \(\sigma\) 个我们用于排序的*桶*中的每一个分配一个输出块。一旦某个块被填满，就从先前消耗的输入块中取出一个新块，并在处理时覆盖输入。排序阶段的结果是 \(\sigma\) 个随机交织的输出块序列，每个序列构成一个输出桶。在每个*排序阶段*之后的*修正阶段*，这些序列被解交织，得到按桶顺序排列的数据块。解交织并不移动数据块——相反，一个置换 \(\pi\) 追踪了处理数据块的顺序，使得数据能按序出现。

与传统的 LSD 基数排序类似，这两个阶段会依次针对每个键位置重复进行。在所有 \(n_t\) 轮排序之后，一个*最终化步骤*会重新排列数据块，在数组 \(A\) 中给出排序结果。虽然 Radsort 分配空间的方式比传统的 LSD 基数排序更复杂，但它遵循相同的原则，在更低的内存开销和更好的缓存局部性下，提供了相同的稳定性保证。

### 2.1 数据结构

除了数组 \(A\) 的 \(\lfloor n/b \rfloor\) 个块之外，还需要存储在临时数组 \(T\) 中的 \(2\sigma\) 个*暂存块*：\(\sigma\) 个*领先*块用于覆盖当前轮次各桶的初始输出块，另外 \(\sigma\) 个块用于覆盖上一轮排序的部分块。该算法将 \(T\) 和 \(A\) 的连接视为一个由 \(n_\pi\) 个块组成的长数组，每个块包含 \(b\) 个元素，通过 `blockat` 函数进行索引。前 \(2\sigma\) 个索引引用 \(T\) 中的块，其余索引引用叠加在 \(A\) 上的块。`blockat(i)` 函数返回索引 \(i\) 处块的起始地址。
\[
\textit{blockat}(i)=
\begin{cases}
\text{address of } T[ib] & \text{if } 0 \leq i < 2\sigma \\
\text{address of } A[(i-2\sigma)b] & \text{if } 2\sigma \leq i < n_\pi
\end{cases}
\quad (1)
\]
然后，数组 \(A\) 末尾的 \(n \bmod b\) 个元素在算法开始时被复制到一个暂存块中。在算法结束之前，它们不会被使用，直到 \(T\) 的内容被复制回 \(A\)。数组 \(\pi\) 描述了这些 \(n_\pi\) 个块的一个置换，使得 \(\pi[i]\) 保存*逻辑索引* \(i\) 处的块索引。图1 (https://arxiv.org/html/2607.05302#S2.F1) 展示了初始化时对 \(n_t=4\) 个字符的字符串的此数据结构示例。

每个块要么是*已分配*的，要么是*未分配*的。已分配的块要么是*满的*，意味着其所有元素都持有数据；要么是*部分的*，意味着其部分元素持有数据。未分配的块不持有任何数据，可以重用。在最初的 \(\sigma\) 个逻辑索引处，存在一个由 \(\sigma\) 个未分配块组成的*领先*区域。它们后面是逻辑索引 \(\sigma \leq i < f\) 处的已分配块，最后是逻辑索引 \(f \leq i < n_\pi\) 处的未分配块。

数组 \(P\) 追踪最多 \(\sigma\) 个部分块。它按逻辑索引顺序保存部分块 \((i, l)\)，其中 \(i\) 是逻辑索引，\(0 \leq l < b\) 是块长度。对于每个这样的部分块，前 \(l\) 个元素被使用，其余 \(b-l\) 个元素未被使用。为简化实现，逻辑索引 \(f-1\) 处的块始终是一个部分块。

算法3：查找下一个块的大小
1. 函数 `sizeOfNextBlock(i, iP)`
2.     \((i_{i_P}, l_{i_P}) \leftarrow P[i_P]\)
3.     如果 \(\pi[i] = i_{i_P}\) 则
   ▷ \(\pi[i]\) 是否是一个部分块？
4.         \(i_P \leftarrow i_P + 1\)
5.         返回 \(l_{i_P}\)
6.     否则
7.         返回 \(b\)
8.     结束判断
9. 函数结束

按照 \(\pi\) 给定的顺序连接，满块和部分块的内容持有待排序的元素，其顺序由到目前为止执行的排序轮次所诱导。我们可以按逻辑顺序遍历已分配的块，同时通过遍历 \(P\) 来发现每个块的长度。这个想法被封装在函数 `sizeOfNextBlock(i, iP)` 中，它给出逻辑索引 \(i\) 处块的长度，假设下一个部分块由 \(P[i_P]\) 追踪。如果该块是部分的，则更新 \(i_P\) 以追踪下一个部分块。

总之，Radsort 算法的状态在以下变量中追踪。单个算法步骤内需要的额外变量在相应步骤中解释。

| 变量 | 描述 | 类型 |
| :--- | :--- | :--- |
| \(A\) | 输入数组 | 字符串数组 |
| \(n\) | 输入长度 | \(n = |A|\) |
| \(\sigma\) | 字母表大小 | \(\sigma = |\Sigma|\) |
| \(n_t\) | 键长度 | 整数 |
| \(b\) | 块大小 | 整数 |
| \(n_\pi\) | 块计数 | \(n_\pi = \lfloor n/b \rfloor + 2\sigma\) |
| \(T\) | 暂存缓冲区 | \(2\sigma b\) 个字符串的数组 |
| \(f\) | 填充级别 | 整数 |
| \(\pi\) | 置换 | \(n_\pi\) 个整数的数组 |
| \(P\) | 部分块 | \(\sigma\) 对整数的数组 |

**图1：** Radsort 初始化时的数据结构。一个包含 \(n_t=4\) 个字符的块，大小为 \(b\)。
[此处省略图1的ASCII图示]

### 2.2 初始化

变量初始化设置如下（参见图1 (https://arxiv.org/html/2607.05302#S2.F1)），使得数组 \(A\) 的内容由逻辑索引 \(\sigma\) 到 \(\sigma + \lceil n/b \rceil\) 的块表示。数组 \(A\) 最后的 \(n \bmod b\) 个元素被复制到块 \(\sigma\) 中。这建立了状态变量的*数据不变性*，如下所示。每一轮排序都假设此不变性在开始时成立，并在结束时重新建立。

1. \(\pi\) 是 \(0 \ldots n_\pi-1\) 的一个置换。
2. \(P\) 按逻辑索引升序追踪部分块。
3. 逻辑索引 \(f-1\) 处的块是一个部分块。
4. 逻辑索引 \(0\) 到 \(\sigma-1\) 以及 \(f\) 到 \(n_\pi-1\) 指向未分配块。
5. 逻辑索引 \(\sigma\) 到 \(f-1\) 指向已分配块，这些块中已使用元素的连接持有根据算法当前进度置换后的输入。

算法4：变量的初始化
1. 过程 `setup(A, n)`
2.     \(f \leftarrow \sigma + \lfloor n/b \rfloor + 1\)
3.     \(\pi[0 \ldots \sigma-1] \leftarrow 0, \ldots, \sigma-1\)
   ▷ 分配领先块
4.     \(\pi[\sigma \ldots f-2] \leftarrow 2\sigma, \ldots, \sigma+f-2\)
   ▷ 分配数组 A 的内容
5.     \(\pi[f-1 \ldots n_\pi-1] \leftarrow \sigma, \ldots, 2\sigma-1\)
   ▷ 分配剩余块
6.     \(T[0 \dots 2\sigma b-1] \leftarrow 0\)
   ▷ 用虚拟值填充 T（可选）
7.     \(\textit{blockat}(\sigma)[0 \dots (n \bmod b)-1] \leftarrow A[b \lfloor n/b \rfloor \ldots n-1]\)
   ▷ 将尾部数据复制到暂存块
8.     \(P[0] \leftarrow (f-1, n \bmod b)\)
   ▷ 将输入尾部追踪为部分块
9.     \(P[1 \ldots \sigma-1] \leftarrow (n_\pi, 0)\)
   ▷ 用虚拟值填充 P（可选）
10. 结束过程

### 2.3 排序

每一轮排序（参见图2 (https://arxiv.org/html/2607.05302#S2.F2)）根据某个键位置 \(t\) 对数组 \(A\) 的元素进行稳定排序。我们用 \(\textit{key}(A[i],t) \in \Sigma\) 表示元素 \(A[i]\) 在键位置 \(t\) 的值。首先，*排序阶段*根据键值将元素分配到桶中。然后，*修正阶段*根据排序阶段的结果计算新的 \(\pi\) 和 \(P\)，恢复数据不变性。排序和修正阶段需要几个额外的变量，这些变量在之后可以丢弃：

| 变量 | 描述 | 类型 |
| :--- | :--- | :--- |
| \(B\) | 桶分配 | \(\sigma\) 对指针的数组 |
| \(C\) | 块计数 | \(\sigma\) 个整数的数组 |
| \(S\) | 桶起始位置 | \(\sigma\) 个整数的数组 |
| \(U\) | 块使用情况 | \(n_\pi\) 个字符的数组 |
| \(\pi'\) | 新置换 | \(n_\pi\) 个整数的数组 |
| \(i_{\text{in}}\) | 下一个输入块 | 整数 |
| \(i_{\text{out}}\) | 下一个输出块 | 整数 |
| \(i_P\) | 下一个部分块 | 整数 |

#### 2.3.1 排序阶段

排序阶段在桶数组 \(B\) 中为每个桶分配一个块。该数组为每个桶 \(c\) 保存一对指针 \((p_{\text{next}}, p_{\text{end}})\)，分别指向下一个空闲元素和块的末尾。然后它通过 \(i_{\text{in}}\) 按逻辑顺序遍历已分配的块，并对每个块的元素进行排序，放入正确的输出块中。如果一个输出块满了，就使用 `newBucketBlock(c)` 过程从下一个逻辑索引 \(i_{\text{out}}\) 中取出一个新块。每个输出块所用于的桶由使用数组 \(U\) 追踪，而分配给每个桶的块数则由 \(C\) 追踪。

随着排序过程在输入数组中进行，已消耗的输入块会被重新分配为输出块。关键在于，在将输入块重新分配为输出块之前，必须完全消耗这些输入块，否则元素可能在排序之前就被覆盖。这通过...