为什么花式作业比花式查找慢?

这是一个复杂的问题,与现代处理器的架构特性和您的直觉密切相关 随机读取比随机写入慢,因为CPU必须等待读取数据 未经验证(大多数情况下)。有几个原因我会详细说明。

1. 现代处理器非常有效地隐藏读取延迟

2. 而内存写入比内存读取更昂贵

3. 尤其是在多核环境中

原因1现代处理器能够有效地隐藏读取延迟。

现代的 superscalar 可以同时执行多条指令,改变指令执行顺序( out of order execution ). 虽然这些特性的第一个原因是为了增加教学思想, 最有趣的结果之一是处理器隐藏内存写入延迟(或复杂的运算符、分支等)的能力。

为了解释这一点,让我们考虑一个将数组复制到另一个数组中的简单代码。

for i in a:
    c[i] = b[i]

一个经过编译的,由处理器执行的代码将以某种方式类似于

#1. (iteration 1) c[0] = b[0]
1a. read memory at b[0] and store result in register c0
1b. write register c0 at memory address c[0]
#2. (iteration 2) c[1] = b[1]
2a. read memory at b[1] and store result in register c1
2b. write register c1 at memory address c[1]
#1. (iteration 2) c[2] = b[2]
3a. read memory at b[2] and store result in register c2
3b. write register c2 at memory address c[2]
# etc

(这过于简单化,实际的代码更复杂,必须处理循环管理、地址计算等,但这种简单化的模型目前已经足够了)。

如问题中所述,对于读取,处理器必须等待实际数据。实际上,1b需要1a获取的数据,并且只要1a没有完成就不能执行。这种约束称为 附属国 我们可以说1b依赖于1a,依赖是现代处理器中的一个主要概念。依赖项表示算法(如我写b到c),必须绝对遵守。但是,如果指令之间不存在依赖关系,处理器将尝试执行其他挂起的指令,以保持操作管道始终处于活动状态。只要遵守依赖关系(类似于as-if规则),这可能导致执行无序。

对于考虑过的代码,这里 不 高级指令2之间的依赖关系。和1。(或在asm说明2a和2b与先前的说明之间)。实际上,最终的结果甚至是相同的是2。在1.之前执行,处理器将尝试在1a和1b完成之前执行2a和2b。2a和2b之间仍然存在依赖关系,但两者都可以发出。同样地,对于3a.和3b.等。这是一个强大的手段 隐藏内存延迟 . 如果出于某种原因2,3。和4。可以在1之前终止。加载它的数据,你甚至可能根本没有注意到任何减速。

此指令级并行性由处理器中的一组“队列”管理。

• 保留站RS中等待的指令队列(在最近的奔腾中键入128μ指令)。只要指令所需的资源可用(例如指令1b的寄存器c1的值),指令就可以执行。

• 一级缓存之前内存顺序缓冲区MOB中挂起的内存访问队列。这需要处理内存别名,并确保在同一地址(typ)的内存写入或加载中的顺序性。64次装载,32次储存)

• 由于类似的原因,在写回时执行顺序性的队列会导致寄存器(重新排序缓冲区或抢夺168个条目)。

• 以及在指令获取、μops生成、缓存中的写缓冲区和未命中缓冲区等方面的一些其他队列

在前一个程序的一个执行点上,会有许多等待存储的RS指令、MOB中的几个加载指令和ROB中等待失效的指令。

一旦数据变为可用(例如,读取终止),依赖于指令就可以执行并释放队列中的位置。但是,如果没有发生终止,并且其中一个队列已满,则与此队列关联的功能单元将暂停(如果处理器缺少寄存器名,则在发出指令时也可能发生这种情况)。暂停是造成性能损失的原因,为了避免这种损失,必须限制队列填充。

这解释了线性内存访问和随机内存访问之间的区别。
在线性访问中,1/由于更好的空间位置性和缓存可以用常规模式预取访问以进一步减少未命中次数,2/每当读取终止时,它将涉及一个完整的缓存线,并可以释放几个限制指令队列填充的挂起加载指令。这样,处理器就永远处于忙碌状态,并且内存延迟被隐藏。
对于随机访问,未命中的次数会更高,并且在数据到达时只能为单个加载提供服务。因此,指令队列将迅速饱和,处理器暂停和内存延迟不能再通过执行其他指令来隐藏。

处理器架构必须在吞吐量方面保持平衡,以避免队列饱和和暂停。实际上,在处理器的某个执行阶段通常有几十条指令,而全局吞吐量(即存储器(或功能单元)为指令请求提供服务的能力)是决定性能的主要因素。事实上,这些挂起的指令中有一些正在等待内存值,这对。。。

…除非你有很长的依赖链。

当指令必须等待前一条指令完成时,就存在依赖关系。使用读取结果是一种依赖关系。当涉及到依赖链时,依赖关系可能是一个问题。

例如,考虑一下 for i in range(1,100000): s += a[i] . 所有的内存读取都是独立的,但是在 s . 在前一个添加终止之前,不能进行任何添加。这些依赖关系将使保留站快速填满,并在管道中创建暂停。

但是读取很少涉及依赖链。仍然可以想象病理代码的所有读取都依赖于前一个(例如 for i in range(1,100000): s = a[s] ),但它们在实际代码中并不常见。这个问题来自于依赖链,而不是它是一个read;这种情况与计算绑定依赖代码类似(甚至可能更糟),比如 for i in range(1,100000): x = 1.0/x+1.0 .

因此,除了在某些情况下,计算时间与吞吐量的关系比与读取依赖性的关系更大,这是因为超标量失效或顺序执行隐藏了延迟。对于吞吐量来说,写操作比读操作更糟糕。

原因2:内存写入(特别是随机写入)比内存读取更昂贵

这与 caches 表现。缓存是存储一部分内存的快速内存(称为 线 )由处理器处理。缓存线目前为64字节,允许利用内存引用的空间位置:存储一行后,该行中的所有数据立即可用。重要的是 缓存和内存之间的所有传输都是行 .

当处理器对数据执行读取时,缓存检查数据所属的行是否在缓存中。否则,该行将从内存中提取,存储在缓存中,并将所需数据发送回处理器。

当处理器将数据写入内存时,缓存也会检查是否存在行。如果该行不存在,则缓存无法将其数据发送到内存(因为 全部的 传输是基于行的),并执行以下步骤:

1. cache从内存中获取行并将其写入缓存行。
2. 数据写入缓存,整行标记为已修改(脏)
3. 当缓存中的某一行被抑制时,它将检查modified标志,如果该行已被修改,则将其写回内存(写回缓存)

因此, 每次内存写入之前都必须有一个内存读取 在缓存中获取行。这增加了一个额外的操作,但对于线性写入来说并不昂贵。对于第一个写入的字,将有一个缓存未命中和一个内存读取,但是连续的写入只涉及缓存和被命中。

但随机写入的情况则大不相同。如果未命中的次数很重要,则每一次缓存未命中都意味着在从缓存中弹出行之前只读取少量的写操作,这将显著增加写入成本。如果一行在一次写入之后被弹出,我们甚至可以认为一次写入的时间开销是一次读取的两倍。

需要注意的是,增加内存访问(读或写)的数量往往会使内存访问路径饱和,并全局地减慢处理器和内存之间的所有传输。

无论哪种情况,写总是比读更贵。多核增强了这一点。

原因3:随机写入在多核中创建缓存未命中

不确定这是否真的适用于问题的情况。虽然numpy BLAS例程是多线程的,但我不认为基本数组拷贝是多线程的。但这是密切相关的,也是写作成本更高的另一个原因。

多核的问题是确保 cache coherence 以这样的方式,由多个处理器共享的数据在每个核心的缓存中被正确地更新。这是通过一个协议来完成的,比如 MESI 在写入缓存线之前更新缓存线,并使其他缓存副本失效(为所有权读取)。

虽然这些数据实际上都不是在问题的核心(或并行版本)之间共享的,但请注意,该协议适用于 缓存线 . 无论何时修改缓存线,都会从保存最新副本、本地更新的缓存中复制该缓存线,并且所有其他副本都将失效。即使核心正在访问缓存线的不同部分。这种情况叫做 false sharing 这是多核程序设计中的一个重要问题。

对于随机写入问题,缓存线是64字节,可以容纳8个int64,如果计算机有8个内核,每个内核平均处理2个值。因此,有一个重要的错误共享会减慢写入速度。

我们做了一些绩效评估。它是在C语言中执行的,以便包括对并行化影响的评估。我们比较了5个 处理大小为N的int64数组的函数。

1. 只是一份b到c的副本( c[i] = b[i] )(由编译器使用 memcpy() )

2. 用线性索引复制 c[i] = b[d[i]] 哪里 d[i]==i ( read_linear )

3. 用随机索引复制 c[i] = b[a[i]] 哪里 a 是随机的 0..N-1的置换( read_random 相当于 fwd 在最初的问题中)

4. 线性写入 c[d[i]] = b[i] 哪里 d[i]==i ( write_linear )

5. 随机写入 c[a[i]] = b[i] 具有 一 随机的 0..N-1的置换( write_random 相当于 inv 在这个问题上)

代码已用 gcc -O3 -funroll-loops -march=native -malign-double 在 天湖处理器。性能用 _rdtsc() 而且是 每次迭代的周期数。函数执行了几次(1000-20000取决于数组大小),执行了10个实验,并保留了最小的时间。

阵列大小从4000到1200000不等。所有代码都是用openmp的顺序和并行版本进行测量的。

这是结果的图表。函数有不同的颜色,顺序版本用粗线条表示,并行版本用细线条表示。

直接拷贝(显然)是最快的,由gcc实现 高度优化的 memcpy() . 这是一种利用内存估计数据吞吐量的方法。它的范围从小矩阵的每次迭代0.8个周期(CPI)到大矩阵的每次迭代2.0cpi。

读取线性性能大约是MEMCPY的两倍,但是有2个读和写,VS 1。 直接拷贝的读写。更多的索引添加了一些依赖项。最小值为1.56 CPI,最大值为3.8 CPI。写线性的稍长(5-10%)。

使用随机索引进行读写是原始问题的目的,需要更长的注释。这是结果。

size    4000    6000    9000    13496   20240   30360   45536   68304   102456  153680  230520  345776  518664  777992  1166984
rd-rand 1.86821 2.52813 2.90533 3.50055 4.69627 5.10521 5.07396 5.57629 6.13607 7.02747 7.80836 10.9471 15.2258 18.5524 21.3811
wr-rand 7.07295 7.21101 7.92307 7.40394 8.92114 9.55323 9.14714 8.94196 8.94335 9.37448 9.60265 11.7665 15.8043 19.1617 22.6785

• 小值(<10k):一级缓存是32k,可以容纳4k的uint64数组。注意,由于索引的随机性,在大约1/8次迭代之后,L1缓存将完全填充随机索引数组的值(因为缓存线是64字节的,可以容纳8个数组元素)。访问其他线性阵列我们将快速生成许多一级未命中,我们必须使用二级缓存。一级缓存访问是5个周期,但它是流水线式的,每个周期可以提供几个值。二级访问时间较长,需要12个周期。对于随机读写,未命中的数量是相似的,但是我们看到,当数组大小很小时,我们完全支付了写操作所需的双重访问费。

• 中等值(10k-100k):二级缓存为256k,可容纳32k int64数组。之后,我们需要转到L3缓存(12Mo)。随着大小的增加,L1和L2中的未命中数和计算时间相应地增加。这两种算法都有相似数量的未命中,主要是由于随机读写(其他访问是线性的,缓存可以非常有效地预取)。我们检索随机读写之间的因子2,在B.M.答案中已经注意到了。写作的双重成本可以部分解释这一点。

• 大值(>100k):方法之间的差异逐渐减小。对于这些大小,很大一部分信息存储在三级缓存中。L3大小足以容纳1.5米的完整阵列,而且管线不太可能被弹出。因此,对于写操作,在初始读取之后,可以在不弹出行的情况下完成大量的写操作,并且减少了写操作与读取操作的相对成本。对于这些大尺寸,还需要考虑许多其他因素。例如,缓存只能提供有限数量的未命中(typ。16) 当未命中的次数很大时,这可能是限制因素。

一个字的并行omp版本的随机读写。除了较小的大小之外,让随机索引数组分布在几个缓存上可能不是一个优势,它们在系统上要快~两倍。对于大容量,我们清楚地看到,由于错误共享,随机读写之间的差距增大。

几乎不可能用当前计算机体系结构的复杂性来进行定量预测,即使对于简单的代码,甚至对行为的定性解释也是困难的,并且必须考虑到许多因素。正如在其他答案中提到的,与python相关的软件方面也会产生影响。但是,虽然在某些情况下可能会发生这种情况,但大多数情况下,由于数据依赖性,人们不能认为读操作更昂贵。

• 首先反驳你的直觉: fwd 跳动 inv 即使没有新墨西哥主义。

是这样的 麻木 版本:

import numba

@numba.njit
def fwd_numba(a,b,c):
    for i in range(N):
        c[a[i]]=b[i]

@numba.njit
def inv_numba(a,b,c):
    for i in range(N):
        c[i]=b[a[i]]

N=10000的计时:

%timeit fwd()
%timeit inv()
%timeit fwd_numba(a,b,c)
%timeit inv_numba(a,b,c)
62.6 µs ± 3.84 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
144 µs ± 2 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
16.6 µs ± 1.52 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
34.9 µs ± 1.57 µs per loop (mean ± std. dev. of 7 runs, 100000 loops each)

• 其次,Numpy必须处理对齐和(缓存)局部性的可怕问题。

它本质上是来自 布拉斯/阿特拉斯/MKL 为之调整。 花式索引是一个很好的高级工具,但对于这些问题的异端者,在低层次上没有直接的概念。

• 第三, numpy dev docs : 细节花哨索引。特别地:

除非在项目获取过程中只有一个索引数组,否则 事先检查指标的有效性。否则会被处理 在内部循环中进行优化。

我们在这件事上。我认为这可以解释差异,以及为什么set比get慢。

这也解释了为什么手工制作 麻木 通常更快:它不检查任何内容,并在不一致的索引上崩溃。

你的两个裸体片段 b[a] 和 c[a] = b 似乎是测量无序/线性读/写速度的合理的启发式方法,我将在下面的第一节中通过查看底层NumPy代码进行论证。

关于哪一个应该更快的问题,似乎有理由认为,随机读线性写通常可以胜出(正如基准显示的那样),但速度的差异可能会受到随机索引的“随机”程度的影响,并且一个或多个:

• CPU缓存读取/更新策略( write-back vs. write-through 等等)。
• CPU如何选择(重新)排序它需要执行的指令(流水线)。
• CPU识别内存访问模式并预取数据。
• 缓存逐出逻辑。

即使假设有哪些策略,这些效果也很难进行建模和分析,因此我不确定是否可以得到适用于所有处理器的通用答案(尽管我不是硬件专家)。

不过,在下面的第二部分中,我将尝试解释为什么在某些假设下,无序读线性写显然更快。

“琐碎”花式索引

本节的目的是浏览NumPy源代码,以确定是否对计时有任何明显的解释,并尽可能清楚地了解当 A[B] 或 A[B] = C 被处决。

迭代程序支持虚拟索引 getitem 和 setitem 这个问题的操作是“ trivial ":

• B 是一个单跨的单索引数组
• A 和 乙 具有相同的内存顺序(两个C-连续或两个Fortran连续)

此外,在我们的案例中 一个 和 乙 是 Uint Aligned :

跨步复制代码:这里使用“uint alignment”。如果数组的itemsize[N]等于1、2、4、8或16个字节,并且数组是uint对齐的,那么[而不是使用缓冲]numpy将执行 *(uintN*)dst) = *(uintN*)src) 以获得适当的N。否则,通过 memcpy(dst, src, N) .

这里的要点是使用内部缓冲区来确保避免对齐。使用 *(单位*)dst)=*(单位*)src) 与“将偏移量src的X字节放入偏移量dst的X字节”一样直接。

编译器可能会很简单地将其转换为 mov 指令(例如在x86上)或类似的。

执行项目获取和设置的核心低级代码在函数中 mapiter_trivial_get 和 mapiter_trivial_set . 这些功能产生于 lowlevel_strided_loops.c.src ,其中模板和宏使阅读变得有些困难(这是一个感谢高级语言的机会)。

坚持不懈,我们最终会发现getitem和setitem之间没有什么区别。这是主循环的简化版本,以供说明。宏行确定运行的是getitem还是setitem:

    while (itersize--) {
        char * self_ptr;
        npy_intp indval = *((npy_intp*)ind_ptr);

#if @isget@
        if (check_and_adjust_index(&indval, fancy_dim, 0, _save) < 0 ) {
            return -1;
        }
#else
        if (indval < 0) {
            indval += fancy_dim;
        }
#endif

        self_ptr = base_ptr + indval * self_stride; /* offset into array being indexed */

#if @isget@
        *(npy_uint64 *)result_ptr = *(npy_uint64 *)self_ptr;
#else
        *(npy_uint64 *)self_ptr = *(npy_uint64 *)result_ptr;
#endif

        ind_ptr += ind_stride;         /* move to next item of index array */
        result_ptr += result_stride;   /* move to next item of result array */

如我们所料,这仅仅相当于一些算法,以获得数组的正确偏移量,然后将字节从一个内存位置复制到另一个内存位置。

setitem的额外索引检查

值得一提的是,对于setitem,索引的有效性(无论它们是否都是目标数组的inbounds)是 checked before copying 开始(通过 check_and_adjust_index ),也用相应的正指数代替了负指数。

在上面的片段中,您可以看到 检查并调整索引 在主循环中调用getitem,同时对setitem进行更简单(可能是多余的)的负索引检查。

这种额外的初步检查可能会对setitem的速度产生微小但负面的影响( A[B]=C ).

缓存未命中

因为这两个代码片段的代码非常相似,所以怀疑的焦点在于CPU以及它如何处理对底层内存数组的访问。

CPU缓存最近访问过的小内存块(缓存线),因为它可能很快需要再次访问该内存区域。

对于上下文,缓存线通常是64字节。我老化的笔记本电脑CPU上的L1(最快)数据缓存为32KB(足以容纳阵列中的500个int64值,但请记住,在执行NumPy代码段时,CPU将执行其他需要内存的操作):

$ cat /sys/devices/system/cpu/cpu0/cache/index0/coherency_line_size
64
$ cat /sys/devices/system/cpu/cpu0/cache/index0/size
32K

正如您可能已经知道的,对于读/写内存,顺序缓存工作得很好,因为64字节的内存块是根据需要获取的,并且存储在靠近CPU的地方。对该内存块的重复访问比从RAM(或更高级别的缓存)中获取更快。实际上,CPU甚至可以在程序请求下一个缓存线之前抢先获取它。

另一方面,随机访问内存可能会导致频繁的缓存未命中。这里,具有所需地址的内存区域不在CPU附近的快速缓存中,而是必须从更高级别的缓存(较慢)或实际内存(较慢得多)访问。

那么,对于CPU来说,哪种处理速度更快:频繁的数据读取未命中,还是数据写入未命中?

假设CPU的写策略是写回,这意味着修改后的内存被写回缓存。缓存被标记为正在修改(或“脏的”),只有当缓存线被逐出时,更改才会被写回主内存(CPU仍然可以从脏的缓存线读取)。

如果我们正在写入一个大数组中的随机点,则预期CPU缓存中的许多缓存线将变脏。当每个内存都被逐出时,需要对主内存进行直写,如果缓存已满,这种情况可能会经常发生。

但是,当按顺序写入数据并随机读取数据时,这种直写应该不那么频繁,因为我们希望更少的缓存线变脏,数据回主内存或较慢的缓存不那么频繁。

如前所述,这是一个简化的模型,可能还有许多其他因素影响CPU的性能。比我更专业的人很可能会改进这种模式。

你的职能 fwd 没有触及全局变量 c . 你没告诉我 global c (仅在 setup ),所以它有自己的局部变量,并使用 STORE_FAST 在cpython中:

>>> import dis
>>> def fwd():
...     c = b[a]
...
>>> dis.dis(fwd)
  2           0 LOAD_GLOBAL              0 (b)
              3 LOAD_GLOBAL              1 (a)
              6 BINARY_SUBSCR
              7 STORE_FAST               0 (c)
             10 LOAD_CONST               0 (None)
             13 RETURN_VALUE

现在,让我们在全球范围内尝试一下:

>>> def fwd2():
...     global c
...     c = b[a]
...
>>> dis.dis(fwd2)
  3           0 LOAD_GLOBAL              0 (b)
              3 LOAD_GLOBAL              1 (a)
              6 BINARY_SUBSCR
              7 STORE_GLOBAL             2 (c)
             10 LOAD_CONST               0 (None)
             13 RETURN_VALUE

即使如此,它可能在时间上与 inv 调用的函数 setitem 为了一个全球性的。

不管怎样,如果你想写的话 进入之内 c类 ,你需要像 c[:] = b[a] 或 c.fill(b[a]) . 赋值将变量(名称)替换为右侧的对象,因此 c 可能会被取消分配而不是新的 b[a] ,而这种内存洗牌可能代价高昂。

至于效果,我认为你想衡量,基本上是正向还是反向排列更昂贵,这将是高度依赖缓存。前向排列(以线性读取的随机顺序索引存储)原则上可能更快,因为它可以使用写掩蔽,并且永远不会获取新数组,前提是缓存系统足够聪明,可以在写缓冲区中保留字节掩蔽。如果数组足够大,则在执行随机读取时向后运行缓存冲突的高风险。

这是我最初的印象;正如你所说,结果恰恰相反。这可能是由于缓存实现没有大的写缓冲区或无法利用小的写操作造成的。如果缓存外访问仍然需要相同的内存总线时间,则读取访问将有机会加载在需要之前不会从缓存中删除的数据。使用多路缓存时,部分写入的缓存线也有可能不被选择用于排除;只有脏缓存线需要内存总线时间来丢弃。使用其他知识编写的低级程序(例如,置换是完整的和不重叠的)可以使用诸如非时态SSE写入之类的提示来改进行为。

下面的实验证实了随机写入比随机读取快。对于较小的数据(当它完全适合缓存时),随机写入代码比随机读取代码慢(可能是由于 numpy ),但随着数据大小的增长,执行时间的最初1.7倍差异几乎完全消除(但是,在 numba 这一趋势最终出现了奇怪的逆转)。

$ cat test.py 
import numpy as np
from timeit import timeit
import numba

def fwd(a,b,c):
    c = b[a]

def inv(a,b,c):
    c[a] = b

@numba.njit
def fwd_numba(a,b,c):
    for i,j in enumerate(a):
        c[i] = b[j]

@numba.njit
def inv_numba(a,b,c):
    for i,j in enumerate(a):
        c[j] = b[i]


for p in range(4, 8):
    N = 10**p
    n = 10**(9-p)
    a = np.random.permutation(N)
    b = np.random.random(N)
    c = np.empty_like(b)
    print('---- N = %d ----' % N)
    for f in 'fwd', 'fwd_numba', 'inv', 'inv_numba':
        print(f, timeit(f+'(a,b,c)', number=n, globals=globals()))

$ python test.py 
---- N = 10000 ----
fwd 1.1199337750003906
fwd_numba 0.9052993479999714
inv 1.929507338001713
inv_numba 1.5510062070025015
---- N = 100000 ----
fwd 1.8672701190007501
fwd_numba 1.5000483989970235
inv 2.509873716000584
inv_numba 2.0653326050014584
---- N = 1000000 ----
fwd 7.639554155000951
fwd_numba 5.673054756000056
inv 7.685382894000213
inv_numba 5.439735023999674
---- N = 10000000 ----
fwd 15.065879136000149
fwd_numba 12.68919651500255
inv 15.433822674000112
inv_numba 14.862108078999881