为什么 Julia 无法系统地求解线性系统?
Why does Julia fails to solve linear system systematically?
方块 A 的 Ax=b 问题由 \ 函数解决。考虑到这一点,我尝试执行以下操作:
A = rand(1:4,3,3)
x = fill(1.0, 3)
b = A * x
A\b
出于某种原因,代码有时似乎可以正常工作。但有时 returns 我会出现以下错误:
LinearAlgebra.SingularException(3)
Stacktrace:
[1] checknonsingular
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/factorization.jl:19 [inlined]
[2] checknonsingular
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/factorization.jl:21 [inlined]
[3] #lu!#136
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:85 [inlined]
[4] #lu#140
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:273 [inlined]
[5] lu (repeats 2 times)
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:272 [inlined]
[6] \(A::Matrix{Int64}, B::Vector{Float64})
@ LinearAlgebra /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/generic.jl:1136
[7] top-level scope
@ In[208]:4
[8] eval
@ ./boot.jl:360 [inlined]
[9] include_string(mapexpr::typeof(REPL.softscope), mod::Module, code::String, filename::String)
@ Base ./loading.jl:1116
所以,我试图了解发生了什么,并执行了代码 10000000 次,发现它有 10% 的执行失败。
using Printf
i = 0
test = 10000000
for x in 1:test
try
A = rand(1:4,3,3)
x = fill(1.0, 3)
b = A * x
A\b
catch
i = i+1
end
end
fail_percentage = (i/test)*100
@printf "this code has failed in %.2f%%" fail_percentage
谁能给我解释一下这是怎么回事?
错误很明显:LinearAlgebra.SingularException
。这不是 Julia 的失败,而是方程组的属性。
如果矩阵 A
是奇异的,则没有单一的解决方案 - 如果系统是齐次的,则有无限数量的解决方案,或者在一般情况下 none。似乎您已经根据经验计算了使用您测试的属性生成奇异系统的概率(A
和 x
的维度,x
填充 1,A
填充 1 和4).
如果像 OP 一样,您希望跳过奇异矩阵,则需要确保 A
的行列式不为 0。您可以使用内置函数检查并跳过此类矩阵,或者通过注意行列式本身是一个方程来生成它们,因此,对于没有 0 个条目的 3x3 示例,您选择 8 个数字,您可以计算第 9 个不能是什么以确保行列式不为零。如果允许 0,则需要检查所有可能性。
方块 A 的 Ax=b 问题由 \ 函数解决。考虑到这一点,我尝试执行以下操作:
A = rand(1:4,3,3)
x = fill(1.0, 3)
b = A * x
A\b
出于某种原因,代码有时似乎可以正常工作。但有时 returns 我会出现以下错误:
LinearAlgebra.SingularException(3)
Stacktrace:
[1] checknonsingular
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/factorization.jl:19 [inlined]
[2] checknonsingular
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/factorization.jl:21 [inlined]
[3] #lu!#136
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:85 [inlined]
[4] #lu#140
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:273 [inlined]
[5] lu (repeats 2 times)
@ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:272 [inlined]
[6] \(A::Matrix{Int64}, B::Vector{Float64})
@ LinearAlgebra /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/generic.jl:1136
[7] top-level scope
@ In[208]:4
[8] eval
@ ./boot.jl:360 [inlined]
[9] include_string(mapexpr::typeof(REPL.softscope), mod::Module, code::String, filename::String)
@ Base ./loading.jl:1116
所以,我试图了解发生了什么,并执行了代码 10000000 次,发现它有 10% 的执行失败。
using Printf
i = 0
test = 10000000
for x in 1:test
try
A = rand(1:4,3,3)
x = fill(1.0, 3)
b = A * x
A\b
catch
i = i+1
end
end
fail_percentage = (i/test)*100
@printf "this code has failed in %.2f%%" fail_percentage
谁能给我解释一下这是怎么回事?
错误很明显:LinearAlgebra.SingularException
。这不是 Julia 的失败,而是方程组的属性。
如果矩阵 A
是奇异的,则没有单一的解决方案 - 如果系统是齐次的,则有无限数量的解决方案,或者在一般情况下 none。似乎您已经根据经验计算了使用您测试的属性生成奇异系统的概率(A
和 x
的维度,x
填充 1,A
填充 1 和4).
如果像 OP 一样,您希望跳过奇异矩阵,则需要确保 A
的行列式不为 0。您可以使用内置函数检查并跳过此类矩阵,或者通过注意行列式本身是一个方程来生成它们,因此,对于没有 0 个条目的 3x3 示例,您选择 8 个数字,您可以计算第 9 个不能是什么以确保行列式不为零。如果允许 0,则需要检查所有可能性。