MA49: Sparse over-determined system: least squares by QR

To compute an orthogonal factorization of a sparse overdetermined matrix \(\mathbf{A}\) and optionally to solve the least squares problem \(\mathrm{ min} ||\mathbf{b} -\mathbf{Ax}||_2\). Given a sparse matrix \(\mathbf{A}\) of order \(m \times n\), \(m \ge n\), of full column rank, this subroutine computes the \(\mathbf{QR}\) factorization \(\mathbf{A} = \mathbf{Q} \left [ \begin{array}{c} \mathbf{R} \\ \mathbf{0} \end{array} \right ]\) where \(\mathbf{Q}\) is an \(m \times m\) orthogonal matrix and \(\mathbf{R}\) is an \(n \times n\) upper triangular matrix.

Given an \(n\)-vector \(\mathbf{b}\), this subroutine may also compute the minimum 2-norm solution of the linear system \(\mathbf{A} ^T \mathbf{x} = \mathbf{b}\), by solving \(\left [ \mathbf{R} ^T \mathbf{0} \right ] \mathbf{z} = \mathbf{b}\) and performing the multiplication \(\mathbf{x} = \mathbf{Q} \mathbf{z}\), or, if the \(\mathbf{Q}\) factor is not stored, by solving \(\mathbf{R} ^T \mathbf{R} \mathbf{z} = \mathbf{b}\) and performing the multiplication \(\mathbf{x} = \mathbf{A} \mathbf{z}\).

The subroutine can also solve systems with the coefficient matrix \(\left [ \begin{array}{c} \mathbf{R} \\ \mathbf{0} \end{array} \right ]\) or \([ {\mathbf{R} ^T \mathbf{0}} ]\), or will compute the product of \(\mathbf{Q}\) or \(\mathbf{Q} ^T\) with a vector.

The method used is based on the multifrontal approach and makes use of Householder transformations. Because an ordering for the columns is chosen using the pattern of the matrix \(\mathbf{A} ^T\)\(\mathbf{A}\), this code is not designed for matrices with full rows.