术语"矢量化"在不同的背景下是否意味着不同的东西？

Question

基于我之前读过的内容,矢量化是一种称为SIMD的并行化形式.它允许处理器同时在阵列上执行相同的指令(例如添加).

然而,在阅读关于Julia和R的矢量化性能的矢量化和非矢量化代码之间的关系时,我感到困惑.该帖子声称,Julia和R开发的Julia代码(通过循环)比矢量化代码更快,因为:

这使一些不熟悉R内部的人感到困惑.因此值得注意的是如何提高R代码的速度.性能改进的过程非常简单:首先从devectorized R代码开始,然后用向量化R代码替换它,然后最终在devectorized C代码中实现这个向量化R代码.遗憾的是,最后一步对于许多R用户是不可见的,因此他们认为向量化本身是提高性能的机制.矢量化本身无助于使代码更快.使R中的矢量化有效的原因在于它提供了一种将计算移动到C中的机制,其中一个隐藏的devectorization层可以发挥其神奇作用.

它声称R将用R编写的矢量化代码转换为C中的devectorized代码.如果矢量化更快(作为一种并行化形式),为什么R会驱动代码,为什么这是一个加号？

Answer 1

R中的"矢量化"是R解释器视图中的矢量处理.cumsum以此功能为例.在进入时,R解释器看到向量x被传递到此函数中.但是,然后将工作传递给R解释器无法分析/跟踪的C语言.当C正在做工作时,R只是在等待.当R的解释器恢复工作时,已经处理了一个向量.所以在R的观点中,它发出了一条指令,但处理了一个向量.这类似于SIMD的概念 - "单指令,多数据".

不只是cumsum采用向量并返回向量的函数在R中被视为"向量化",类似于sum采用向量并返回标量的函数也是"向量化".

简单地说:每当R为循环调用一些编译代码时,它就是"向量化".如果你想知道为什么这种"矢量化"是有用的,那是因为编译语言编写的循环比用解释语言编写的循环更快.C循环被转换为CPU可以理解的机器语言.但是,如果CPU想要执行R循环,则需要R的解释器帮助读取它,迭代迭代.这就像,如果你懂中文(最难的人类语言),你可以更快地回复向你说中文的人; 否则,你需要一名翻译人员用英语逐句翻译中文,然后用英语回答,翻译人员会逐句回复中文.沟通的有效性大大降低.

x <- runif(1e+7)

## R loop
system.time({
  sumx <- 0
  for (x0 in x) sumx <- sumx + x0
  sumx
  })
#   user  system elapsed 
#  1.388   0.000   1.347 

## C loop
system.time(sum(x))
#   user  system elapsed 
#  0.032   0.000   0.030 

请注意,R中的"矢量化"只是SIMD的一个类比,但不是真实的.真正的SIMD使用CPU的向量寄存器进行计算,因此通过数据并行是真正的并行计算.R不是可以编程CPU寄存器的语言; 你必须为此目的编写编译代码或汇编代码.

R的"向量化"并不关心用编译语言编写的循环是如何真正执行的; 毕竟这超出了R的翻译知识.关于这些编译的代码是否将使用SIMD执行,读取R在进行矢量化计算时是否利用SIMD？

更多关于R中的"矢量化"

我不是朱莉娅的用户,但BogumiłKamiński展示了该语言的一个令人印象深刻的特征:循环融合.朱莉娅可以做到这一点,因为正如他所指出的那样,"朱莉娅的矢量化是在朱莉娅实施的",而不是在语言之外.

这揭示了R的矢量化的缺点:速度通常以内存使用的代价为代价.我不是说朱莉娅不会有这个问题(因为我不使用它,我不知道),但对于R来说这绝对是正确的.

下面是一个示例:计算R中两个高高矩阵之间的行方点积的最快方法.rowSums(A * B)是R A"矢量化",因为这两个"*"和rowSums在C语言作为一个循环进行编码.但是,R不能将它们融合到单个C循环中,以避免将临时矩阵生成C = A * B到RAM中.

另一个例子是R的回收规则或依赖这种规则的任何计算.例如,当您添加一个标量a的矩阵A通过A + a,真正的情况是,a首先复制是一个矩阵B具有与相同尺寸A,即B <- matrix(a, nrow(A), ncol(A)),然后计算两个矩阵之间的加法:A + B.显然,临时矩阵的产生B是不需要的,但对不起,你不能做的更好,除非你写你自己的C函数A + a,这被描述为调用它R. "这样的融合才有可能在明确实施"中BogumiłKamiński的回答.

为了处理许多临时结果的记忆效应,R有一种称为"垃圾收集"的复杂机制.它有所帮助,但是如果你在代码中的某个地方产生一些非常大的临时结果,内存仍然会爆炸.一个很好的例子就是这个功能outer.我用这个函数写了很多答案,但它特别对内存不友好.

当我开始讨论"矢量化"的副作用时,我可能在这个编辑中偏离了主题.小心使用.

• 记住内存使用; 可能存在更高内存效率的矢量化实现.例如,正如在两个矩阵之间的行方点产品的链接线程中所提到的,c(crossprod(x, y))优于sum(x * y).
• 准备使用已编译代码的CRAN R包.如果您发现R中的现有矢量化函数仅限于执行您的任务,请探索可以执行此操作的可能R包的CRAN.您可以在Stack Overflow上询问有关编码瓶颈的问题,有人可能会在正确的包中指出您正确的功能.
• 很乐意编写自己编译的代码.

Answer 2

我认为值得注意的是,您所指的帖子并未涵盖Julia中所有当前矢量化的功能.

重要的是Julia中的矢量化是在Julia中实现的,而不是R,它是在语言之外实现的.这篇文章对此进行了详细解释:https://julialang.org/blog/2017/01/moredots.

Julia可以将任何广播操作序列融合到单个循环中的结果.在提供矢量化的其他语言中,只有在明确实现的情况下才可能进行融合.

综上所述:

1. 在Julia中,您可以期望矢量化代码与循环一样快.
2. 如果执行矢量化操作序列,那么通常可以预期Julia比R更快,因为它可以避免分配计算的中间结果.

编辑:

根据李哲源的评论,这里有一个例子,表明如果你想通过以下方式增加向量x的所有元素,Julia能够避免任何分配1:

julia> using BenchmarkTools

julia> x = rand(10^6);

julia> @benchmark ($x .+= 1)
BenchmarkTools.Trial:
  memory estimate:  0 bytes
  allocs estimate:  0
  --------------
  minimum time:     819.230 ?s (0.00% GC)
  median time:      890.610 ?s (0.00% GC)
  mean time:        929.659 ?s (0.00% GC)
  maximum time:     2.802 ms (0.00% GC)
  --------------
  samples:          5300
  evals/sample:     1

在代码中.+=执行加法($在表达式前面添加只需要进行基准测试,在普通代码中就是这样x .+= 1).我们看到没有完成内存分配.

如果我们将其与R中的可能实现进行比较:

> library(microbenchmark)
> x <- runif(10^6)
> microbenchmark(x <- x + 1)
Unit: milliseconds
       expr      min       lq     mean   median       uq      max neval
 x <- x + 1 2.205764 2.391911 3.999179 2.599051 5.061874 30.91569   100

我们可以看到它不仅可以节省内存,还可以加快代码的执行速度.