Presently, nullspace is essentially a simple wrapper over svd, picking out the relevant rows of Vt:
|
function nullspace(A::AbstractVecOrMat; atol::Real=0, rtol::Real = (min(size(A, 1), size(A, 2))*eps(real(float(oneunit(eltype(A))))))*iszero(atol)) |
|
m, n = size(A, 1), size(A, 2) |
|
(m == 0 || n == 0) && return Matrix{eigtype(eltype(A))}(I, n, n) |
|
SVD = svd(A; full=true) |
|
tol = max(atol, SVD.S[1]*rtol) |
|
indstart = sum(s -> s .> tol, SVD.S) + 1 |
|
return copy((@view SVD.Vt[indstart:end,:])') |
|
end |
By default, svd uses a divide-and-conquer algorithm (LinearAlgebra.DivideAndConquer() -> gesdd), but it would be nice to be able to choose other algorithms - namely the QR iteration algorithm (LinearAlgebra.QRIteration() -> gesvd).
The purpose of this issue is therefore to suggest that an alg or svd_alg argument is added to nullspace, which is simply forwarded to svd.
The motivation is that the divide-and-conquer approach sometimes fails to converge the SVD, especially for very large and/or nearly rank-deficient problems.
Presently,
nullspaceis essentially a simple wrapper oversvd, picking out the relevant rows ofVt:LinearAlgebra.jl/src/dense.jl
Lines 1881 to 1888 in 17325b7
By default,
svduses a divide-and-conquer algorithm (LinearAlgebra.DivideAndConquer()->gesdd), but it would be nice to be able to choose other algorithms - namely the QR iteration algorithm (LinearAlgebra.QRIteration()->gesvd).The purpose of this issue is therefore to suggest that an
algorsvd_algargument is added tonullspace, which is simply forwarded tosvd.The motivation is that the divide-and-conquer approach sometimes fails to converge the SVD, especially for very large and/or nearly rank-deficient problems.