JuliaLinearAlgebra / IterativeSolvers.jl

Iterative algorithms for solving linear systems, eigensystems, and singular value problems
MIT License
394 stars 106 forks source link

lsqr! not minimizing norm #327

Open FHell opened 1 year ago

FHell commented 1 year ago

According to the documentation lsqr! should return the minimum norm solution if the matrix is ambiguous. This doesn't seem to work:

using IterativeSolvers

#-

A = zeros(5,5)
b = ones(5)
x = ones(5)

#-

lsqr!(x, A, b; damp = 1.)

#-

isapprox(x, ones(5)) # true

lsqr is working as expected:

y = lsqr(A, b; damp = 1.)

isapprox(y, zeros(5)) # true

Julia 1.8, IterativeSolvers v0.9.2

FHell commented 1 year ago

Incidentally: lsmr! and lsmr show the same behavior.

Edit:

And a nontrivial example:

x = ones(4)
b = ones(4)

v1 = rand(4)
v2 = rand(4)
v3 = rand(4)

w1 = rand(4)
w2 = rand(4)
w3 = rand(4)

A = v1 * transpose(w1) + v2 * transpose(w2) + v3 * transpose(w3)

#-

lsqr!(x, A, b)
y = lsqr(A, b)

#-

isapprox(x, y) # false