“向量化”一词在不同的语境中是否意味着不同的事物？

Question

基于我之前所读到的，矢量化是一种称为SIMD的并行化形式。它允许处理器在数组上同时执行相同的指令(例如加法)。

然而，当我阅读关于朱莉娅和R的向量化性能的矢量编码与去卷积码的关系时，我感到很困惑。这篇文章声称，开发Julia代码(via循环)比Julia和R中的矢量化代码更快，因为：

这使一些不熟悉R的内部结构的人感到困惑，因此值得注意的是如何提高R代码的速度。性能改进的过程非常简单:先从去卷积R码开始，再用矢量化R码代替R码，最后用去vectorized代码实现这个矢量化R码。遗憾的是，这最后一步对于许多R用户来说是不可见的，因此他们认为向量化本身就是一种提高性能的机制。矢量化本身无助于使代码更快。使R中的矢量化有效的是，它提供了一种将计算迁移到C中的机制，在C中，隐藏层可以发挥其魔力。

它声称R将用R编写的矢量化代码转换成C中的去中心化代码。如果矢量化速度更快(作为并行化的一种形式)，为什么R要去代码，为什么这是一个加号？

Answer 1

R中的“向量化”，是R的解释器观点中的向量处理。以函数cumsum为例。在输入时，R解释器会看到向量x被传递到这个函数中。然而，工作被传递到C语言，R解释器无法分析/跟踪。当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ń滑雪展示了该语言的一个令人印象深刻的特性:循环融合。朱莉娅能够做到这一点，因为，正如他所指出的，“朱莉娅的向量化是在朱莉娅实现的”，而不是在语言之外。

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

下面是一个例子：计算R中两个瘦高矩阵间逐行点积的最快方法。rowSums(A * B)是R中的“向量化”，因为"*"和rowSums都是用C语言作为循环编码的。然而，R不能将它们合并成一个C循环，以避免将临时矩阵C = A * B生成到内存中。

另一个例子是R的回收规则或依赖于这种规则的任何计算。例如，当您通过a将标量A添加到矩阵A中时，真正发生的情况是，a首先被复制为与A具有相同维数的矩阵B，即B <- matrix(a, nrow(A), ncol(A))，然后计算两个矩阵之间的相加：A + B。显然，临时矩阵B的生成是不理想的，但遗憾的是，除非您为A + a编写自己的C函数并在R中调用它，否则就无法更好地实现它。

为了处理许多临时结果的内存效应，R有一种称为“垃圾收集”的复杂机制。这很有帮助，但是如果您在代码中的某个地方生成一些非常大的临时结果，那么内存仍然会爆炸。一个很好的例子是函数outer。我用这个函数写了许多答案，但它特别是记忆-不友好。

当我开始讨论“向量化”的副作用时，我可能在这个编辑中偏离了主题。小心使用。

• 将内存使用放在心上；可能会有一个更高效的向量化实现。例如，正如在两个矩阵之间的逐行点积上的链接线程中所提到的，c(crossprod(x, y))优于sum(x * y)。
• 准备使用已编译代码的CRAN R包。如果您发现R中现有的矢量化函数仅用于执行您的任务，请探索CRAN中可能完成任务的R包。您可以使用堆栈溢出的编码瓶颈提出一个问题，有人可能会指出正确的包中的正确函数。
• 很高兴自己编写自己的编译代码。

Answer 2

我认为值得注意的是，您所指的帖子并没有涵盖朱莉娅目前矢量化的所有功能。

重要的是，Julia中的向量化是在Julia实现的，而不是R，后者是在语言之外实现的。这一点在这篇文章中有详细的解释：https://julialang.org/blog/2017/01/moredots。

朱莉娅可以将任何广播运算序列融合成一个循环的结果。在提供矢量化的其他语言中，这种融合只有在显式实现的情况下才有可能。

总结如下：

1. 在Julia中，您可以预期向量化代码的速度与循环一样快。
2. 如果您执行一系列向量化操作，那么一般来说，您可以期望Julia比R更快，因为它可以避免分配计算的中间结果。

编辑：

在李哲源的注释之后，有一个例子显示，如果您想通过1增加向量x的所有元素，朱莉娅可以避免任何分配。

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

我们可以看到，它不仅节省内存，而且导致更快的代码执行。