When doing linear algebra procedures such as
A\b
or pinv(A)
,
where A
is a square matrix, is there a way to check for singularity of A
? Otherwise I keep getting errors such as SingularException(3)
or LAPACKException(1)
.
pretty new to julia, but I think this would work:
is_singular = A → rank(A) == max(size(A)…)
#only true for matrix A when square and no dependent columns (full rank) = singular
Numerically, this question is a bit ill-defined, because roundoff errors make it hard to distinguish between a matrix that is exactly singular and a matrix that is nearly singular. Moreover, as a practical matter it is rarely a good idea to try to make this distinction — matrices that are nearly singular are in many ways just as “bad” as matrices that are singular.
If you find yourself solving a lot of systems that are nearly singular, that it is a good sign that you need to re-think what you are doing. (e.g. perhaps you need some regularization in your equations)
Where are these matrices coming from? Why are they sometimes singular or nearly so?
Just adding that checking the condition number of a matrix can be done with:
julia> using LinearAlgebra
julia> M = rand(4,4)
4×4 Matrix{Float64}:
0.596334 0.699729 0.632193 0.931758
0.910922 0.150834 0.605005 0.90234
0.640964 0.0686345 0.0811358 0.372114
0.67283 0.458002 0.346649 0.485548
julia> cond(M)
23.45879906432043
julia> cond(zeros(4,4)) # zero matrix is singular
Inf
and for any singular matrix, the condition number would be Inf
. A smaller condition number is better.
Very nice point. I shall follow it.
Thank you, that is quite handy