Julia：三维数组（性能）

Question

想到 Julia 的 performance tips 我还没有找到任何关于如何使用三维数组加速代码的建议。

据我了解，d-element Array{Array{Float64,2},1} 在d（第三维度）较小时表现最佳。但是，我不确定d 很大时是否是这种情况。

对于 Julia 有没有关于这个主题的教程？

---

示例 1a (d=50)

x = [zeros(100, 10) for d=1:50];

@time for d=1:50
    x[d] = rand(100,10);
end

0.000100 seconds (50 allocations: 396.875 KB)

示例 1b (d=50)

y=zeros(100, 10, 50);

@time for d=1:50
    y[:,:,d] = rand(100,10);
end

0.000257 seconds (200 allocations: 400.781 KB)

---

示例 2a (d=50000)

x = [zeros(100, 10) for d=1:50000];

@time for d=1:50000
    x[d] = rand(100,10);
end

0.410813 seconds (99.49 k allocations: 388.328 MB, 81.88% gc time)

示例 2b (d=50000)

y=zeros(100, 10, 50000);

@time for d=1:50000
    y[:,:,d] = rand(100,10);
end

0.185929 seconds (298.98 k allocations: 392.898 MB, 6.83% gc time)

Answer 1

据我了解，当 d（第三维）较小时，d-element Array{Array{Float64,2},1} 表现最佳。但是，我不确定当 d 很大时是否是这种情况。

不，更重要的是你如何使用它。 A = Array{Array{Float64,2},1} 是一个矩阵指针数组。数组的值是指针或引用。因此A[i] 返回一个引用，即它很便宜。 A2 = Array{Float64,3} 是一个连续的浮点数组。它实际上只是在线性内存块上的索引设置（并且有一个线性索引A2[i]，它使用该线性形式贯穿整个事物）。

后者有一些优势，因为它是连续的。没有间接性，因此循环所有 A2s 值会更快。 A 必须服从两个指针才能获得一个值，因此如果您不知道只服从每个内部矩阵一次，那么简单的 3D 循环会更慢。此外，您可以通过@view A2[:,:,1] 等获取矩阵的视图，但您必须注意A2[:,:,1] 本身会制作矩阵的副本。 A[1] 是一个自然视图，因为它返回对 matirx 的引用，如果你想复制你必须明确地做 copy(A[1])。因为A 只是一个线性指针数组，所以push! 在其上添加一个新矩阵很便宜，因为它只是增加一个相对较小的数组（并且push! 会自动摊销）以在末尾添加一个新指针（这这就是为什么像DifferentialEqautions.jl 这样的东西使用数组的数组而不是更传统的矩阵来构建时间序列）。

因此它们是具有不同优点和缺点的不同工具。

至于你的时间安排，你正在做两件不同的事情。 x[d] = rand(100,10) 正在创建一个新矩阵并将其引用添加到 x。 y[:,:,d] = rand(100,10) 正在创建一个新矩阵并循环遍历 y 的值以更改 y 的值。你可以看到为什么它更慢。但是您遗漏的是免分配案例。

function f2()
    y=zeros(100, 10, 50);

    @time for i in eachindex(y)
        y[i] = rand()
    end
    y
end

在较小的情况下，这与数组创建相匹配。在第一种情况下你不能天真地这样做，但正如我所说，如果你在做得很好后取消引用矩阵的指针：

function f()
    x = [zeros(100, 10) for d=1:5000];

    @time @inbounds for d=1:50
        xd = x[d]
        for i in eachindex(xd)
            xd[i] = rand()
        end
    end
    x
end

因此，在适当的情况下，数组的数组可以是很好的数据结构。创建库RecursiveArrayTools.jl 是为了更好地利用它。例如，A3 = VectorOfArrays(A) 通过将A[i,j,k] 延迟转换为A[k][i,j]，为A3 提供与A2 相同的索引结构。但是，它保留了A 的优点，但会自动确保以正确的方式广播，如f。像这样的另一个工具是ArrayPartition，它允许以广播性能的方式进行异构输入。

是的，它并不总是正确的工具，但如果使用得当，这些异构和递归数组是很棒的工具。