試圖lu_nopivot從這個答案https://stackoverflow.com/a/41151228重寫為 JULIA 并僅使用一個回圈。我寫的朱莉婭代碼
using LinearAlgebra
function lu_nopivot(A)
n = size(A, 1)
L = Matrix{eltype(A)}(I, n, n)
U = copy(A)
for k = 1:n
L[k 1:n,k] = U[k 1:n,k] / U[k,k]
U[k 1:n,:] = U[k 1:n,:] - L[k 1:n,k]*U[k,:]
end
return L, U
end
但是呼叫該函式
L, U = lu_nopivot(A)會MethodError: no method matching *(::Vector{Float64}, ::Vector{Float64})在L[k 1:n,k]*U[k,:]. 這是什么原因？
uj5u.com熱心網友回復：
問題是當你這樣做時U[k, :]，即使你提取了一行，你得到的是一個列向量。因此L[k 1:n,k]*U[k,:]嘗試將列向量乘以列向量。
在 Julia中獲取行向量又名1 x N矩陣的一種方法（盡管我不知道這是否是慣用方式）是U[k:k, :]：
U[k 1:n,:] = U[k 1:n,:] - L[k 1:n,k] * U[k:k,:]
但是請注意，Julia 的實作lu已經有一個 no-pivot 選項：
可以通過傳遞 pivot = NoPivot() 來關閉旋轉。
julia> M = rand(1.0:100.0, 3, 3)
3×3 Matrix{Float64}:
71.0 50.0 23.0
82.0 63.0 37.0
96.0 28.0 5.0
julia> L, U = lu_nopivot(M); # after the above k:k change has been made
julia> L
3×3 Matrix{Float64}:
1.0 0.0 0.0
1.15493 1.0 0.0
1.35211 -7.53887 1.0
julia> U
3×3 Matrix{Float64}:
71.0 50.0 23.0
0.0 5.25352 10.4366
0.0 0.0 52.5818
julia> lu(M, NoPivot())
LU{Float64, Matrix{Float64}}
L factor:
3×3 Matrix{Float64}:
1.0 0.0 0.0
1.15493 1.0 0.0
1.35211 -7.53887 1.0
U factor:
3×3 Matrix{Float64}:
71.0 50.0 23.0
0.0 5.25352 10.4366
0.0 0.0 52.5818
使用它可能更高效，也更健壯（例如，您當前的實作無法處理 eltype 的矩陣Int，因為 L 和 U 的型別與輸入相同，但需要包含Floats）。
轉載請註明出處，本文鏈接：https://www.uj5u.com/yidong/440851.html