3
<title>QR Decomposition with Column Pivoting - GNU Scientific Library -- Reference Manual</title>
4
<meta http-equiv="Content-Type" content="text/html">
5
<meta name="description" content="GNU Scientific Library -- Reference Manual">
6
<meta name="generator" content="makeinfo 4.8">
7
<link title="Top" rel="start" href="index.html#Top">
8
<link rel="up" href="Linear-Algebra.html#Linear-Algebra" title="Linear Algebra">
9
<link rel="prev" href="QR-Decomposition.html#QR-Decomposition" title="QR Decomposition">
10
<link rel="next" href="Singular-Value-Decomposition.html#Singular-Value-Decomposition" title="Singular Value Decomposition">
11
<link href="http://www.gnu.org/software/texinfo/" rel="generator-home" title="Texinfo Homepage">
13
Copyright (C) 1996, 1997, 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 The GSL Team.
15
Permission is granted to copy, distribute and/or modify this document
16
under the terms of the GNU Free Documentation License, Version 1.2 or
17
any later version published by the Free Software Foundation; with the
18
Invariant Sections being ``GNU General Public License'' and ``Free Software
19
Needs Free Documentation'', the Front-Cover text being ``A GNU Manual'',
20
and with the Back-Cover Text being (a) (see below). A copy of the
21
license is included in the section entitled ``GNU Free Documentation
24
(a) The Back-Cover Text is: ``You have freedom to copy and modify this
25
GNU Manual, like GNU software.''-->
26
<meta http-equiv="Content-Style-Type" content="text/css">
27
<style type="text/css"><!--
28
pre.display { font-family:inherit }
29
pre.format { font-family:inherit }
30
pre.smalldisplay { font-family:inherit; font-size:smaller }
31
pre.smallformat { font-family:inherit; font-size:smaller }
32
pre.smallexample { font-size:smaller }
33
pre.smalllisp { font-size:smaller }
34
span.sc { font-variant:small-caps }
35
span.roman { font-family:serif; font-weight:normal; }
36
span.sansserif { font-family:sans-serif; font-weight:normal; }
42
<a name="QR-Decomposition-with-Column-Pivoting"></a>
43
Next: <a rel="next" accesskey="n" href="Singular-Value-Decomposition.html#Singular-Value-Decomposition">Singular Value Decomposition</a>,
44
Previous: <a rel="previous" accesskey="p" href="QR-Decomposition.html#QR-Decomposition">QR Decomposition</a>,
45
Up: <a rel="up" accesskey="u" href="Linear-Algebra.html#Linear-Algebra">Linear Algebra</a>
49
<h3 class="section">13.3 QR Decomposition with Column Pivoting</h3>
51
<p><a name="index-QR-decomposition-with-column-pivoting-1223"></a>
52
The QR decomposition can be extended to the rank deficient case
53
by introducing a column permutation P,
54
The first r columns of Q form an orthonormal basis
55
for the range of A for a matrix with column rank r. This
56
decomposition can also be used to convert the linear system A x =
57
b into the triangular system R y = Q^T b, x = P y, which can be
58
solved by back-substitution and permutation. We denote the QR
59
decomposition with column pivoting by QRP^T since A = Q R
63
— Function: int <b>gsl_linalg_QRPT_decomp</b> (<var>gsl_matrix * A, gsl_vector * tau, gsl_permutation * p, int * signum, gsl_vector * norm</var>)<var><a name="index-gsl_005flinalg_005fQRPT_005fdecomp-1224"></a></var><br>
64
<blockquote><p>This function factorizes the M-by-N matrix <var>A</var> into
65
the QRP^T decomposition A = Q R P^T. On output the
66
diagonal and upper triangular part of the input matrix contain the
67
matrix R. The permutation matrix P is stored in the
68
permutation <var>p</var>. The sign of the permutation is given by
69
<var>signum</var>. It has the value (-1)^n, where n is the
70
number of interchanges in the permutation. The vector <var>tau</var> and the
71
columns of the lower triangular part of the matrix <var>A</var> contain the
72
Householder coefficients and vectors which encode the orthogonal matrix
73
<var>Q</var>. The vector <var>tau</var> must be of length k=\min(M,N). The
74
matrix Q is related to these components by, Q = Q_k ... Q_2
75
Q_1 where Q_i = I - \tau_i v_i v_i^T and v_i is the
76
Householder vector v_i =
77
(0,...,1,A(i+1,i),A(i+2,i),...,A(m,i)). This is the same storage scheme
78
as used by <span class="sc">lapack</span>. The vector <var>norm</var> is a workspace of length
79
<var>N</var> used for column pivoting.
81
<p>The algorithm used to perform the decomposition is Householder QR with
82
column pivoting (Golub & Van Loan, <cite>Matrix Computations</cite>, Algorithm
84
</p></blockquote></div>
87
— Function: int <b>gsl_linalg_QRPT_decomp2</b> (<var>const gsl_matrix * A, gsl_matrix * q, gsl_matrix * r, gsl_vector * tau, gsl_permutation * p, int * signum, gsl_vector * norm</var>)<var><a name="index-gsl_005flinalg_005fQRPT_005fdecomp2-1225"></a></var><br>
88
<blockquote><p>This function factorizes the matrix <var>A</var> into the decomposition
89
A = Q R P^T without modifying <var>A</var> itself and storing the
90
output in the separate matrices <var>q</var> and <var>r</var>.
91
</p></blockquote></div>
94
— Function: int <b>gsl_linalg_QRPT_solve</b> (<var>const gsl_matrix * QR, const gsl_vector * tau, const gsl_permutation * p, const gsl_vector * b, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQRPT_005fsolve-1226"></a></var><br>
95
<blockquote><p>This function solves the square system A x = b using the QRP^T
96
decomposition of A into (<var>QR</var>, <var>tau</var>, <var>p</var>) given by
97
<code>gsl_linalg_QRPT_decomp</code>.
98
</p></blockquote></div>
101
— Function: int <b>gsl_linalg_QRPT_svx</b> (<var>const gsl_matrix * QR, const gsl_vector * tau, const gsl_permutation * p, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQRPT_005fsvx-1227"></a></var><br>
102
<blockquote><p>This function solves the square system A x = b in-place using the
103
QRP^T decomposition of A into
104
(<var>QR</var>,<var>tau</var>,<var>p</var>). On input <var>x</var> should contain the
105
right-hand side b, which is replaced by the solution on output.
106
</p></blockquote></div>
109
— Function: int <b>gsl_linalg_QRPT_QRsolve</b> (<var>const gsl_matrix * Q, const gsl_matrix * R, const gsl_permutation * p, const gsl_vector * b, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQRPT_005fQRsolve-1228"></a></var><br>
110
<blockquote><p>This function solves the square system R P^T x = Q^T b for
111
<var>x</var>. It can be used when the QR decomposition of a matrix is
112
available in unpacked form as (<var>Q</var>, <var>R</var>).
113
</p></blockquote></div>
116
— Function: int <b>gsl_linalg_QRPT_update</b> (<var>gsl_matrix * Q, gsl_matrix * R, const gsl_permutation * p, gsl_vector * u, const gsl_vector * v</var>)<var><a name="index-gsl_005flinalg_005fQRPT_005fupdate-1229"></a></var><br>
117
<blockquote><p>This function performs a rank-1 update w v^T of the QRP^T
118
decomposition (<var>Q</var>, <var>R</var>, <var>p</var>). The update is given by
119
Q'R' = Q R + w v^T where the output matrices Q' and
120
R' are also orthogonal and right triangular. Note that <var>w</var> is
121
destroyed by the update. The permutation <var>p</var> is not changed.
122
</p></blockquote></div>
125
— Function: int <b>gsl_linalg_QRPT_Rsolve</b> (<var>const gsl_matrix * QR, const gsl_permutation * p, const gsl_vector * b, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQRPT_005fRsolve-1230"></a></var><br>
126
<blockquote><p>This function solves the triangular system R P^T x = b for the
127
N-by-N matrix R contained in <var>QR</var>.
128
</p></blockquote></div>
131
— Function: int <b>gsl_linalg_QRPT_Rsvx</b> (<var>const gsl_matrix * QR, const gsl_permutation * p, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQRPT_005fRsvx-1231"></a></var><br>
132
<blockquote><p>This function solves the triangular system R P^T x = b in-place
133
for the N-by-N matrix R contained in <var>QR</var>. On
134
input <var>x</var> should contain the right-hand side b, which is
135
replaced by the solution on output.
136
</p></blockquote></div>