MA75: Sparse over-determined system: weighted least squares

This subroutine solves weighted sparse least-squares problems. Given an \(m \times n\) (\(m \ge n\)) sparse matrix \(\mathbf{A } = \{ a _{ij} \}\) of rank \(n\), an \(m \times m\) diagonal matrix \(\mathbf{W}\) of weights, and an \(m\)-vector \(\mathbf{ b}\), the routine calculates the solution vector \(\mathbf{ x}\) that minimizes the Euclidean norm of the weighted residual vector \(\mathbf{ r } = \mathbf{W} (\mathbf{A} \mathbf{ x } - \mathbf{ b})\) by solving the normal equations \(\mathbf{A} ^ T \mathbf{W} ^2 \mathbf{A}\mathbf{x} = \mathbf{A} ^ T \mathbf{W} ^2 \mathbf{b}\).

Three forms of data storage are permitted for the input matrix: storage by columns, where row indices and column pointers describe the matrix; storage by rows, where column indices and row pointers describe the matrix; and the coordinate scheme, where both row and column indices describe the position of entries in the matrix. For the statistical analysis of the weighted least-squares problem, there are two entries: one to obtain a column and one to obtain the diagonal of the covariance matrix \(( \mathbf{A} ^ T \mathbf{W} ^2 \mathbf{A} ) ^{-1}\).