No. There is absolutely no presumption that lsqr (or ANY such tool) will give an exact solution in terms of essentially a zero residual. Consider this example:
A = rand(3,2)
A =
0.81472 0.91338
0.90579 0.63236
0.12699 0.09754
b = rand(3,1)
b =
0.2785
0.54688
0.95751
x = lsqr(A,b)
lsqr converged at iteration 2 to a solution with relative residual 0.77.
x =
1.2896
-0.82073
As you can see, it agrees with that which backslash returns.
A\b
ans =
1.2896
-0.82073
Is it exact? OF COURSE NOT!
A*x-b
ans =
0.022536
0.10223
-0.8738
No exact solution exists. This is what you have, a problem with no exact solution.
Merely setting a small tolerance cannot enforce a solution to exist where no solution exists. The tolerance only tells lsqr when to stop iterating. So lets see what you THINK you read in the documentation.
[X,FLAG,RELRES] = lsqr(A,B,...) also returns estimates of the relative
residual NORM(B-A*X)/NORM(B). If RELRES <= TOL, then X is a
consistent solution to A*X=B. If FLAG is 0 but RELRES > TOL, then X is
the least squares solution which minimizes norm(B-A*X).
LSQR returns a relative residual, but as you can read above, there is no presumption that the solution is exact.
