Iterative Solution of Sparse Linear Least Squares using LU Factorization

In this paper, we are interested in computing the solution of an overdetermined sparse linear least squares problem Ax=b via the normal equations method. Transforming the normal equations using the L factor from a rectangular LU decomposition of A usually leads to a better conditioned problem. Here we explore a further preconditioning by inv(L1) where L1 is the n × n upper part of the lower trapezoidal m × n factor L. Since the condition number of the iteration matrix can be easily bounded, we can determine whether the iteration will be effective, and whether further pre-conditioning is required. Numerical experiments are performed with the Julia programming language. When the upper triangular matrix U has no near zero diagonal elements, the algorithm is observed to be reliable. When A has only a few more rows than columns, convergence requires relatively few iterations and the algorithm usually requires less storage than the Cholesky factor of AtA or the R factor of the QR factorization of A.