JULIA 中没有旋转的 LU 分解

Question

尝试将此答案中的 lu_nopivot 重写为 JULIA 并仅使用一个循环。我写的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) 在 L[k+1:n,k]*U[k,:] 上给出错误 MethodError: no method matching *(::Vector{Float64}, ::Vector{Float64})。这是什么原因？

Answer 1

问题是，当您执行 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 选项：

Pivoting can be turned off by passing 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).

JULIA 中没有旋转的 LU 分解

LU decomposition without pivoting in JULIA

matlab

linear-algebra

numerical-methods

julia