Caiwen的博客

软件性能工程-缓存高效算法

2025-09-28 00:59

1. Ideal-Cache Model

1.1 Overview

后续我们对缓存高效算法的分析就基于这样一个 Ideal-Cache Model 上:只有一个 cache set,缓存的 大小为 MM,cache line 存储的块大小为 BB,cache line 的数量总共是 MB\frac{M}{B}。只有这一层缓存。并且我们假设如果出现了 capacity miss 的时候,我们永远能做出最优的替换选择(也就是假设我们知道后续的内存访问情况)。

注意,由于我们这里只有一个 cache set,所以 ideal-cache model 并没有考虑 conflict miss 的情况,这和实际还是有差别的。

image-20250928091614831

然后我们有 LRU Lemma:假设一个程序在大小为 MM 的 ideal-cache 上发生了 QQ 次的缓存不命中,那么在一个大小为 2M2M 全相连缓存(即只有一个 cache set)的缓存上,使用 LRU 的缓存替换策略,最多只会发生 2Q2Q 次的缓存不命中。这也就意味着在全相连缓存上,我们可以近似认为 LRU 就是最优策略。

Tall-cache assumption:我们一律假设我们的缓存是比较高的,也就是 MB>B\frac{M}{B}>B,也就是 M>B2M>B^2。这不仅符合现实情况,而且理论上也应该是这样。感性地理解为我们缓存肯定是少量多次地去覆盖工作集。考虑说如果我们有一个 short cache,那么可能一个 n2<Mn^2<M 的 submatrix 都无法放进缓存中。

image-20250928094949406

1.2 并行

下面的分析基于使用 work stealing 的调度策略的。

我们令 QPQ_P 为在 PP 个核心上并行执行时总共产生的缓存不命中次数,SpS_p 表示 work stealing 策略成功 steal 的次数。那么我们有:

QP=Q1+O(SPMB)Q_P = Q_1 + O(\frac{S_PM}{B})

因为每次 steal 得来的任务最多需要 MB\frac{M}{B} 次缓存不命中使得该核上的缓存变为任务被 steal 前的缓存状态。

我们之前分析过,SPS_P 在期望上大概是 O(PT)O(PT_{\infin}) 的。

这意味着,减少并行时的缓存不命中次数,就等同于减少顺序执行时的缓存不命中次数。

1.3 其他

Cache-Miss Lemma:假设一个程序从内存中读取了 rr 段,第 ii 段包含 sis_i 字节。总读取量 N=i=1rsiN=\sum_{i=1}^r s_i。程序使用的是总大小为 MM,cache line 大小为 BB 的 ideal-cache。

如果 N<M3N < \frac{M}{3}NrB\frac{N}{r} \ge B,那么我们说这 rr 段数据都能放进缓存中,且整个读取过程中产生的缓存不命中次数为 3NB\frac{3N}{B}

证明:直接顺序读取大小为 sis_i 的数据,产生的缓存不命中次数最多是 siB+2\frac{s_i}{B}+2,如图:

image-20250928095054152

于是:

i=1rsiB+2=NB+2r=NB+2rBBNB+2NB=3NB\begin{align} \sum_{i=1}^r \frac{s_i}{B}+2 &=\frac{N}{B}+2r\\ &=\frac{N}{B}+\frac{2rB}{B}\\ &\le \frac{N}{B}+\frac{2N}{B}\\ &=\frac{3N}{B} \end{align}

Submatrix Caching Lemma:如果一个 n×nn\times n 的 submatrix 满足 cMn2<M3cM\le n^2 < \frac{M}{3}c1c\le 1),那么对于这个 submatrix 是能够被放进 cache 里的,并且遍历一遍后缓存不命中次数最多是 3n2B\frac{3n^2}{B}。应用 Cache-Miss Lemma 即可证明。

2. 算法

2.1 矩阵乘法

2.1.1 朴素算法

image-20250928101200588

访问情况:

image-20250928101228704

有如下几种情况:(注意下面的分析的前提是,matrix 都是阵,cache 是的)

  • n>cMBn>c\frac{M}{B}:由于 BB 的访问是沿着列的,且矩阵的一行放不到一个 cache line 中,并且整个缓存也放不下 nn 行,这就导致每访问一个元素就会带来一次缓存不命中。因此 Q(n)=O(n3)Q(n)=O(n^3)
  • cM12<n<cMBc'M^{\frac{1}{2}}<n<c\frac{M}{B}:此时矩阵的一行仍放不到一个 cache line,但是整个缓存能放得下 nn 行,所以会呈现一定的空间局部性。BB 每全局访问一遍都会带来 O(n2B)O(\frac{n^2}{B}) 个缓存不命中。我们会全局访问 nnBB,因此 Q(n)=O(n3B)Q(n)=O(\frac{n^3}{B})
  • n<cM12n<c'M^{\frac{1}{2}}:此时整个矩阵都能放入 cache 中,所以我们只需要承担 cold miss 即可,Q(n)=O(n2B)Q(n)=O(\frac{n^2}{B})

2.1.2 调整访问顺序

如果我们交换一下循环:

image-20250928102517206

image-20250928102615410

BB 的访问情况大概是按行访问,然后完整访问 nnBB,于是 Q(n)=nO(n2B)=O(n3B)Q(n)=n\cdot O(\frac{n^2}{B})=O(\frac{n^3}{B})

2.1.3 分块

我们可以将矩阵分成若干个 s×ss\times s 的小矩阵,然后进行运算。通过调整 s×ss\times s 使得小矩阵能够放入缓存中,也就是我们尽可能满足 s=O(M12)s=O(M^{\frac{1}{2}})

image-20250928104253505

image-20250928104317126

然后一个小矩阵会带来 O(s2B)O(\frac{s^2}{B}) 个缓存不命中。我们要进行 O((ns)3)O((\frac{n}{s})^3) 次对小矩阵的操作,于是 Q(n)=O(n3BM12)Q(n)=O(\frac{n^3}{BM^{\frac{1}{2}}})

当然我们还可以对小矩阵再继续分块,使得能够利用二级缓存。以此类推,还能继续分块...

这个方法有个很大的缺点,就是我们需要对分块大小进行调整。但是受到外部因素的影响,我们很难非常精准地调整好。

2.1.4 分治

和分块差不多:

image-20250928105721115

image-20250928105840470

缓存不命中次数:

Q(n)={O(n2B)n2<cM, n 比较小,足以存放在缓存中8Q(n2)+O(1) otherwise Q(n) = \begin{cases} O(\frac{n^2}{B}) & n^2 < cM, \space n \text{ 比较小,足以存放在缓存中} \\ 8Q(\frac{n}{2}) + O(1) & \text{ otherwise } \end{cases}

其中 otherwise 情况中的 O(1)O(1) 是因为,我们调用函数进行分治的过程中可能会发生缓存不命中的情况,但是不命中次数是和数据规模无关的。

我们使用 recursion tree 的方法来分析复杂度:

image-20250928110746049

我们的算法仅仅是不断地把问题的规模缩小,那么缩小到一定程度后就很容易吃到缓存了。这样做的话我们无需知道缓存的大小是多少,是否有多级缓存等,问题被缩小到一定程度之后自然就会利用好缓存了。我们将这种算法称为 cache-oblivious 算法。

2.2 二维 DP

2.2.1 Cache oblivious

比如现在有如下的 DP 转移方程:

dpi,j=f(dpi1,j1,dpi1,j,dpi1,j+1)dp_{i,j}=f(dp_{i-1,j-1},dp_{i-1,j},dp_{i-1,j+1})

我们令 iT,jNi\le T,j\le N

如果直接朴素的去做的话,那么缓存不命中次数大概是 Q(n)=O(NTB)Q(n)=O(\frac{NT}{B})。现在我们考虑将朴素算法转为缓存高效的算法。首先是我们可以划分出如图所示的一个梯形区域:

image-20250929085532696

这个梯形区域中,如果最底下的点已经被计算完成的话,那么剩余的其他点也都可以被计算出来,且不依赖于梯形区域以外的点。

然后我们有如下两种的划分方案,如果这个梯形太“高”的话,或者说是 width<2height\text{width} < 2\cdot \text{height} ,那么我们可以直接从中间分开:

image-20250929085859454

如果这个梯形太“宽”,或者说是 width2height\text{width} \ge 2\cdot \text{height},那么我们就可以从中间用一条斜率为 1-1 的线去分成两半:

image-20250929090253246

当我们划分到区域高度为 1 的时候就到达 base case 了,不用再划分。

总之我们就是不断地去划分这些区域,使得这些区域变得越来越小。这一部分我们类似上面的分治的方法来实现。同样地,当区域小到一定程度的时候整个区域都能被完整地放入缓存中,从而降低了 cache miss 数量。

image-20250929090648521

使用 recursion tree 来分析缓存不命中的复杂度:

image-20250929090823392

叶子节点大概是 O(hw)O(hw) 个点,其中由于我们划分方法的特性,h=O(w)h=O(w)

每个叶子节点大概会造成 O(wB)O(\frac{w}{B}) 个 misses(因为 base case 只有一层的点需要求),其中 w=O(M)w=O(M),因为到叶子节点说明能被放入缓存中。

叶子节点点数总和应为总点数,因此我们有 O(NThw)O(\frac{NT}{hw}) 个叶子节点。

内部节点也是 O(NThw)O(\frac{NT}{hw}) 的,但每个内部节点只贡献 O(1)O(1) 的缓存不命中,总起来不影响渐进复杂度。

于是 Q(n)=O(NThw)O(wB)=O(NTM2)O(MB)=O(NTMB)Q(n)=O(\frac{NT}{hw})\cdot O(\frac{w}{B})=O(\frac{NT}{M^2})\cdot O(\frac{M}{B})=O(\frac{NT}{MB})

我们大幅度降低了缓存不命中次数,但是实际性能提升不明显。一个主要原因就是,我们这个问题的访问模式比较有规律,硬件可能会启用 prefetch 机制。不过我们这个算法又并非是朴素的做法,prefetch 会不准确,导致性能下降。

2.2.2 并行

上述算法并不能直接用于并行,因为划分出来的两个部分存在依赖关系。为了使用并行优化,我们可以像下面这样多划分一部分:

image-20250929093912876

image-20250929093931406

同色部分可以并行执行。

2.3 归并排序

朴素归并排序的缓存不命中情况:

Q(n)={O(nB) if ncM , constant c12Q(n2)+O(nB) otherwise Q(n) = \begin{cases} O(\frac{n}{B}) & \text{ if } n\le cM \text{ , constant } c \le 1 \\ 2Q(\frac{n}{2})+O(\frac{n}{B}) & \text{ otherwise } \end{cases}

缓存不命中复杂度分析

cache oblivious 的思路是,我们不再分成两部分再合并,而是分成 RR 部分合并。

R路归并排序

每个元素要经过 logn\log n 层来输出,所以合并的时间复杂度是 O(nlogn)O(n\log n) 的,于是我们有:

W(n)={O(1) if n=1RW(nR)+O(nlogR) otherwise W(n) = \begin{cases} O(1) &\text{ if } n = 1 \\ RW(\frac{n}{R}) + O(n\log R) &\text{ otherwise } \end{cases}

总时间复杂度分析

然后我们考虑一下缓存不命中次数。当 n<cMn<cM 时,如果取 R=cMBR=\frac{cM}{B} ,那么数据能够放入缓存中。于是我们有:

Q(n)={O(nB) if n<cMRQ(nR)+O(nB) otherwise Q(n) = \begin{cases} O(\frac{n}{B}) & \text{ if } n < cM \\ R\cdot Q(\frac{n}{R}) + O(\frac{n}{B}) & \text{ otherwise } \end{cases}

缓存不命中复杂度分析

一般我们有 nMn\gg M,于是:O(lognM)=O(logn)O(\log \frac{n}{M}) = O(\log n)

Qmultiway(n)=O(nBlogRnM)=O(nBlogMBnM)=O(nBlogMnM)=O(nlognMBlogM)=O(nlognBlogM)\begin{align} Q_{\text{multiway}}(n) &= O(\frac{n}{B}\log_R\frac{n}{M}) \\ &= O(\frac{n}{B}\log_{\frac{M}{B}}\frac{n}{M}) \\ &= O(\frac{n}{B}\log_M\frac{n}{M}) \\ &= O(\frac{n\log \frac{n}{M}}{B\log M}) \\ &= O(\frac{n\log n}{B\log M}) \end{align}

对比一下之前的:

Qbinary(n)=O(nBlognM)=O(nlognB)\begin{align} Q_{\text{binary}}(n) &= O(\frac{n}{B}\log \frac{n}{M}) \\ &= O(\frac{n\log n}{B}) \end{align}

相当于减少了 logM\log M 倍。当 L1-cache 大小为 2152^{15} 时,大概能有 1515 倍的加速。

2.4 BFS

在一个 CSR 形式的图上进行 BFS:

image-20251007130024981

其中箭头指向的那一行是最耗时的代码,因为对于每一条边,都会尝试判断其连向的点是否已经遍历过,这会带来 mm 次的随机访问,使得缓存基本上发挥不了什么作用。

一个优化方法是,我们开一个 bitset 来记录某个点是否遍历过。这样如果有 nn 个点,则只需要 nn 个位来存储。这样的 bitset 大小会比较小,容易放入到缓存中:

image-20251007130717314