在Julia中有效地将标量添加到矩阵中

Question

我需要为一个巨大矩阵的所有元素添加一个标量.矩阵将尽可能大.在示例中,我将使用2 GiB的大小,但在我的实际计算中它将更大.

A = rand(2^14, 2^14)

如果我执行

A += 1.0

Julia分配额外的2 GiB内存.该操作大约需要1秒.我可以使用一个for循环:

for jj = 1:size(A, 2), ii = 1:size(A, 1)
  A[ii, jj] = A[ii, jj] + 1.0
end

这不会分配任何内存,但需要一分钟.这两种方法对我来说都不可行,因为第一种方法违反了内存约束,第二种方法明显效率低下.对于元素乘法,有scal!使用BLAS.有没有什么方法可以像乘法一样有效地执行加法scal!？

Answer 1

@ DSM的答案很好.不过,我还有很多事情需要解决.你的for循环很慢的原因是因为它A是一个非常量的全局变量,你的代码直接改变了那个全局变量.由于A是非常量的,因此代码必须防止A在循环执行期间在任何时刻变为具有不同类型的不同值的可能性.代码必须A在循环的每次迭代中查找类型和位置,并在表达式中动态调度方法调用A[ii, jj] = A[ii, jj] + 1.0- 这是一个调用getindex,+并且setindex!所有这些都依赖于静态未知类型A.只需在函数中执行此工作,即可立即获得更好的性能:

julia> A = rand(2^10, 2^10);

julia> @time for jj = 1:size(A, 2), ii = 1:size(A, 1)
           A[ii, jj] += 1
       end
elapsed time: 0.288340785 seconds (84048040 bytes allocated, 15.59% gc time)

julia> function inc!(A)
           for jj = 1:size(A, 2), ii = 1:size(A, 1)
               A[ii, jj] += 1
           end
       end
inc! (generic function with 1 method)

julia> @time inc!(A)
elapsed time: 0.006076414 seconds (171336 bytes allocated)

julia> @time inc!(A)
elapsed time: 0.000888457 seconds (80 bytes allocated)

避免像这样的非常数全局变量是本手册" 性能提示"部分的第一个建议.您可能也想仔细阅读本章的其余部分.

我们可以inc!使用@inbounds注释进一步提高函数的性能,以指示此代码不需要边界检查,并使用线性索引而不是二维索引:

julia> function inc!(A)
           @inbounds for i = 1:length(A)
               A[i] += 1
           end
       end
inc! (generic function with 1 method)

julia> @time inc!(A)
elapsed time: 0.000637934 seconds (80 bytes allocated)

大多数加速来自@inbounds注释而不是线性索引,尽管这确实提供了一点速度提升.该@inbounds注释应节制地使用,然而,只有其中一个是确定的,该索引不能出界和性能是至关重要的.正如您所看到的,现有的额外性能改进并非压倒性的.大部分好处来自于不直接改变全局变量.

Answer 2

你可以做一个就地广播操作:

julia> A = rand(2^14, 2^14); A[1:5, 1:5]
5x5 Array{Float64,2}:
 0.229662  0.680236    0.131202  0.111664   0.802698
 0.500575  0.580994    0.385844  0.983806   0.324382
 0.701694  0.577749    0.532591  0.0508955  0.94325 
 0.592929  0.00319653  0.759241  0.448704   0.706204
 0.867945  0.0413606   0.586151  0.82561    0.679233

julia> @time broadcast!(.+, A, A, 100);
elapsed time: 0.382669486 seconds (11490976 bytes allocated)

julia> A[1:5, 1:5]
5x5 Array{Float64,2}:
 100.23   100.68   100.131  100.112  100.803
 100.501  100.581  100.386  100.984  100.324
 100.702  100.578  100.533  100.051  100.943
 100.593  100.003  100.759  100.449  100.706
 100.868  100.041  100.586  100.826  100.679

它只使用了总共​​~2G的内存.