Julia 是否需要矢量化来加快计算速度？

Question

在 Python 中，通常建议将代码向量化以加快计算速度。例如，如果你想计算两个向量的内积，比如a 和b，通常是

代码 A

c = np.dot(a, b)

优于

代码 B

c = 0
for i in range(len(a)):
    c += a[i] * b[i]

但在 Julia 中，有时矢量化似乎没那么有用。我将'* 和dot 视为矢量化版本，将显式for 循环视为非矢量化版本，并得到以下结果。

using Random
using LinearAlgebra

len = 1000000
rng1 = MersenneTwister(1234)
a = rand(rng1, len)
rng2 = MersenneTwister(12345)
b = rand(rng2, len)


function inner_product(a, b)
    result = 0
    for i in 1: length(a) 
        result += a[i] * b[i] 
    end
    return result
end


@time a' * b
@time dot(a, b)
@time inner_product(a, b);

  0.013442 seconds (56.76 k allocations: 3.167 MiB)
  0.003484 seconds (106 allocations: 6.688 KiB)
  0.008489 seconds (17.52 k allocations: 976.752 KiB)

（我知道使用BenchmarkTools.jl 是衡量性能的更标准方法。）

从输出结果来看，dot 的运行速度比 for 快于 '*，这与假设的情况相矛盾。

所以我的问题是，

1. Julia 是否需要（或有时需要）矢量化来加快计算速度？
2. 如果是，那么何时使用矢量化以及哪种方法更好（考虑dot 和'*）？
3. 如果不是，那么在矢量化和非矢量化代码的机制方面，Julia 和 Python 有什么区别？

Answer 1

让我将我的实践经验作为评论添加（对于标准评论来说太长了）：

Julia 是否需要（或有时需要）矢量化来加快计算速度？

Julia 不需要像 Python 那样进行矢量化（请参阅 Przemysław 的答案），但在实践中，如果您有一个编写良好的矢量化函数（如 dot），那么尽可能使用它，但有时编写起来可能会很棘手自己作为高性能函数（人们可能花了几天时间优化dot，尤其是优化使用多线程）。

如果是，那么何时使用矢量化以及哪种方式更好用（考虑点和'*）？

当您使用矢量化代码时，这一切都取决于您要使用的函数的实现。在这种情况下dot(a, b) 和a' * b 与您在这种情况下编写@edit a' * b 给您的情况完全相同：

*(u::AdjointAbsVec{<:Number}, v::AbstractVector{<:Number}) = dot(u.parent, v)

你看是一样的。

如果不是，那么在矢量化和非矢量化代码的机制方面，Julia 和 Python 有什么区别？

Julia 是一种编译语言，而 Python 是一种解释语言。在某些情况下，Python 解释器可以提供快速的执行速度，但在其他情况下，它目前无法做到这一点（但这并不意味着将来不会改进）。特别是矢量化函数（如您的问题中的dot）很可能是用某种编译语言编写的，因此 Julia 和 Python 在典型情况下不会有太大差异，因为它们只是调用此编译函数。但是，当您使用循环（非向量化代码）时，当前 Python 会比 Julia 慢。

Answer 2

您没有正确地进行基准测试，并且您的功能的实现不是最理想的。

julia> using BenchmarkTools                   
                                              
julia> @btime $a' * $b                        
  429.400 μs (0 allocations: 0 bytes)         
249985.3680190253                             
                                              
julia> @btime dot($a,$b)                      
  426.299 μs (0 allocations: 0 bytes)         
249985.3680190253                             
                                              
julia> @btime inner_product($a, $b)           
  970.500 μs (0 allocations: 0 bytes)         
249985.36801903677                            

正确的实现方式：

function inner_product_correct(a, b)
    result = 0.0 #use the same type as elements in the args
    @simd for i in 1: length(a)
        @inbounds result += a[i] * b[i]
    end
    return result
end

julia> @btime inner_product_correct($a, $b)
  530.499 μs (0 allocations: 0 bytes)
249985.36801902478

仍然存在差异（但不那么重要），因为 dot 使用的是优化的多线程 BLAS 实现。你可以（按照 Bogumil 的评论集 OPENBLAS_NUM_THREADS=1 然后你会发现 BLAS 函数的时间将与 Julia 实现相同。

另请注意，使用浮点数在许多方面都很棘手：

julia> inner_product_correct(a, b)==dot(a,b)
false

julia> inner_product_correct(a, b) ≈ dot(a,b)
true 

最后，在 Julia 中，决定是使用矢量化还是自己编写循环取决于您的两倍 - 没有性能损失（只要您编写类型稳定的代码并在需要时使用 @simd 和 @inbounds）。但是，在您的代码中，您没有测试矢量化，而是将调用 BLAS 与自己编写循环进行比较。这是了解正在发生的事情的必读内容https://docs.julialang.org/en/v1/manual/performance-tips/