|
|
Ublas : |
From: Paul C. Leopardi (paul.leopardi_at_[hidden])
Date: 2007-02-13 17:05:11
Hi all,
Perhaps cs_qr.c from
http://www.cise.ufl.edu/research/sparse/CXSparse/CXSparse/Source/cs_qr.c
could be adapted? It is LGPL:
http://www.cise.ufl.edu/research/sparse/CXSparse/CXSparse/License.txt
See "Direct Methods for Sparse Linear Systems: the CXSparse package"
Tim Davis, http://www.cise.ufl.edu/research/sparse/CXSparse/
and
Direct Methods for Sparse Linear Systems, T. A. Davis, SIAM, Philadelphia,
Sept. 2006.
I have not tried cs_qr.c myself.
See also
MR1438098 (98b:65049) Matstoms, Pontus
"Sparse linear least squares problems in optimization."
Computational issues in high performance software for nonlinear optimization
(Capri, 1995). Comput. Optim. Appl. 7 (1997), no. 1, 89--110.
doi:10.1023/A:1008680131271
http://www.springerlink.com/content/n6j160823847l397/
and
Pontus Matstoms,
"Sparse QR factorization in MATLAB", 1994,
http://portal.acm.org/citation.cfm?id=174603.174408
Best, Paul
On Wed, 14 Feb 2007, Karl Meerbergen wrote:
> Gunter Winkler wrote:
> >Am Dienstag, 13. Februar 2007 00:28 schrieb
> >
> >Richard_Harding_at_[hidden]:
> >>Hello,
> >>
> >>I have a very sparse matrix (<= 6 nnz entries per row with
> >>potentially thousands of columns) which I need to pre-multiply by its
> >>transpose in order to solve the normal equations associated with a
> >>least-squares computation: (A'A)x = A'b. I may have missed them, but
> >
> >Don't do this - trust me ;-)
> >
> >in order to solve the least squares problem you do a QR-decomposition of
> >A and get the cholesky decomposition of A'A for free:
> >
> >A = QR -> A'A = (QR)'QR = R'Q'QR= R'R
>
> True. But you need a sparse QR factorization. It exists, but I do not
> recall precisely where you can find code. A sparse QR can be done with a
> sparse Gram-Schmidt routine. Needless to say that fill-in is going to
> grow. Developing an efficient code using BLAS 3 etc is quite a job.
>
> Karl
>
> Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm