julia - 是否可以为矩阵分解预先分配数组？

Question

我的问题是而不是 F = svd(A) , 是否可以先为 SVD 结构分配适当的内存，然后执行 F .= svd(A) ?
我想到的是如下内容:

function main()
    F = Vector{SVD}(undef,10)
    # how to preallocate F?
    test(F)
end

function test(F::Vector{SVD})
    for i in 1:10
        F .= svd(rand(3,3))
    end
end

Answer 1

您的代码几乎可以正常工作。但你可能想要的是这样的:

using LinearAlgebra

function main()
    F = Vector{SVD}(undef, 10)
    test(F)
end

function test(F::Vector{SVD})
    for i in 1:10
        F[i] = svd(rand(3, 3))
    end
    return F
end

您在 for 中的行循环是这样的:

F .= svd(rand(3,3))

它对每个循环执行相同的操作，因为您没有索引到 F .特别是，此操作试图广播单个 SVD对象到 F 的所有字段中在循环的每次迭代中。 (并且该广播操作失败，因为默认情况下 struct s 被视为具有 length 方法的可迭代对象，但 SVD 没有 length 方法。)
但是，我建议不要在这种情况下预先分配一个向量。首先我们看一下 F的类型:

julia> typeof(Vector{SVD}(undef, 10))
Array{SVD,1}

这个向量的问题在于它是由抽象类型参数化的。有一个 section在手册的性能提示一章中建议不要这样做。 SVD是一个抽象类型，因为它的参数类型没有被指定。为了使其具体化，您需要指定参数的类型，如下所示:

julia> SVD{Float64,Float64,Array{Float64,2}}
SVD{Float64,Float64,Array{Float64,2}}

julia> Vector{SVD{Float64,Float64,Array{Float64,2}}}(undef, 2)
2-element Array{SVD{Float64,Float64,Array{Float64,2}},1}:
 #undef
 #undef

如您所见，当您处理复杂类型(如 SVD)时，很难正确指定具体类型。 .此外，如果您这样做，您的代码将不会像应有的那样通用。
解决此类问题的更好方法是使用映射、广播或列表推导。然后会自动推断出正确的输出类型。这里有些例子:
列表理解

julia> [svd(rand(3, 3)) for _ in 1:2]
2-element Array{SVD{Float64,Float64,Array{Float64,2}},1}:
 SVD{Float64,Float64,Array{Float64,2}}([-0.6357040496635746 -0.2941425771794837 -0.7136949667270628; -0.45459999623274916 -0.6045700314848496 0.654090147040599; -0.6238743500629883 0.7402534845042064 0.2506104028424691], [1.4535849689665463, 0.7212190827260345, 0.05010669163393896], [-0.5975505057447164 -0.588792736048385 -0.5442945039782142; 0.7619724725128861 -0.6283345569895092 -0.15682358121595258; -0.2496624605679292 -0.5084474392397449 0.8241054891903787])
 SVD{Float64,Float64,Array{Float64,2}}([-0.5593632049776268 0.654338345992878 -0.5088753618327984; -0.6687620652652163 -0.7189576326033171 -0.18936003428293915; -0.4897653570633183 0.23439550227070827 0.8397551092645418], [1.8461274187259178, 0.21226179692488983, 0.14194607536315287], [-0.29089551972856004 -0.7086270946133293 -0.6428276887173754; -0.9203610429640889 0.023709029028269546 0.390350397126212; 0.2613720474647311 -0.7051847436823973 0.6590896221923739])

map

julia> map(_ -> svd(rand(3, 3)), 1:2)
2-element Array{SVD{Float64,Float64,Array{Float64,2}},1}:
 SVD{Float64,Float64,Array{Float64,2}}([-0.5807809149601634 0.5635242755434755 0.5874809951745127; -0.6884131975465821 0.0451903888051729 -0.7239095925620322; -0.43448912329507794 -0.8248625459025509 0.3616918330643316], [1.488618654040125, 0.4122166626927311, 0.004235624485479941], [-0.6721098925787947 -0.2684664121709399 -0.6900681689759235; -0.7384292974335966 0.31185073633575333 0.5978890289498324; -0.05468514413847799 -0.9114136842196914 0.4078414290231468])
 SVD{Float64,Float64,Array{Float64,2}}([-0.3677873424759118 0.8090638526628051 -0.4584191892023337; -0.43071684640222546 -0.5851169278783189 -0.6871107472129654; -0.8241452960126802 -0.055261768200600137 0.5636760310989947], [1.6862363968739773, 0.5899255050748418, 0.24246688716190598], [-0.3751742784957875 -0.7172409091515735 -0.5872050229643736; 0.8600668700980193 -0.505618838823938 0.06807766730822862; -0.3457300098559026 -0.4794945964927631 0.8065703268899])

广播

julia> g = (rand(3, 3) for _ in 1:2)
Base.Generator{UnitRange{Int64},var"#17#18"}(var"#17#18"(), 1:2)

julia> svd.(g)
2-element Array{SVD{Float64,Float64,Array{Float64,2}},1}:
 SVD{Float64,Float64,Array{Float64,2}}([-0.7988295268840152 0.5443221484534134 -0.256095266807727; -0.5436890668169485 -0.8354777569473182 -0.0798693700362902; -0.257436566171119 0.07543418554831638 0.963346302244777], [1.8188722412547844, 0.3934389096422389, 0.2020398396772306], [-0.7147404794808727 -0.37763644211761316 -0.5886737335538281; -0.6944558966482991 0.4830041206449164 0.5333273169925189; -0.08292800854873916 -0.7899985677359054 0.607474450798845])
 SVD{Float64,Float64,Array{Float64,2}}([-0.5910620103531503 0.3599866268397522 0.7218416228050514; -0.7367495542691711 0.12340124384185132 -0.664809918173956; -0.3283988340440176 -0.9247603805931685 0.1922821996018057], [1.826019614357666, 0.5333148215847028, 0.11639139812894106], [-0.6415954756495915 -0.6888196183142843 -0.33746522643279503; -0.5845558664639438 0.7239484700883465 -0.3663236978948133; -0.4966383841474222 0.037764349353666515 0.8671356118331964])

此外，映射、广播和列表推导应该与预分配向量一样有效。如果您正在做一个简单的映射，那么使用映射、广播或列表推导通常更容易且更具可读性。预分配向量是我为从头开始编写自定义算法而保留的工具。
最后一点。在大多数情况下，类型参数被视为实现细节，而不是类型公共(public) API 的一部分。因此，最好使用不依赖于固定类型参数类型的泛型编程方法。当然，这条经验法则也有一些异常(exception)，例如 Array{T,N}和 Dict{K,V} .