QR Decomposition with Column Pivoting
The QR decomposition can be extended to the rank deficient case
by introducing a column permutation P,
A P = Q R
The first r columns of this Q form an orthonormal basis
for the range of A for a matrix with column rank r. This
decomposition can also be used to convert the linear system A x =
b into the triangular system R y = Q^T b, x = P y, which can be
solved by back-substitution and permutation. We denote the QR
decomposition with column pivoting by QRP^T since A = Q R
P^T.
- Function: int gsl_linalg_QRPT_decomp (gsl_matrix * A, gsl_vector * tau, gsl_permutation * p, int *signum, gsl_vector * norm)
-
This function factorizes the M-by-N matrix A into
the QRP^T decomposition A = Q R P^T. On output the
diagonal and upper triangular part of the input matrix contain the
matrix R. The permutation matrix P is stored in the
permutation p. The sign of the permutation is given by
signum. It has the value (-1)^n, where n is the
number of interchanges in the permutation. The vector tau and the
columns of the lower triangular part of the matrix A contain the
Householder coefficients and vectors which encode the orthogonal matrix
Q. The vector tau must be of length k=\min(M,N). The
matrix Q is related to these components by, Q = Q_k ... Q_2
Q_1 where Q_i = I - \tau_i v_i v_i^T and v_i is the
Householder vector v_i =
(0,...,1,A(i+1,i),A(i+2,i),...,A(m,i)). This is the same storage scheme
as used by LAPACK. On output the norms of each column of R
are stored in the vector norm.
The algorithm used to perform the decomposition is Householder QR with
column pivoting (Golub & Van Loan, Matrix Computations, Algorithm
5.4.1).
- Function: int gsl_linalg_QRPT_decomp2 (const gsl_matrix * A, gsl_matrix * q, gsl_matrix * r, gsl_vector * tau, gsl_permutation * p, int *signum, gsl_vector * norm)
-
This function factorizes the matrix A into the decomposition
A = Q R P^T without modifying A itself and storing the
output in the separate matrices q and r.
- Function: int gsl_linalg_QRPT_solve (const gsl_matrix * QR, const gsl_vector * tau, const gsl_permutation * p, const gsl_vector * b, gsl_vector * x)
-
This function solves the system A x = b using the QRP^T
decomposition of A into (QR, tau, p) given by
gsl_linalg_QRPT_decomp
.
- Function: int gsl_linalg_QRPT_svx (const gsl_matrix * QR, const gsl_vector * tau, const gsl_permutation * p, gsl_vector * x)
-
This function solves the system A x = b in-place using the
QRP^T decomposition of A into
(QR,tau,p). On input x should contain the
right-hand side b, which is replaced by the solution on output.
- Function: int gsl_linalg_QRPT_QRsolve (const gsl_matrix * Q, const gsl_matrix * R, const gsl_permutation * p, const gsl_vector * b, gsl_vector * x)
-
This function solves the system R P^T x = Q^T b for x. It can
be used when the QR decomposition of a matrix is available in
unpacked form as (Q,R).
- Function: int gsl_linalg_QRPT_update (gsl_matrix * Q, gsl_matrix * R, const gsl_permutation * p, gsl_vector * u, const gsl_vector * v)
-
This function performs a rank-1 update w v^T of the QRP^T
decomposition (Q, R,p). The update is given by
Q'R' = Q R + w v^T where the output matrices Q' and
R' are also orthogonal and right triangular. Note that w is
destroyed by the update. The permutation p is not changed.
- Function: int gsl_linalg_QRPT_Rsolve (const gsl_matrix * QR, const gsl_permutation * p, const gsl_vector * b, gsl_vector * x)
-
This function solves the triangular system R P^T x = b for the
N-by-N matrix R contained in QR.
- Function: int gsl_linalg_QRPT_Rsvx (const gsl_matrix * QR, const gsl_permutation * p, gsl_vector * x)
-
This function solves the triangular system R P^T x = b in-place
for the N-by-N matrix R contained in QR. On
input x should contain the right-hand side b, which is
replaced by the solution on output.