1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
|
<!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/html4/loose.dtd">
<html>
<!-- Copyright (C) 1996, 1997, 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012, 2013, 2014, 2015, 2016 The GSL Team.
Permission is granted to copy, distribute and/or modify this document
under the terms of the GNU Free Documentation License, Version 1.3 or
any later version published by the Free Software Foundation; with the
Invariant Sections being "GNU General Public License" and "Free Software
Needs Free Documentation", the Front-Cover text being "A GNU Manual",
and with the Back-Cover Text being (a) (see below). A copy of the
license is included in the section entitled "GNU Free Documentation
License".
(a) The Back-Cover Text is: "You have the freedom to copy and modify this
GNU Manual." -->
<!-- Created by GNU Texinfo 5.1, http://www.gnu.org/software/texinfo/ -->
<head>
<title>GNU Scientific Library – Reference Manual: QR Decomposition with Column Pivoting</title>
<meta name="description" content="GNU Scientific Library – Reference Manual: QR Decomposition with Column Pivoting">
<meta name="keywords" content="GNU Scientific Library – Reference Manual: QR Decomposition with Column Pivoting">
<meta name="resource-type" content="document">
<meta name="distribution" content="global">
<meta name="Generator" content="makeinfo">
<meta http-equiv="Content-Type" content="text/html; charset=utf-8">
<link href="index.html#Top" rel="start" title="Top">
<link href="Function-Index.html#Function-Index" rel="index" title="Function Index">
<link href="Linear-Algebra.html#Linear-Algebra" rel="up" title="Linear Algebra">
<link href="Complete-Orthogonal-Decomposition.html#Complete-Orthogonal-Decomposition" rel="next" title="Complete Orthogonal Decomposition">
<link href="QR-Decomposition.html#QR-Decomposition" rel="previous" title="QR Decomposition">
<style type="text/css">
<!--
a.summary-letter {text-decoration: none}
blockquote.smallquotation {font-size: smaller}
div.display {margin-left: 3.2em}
div.example {margin-left: 3.2em}
div.indentedblock {margin-left: 3.2em}
div.lisp {margin-left: 3.2em}
div.smalldisplay {margin-left: 3.2em}
div.smallexample {margin-left: 3.2em}
div.smallindentedblock {margin-left: 3.2em; font-size: smaller}
div.smalllisp {margin-left: 3.2em}
kbd {font-style:oblique}
pre.display {font-family: inherit}
pre.format {font-family: inherit}
pre.menu-comment {font-family: serif}
pre.menu-preformatted {font-family: serif}
pre.smalldisplay {font-family: inherit; font-size: smaller}
pre.smallexample {font-size: smaller}
pre.smallformat {font-family: inherit; font-size: smaller}
pre.smalllisp {font-size: smaller}
span.nocodebreak {white-space:nowrap}
span.nolinebreak {white-space:nowrap}
span.roman {font-family:serif; font-weight:normal}
span.sansserif {font-family:sans-serif; font-weight:normal}
ul.no-bullet {list-style: none}
-->
</style>
</head>
<body lang="en" bgcolor="#FFFFFF" text="#000000" link="#0000FF" vlink="#800080" alink="#FF0000">
<a name="QR-Decomposition-with-Column-Pivoting"></a>
<div class="header">
<p>
Next: <a href="Complete-Orthogonal-Decomposition.html#Complete-Orthogonal-Decomposition" accesskey="n" rel="next">Complete Orthogonal Decomposition</a>, Previous: <a href="QR-Decomposition.html#QR-Decomposition" accesskey="p" rel="previous">QR Decomposition</a>, Up: <a href="Linear-Algebra.html#Linear-Algebra" accesskey="u" rel="up">Linear Algebra</a> [<a href="Function-Index.html#Function-Index" title="Index" rel="index">Index</a>]</p>
</div>
<hr>
<a name="QR-Decomposition-with-Column-Pivoting-1"></a>
<h3 class="section">14.3 QR Decomposition with Column Pivoting</h3>
<a name="index-QR-decomposition-with-column-pivoting"></a>
<p>The <em>QR</em> decomposition of an <em>M</em>-by-<em>N</em> matrix <em>A</em>
can be extended to the rank deficient case by introducing a column permutation <em>P</em>,
</p>
<div class="example">
<pre class="example">A P = Q R
</pre></div>
<p>The first <em>r</em> columns of <em>Q</em> form an orthonormal basis
for the range of <em>A</em> for a matrix with column rank <em>r</em>. This
decomposition can also be used to convert the linear system <em>A x =
b</em> into the triangular system <em>R y = Q^T b, x = P y</em>, which can be
solved by back-substitution and permutation. We denote the <em>QR</em>
decomposition with column pivoting by <em>QRP^T</em> since <em>A = Q R
P^T</em>. When <em>A</em> is rank deficient with <em>r = {\rm rank}(A)</em>, the matrix
<em>R</em> can be partitioned as
</p>
<div class="example">
<pre class="example">R = [ R11 R12; 0 R22 ] =~ [ R11 R12; 0 0 ]
</pre></div>
<p>where <em>R_{11}</em> is <em>r</em>-by-<em>r</em> and nonsingular. In this case,
a “basic” least squares solution for the overdetermined system <em>A x = b</em>
can be obtained as
</p>
<div class="example">
<pre class="example">x = P [ R11^-1 c1 ; 0 ]
</pre></div>
<p>where <em>c_1</em> consists of the first <em>r</em> elements of <em>Q^T b</em>.
This basic solution is not guaranteed to be the minimum norm solution unless
<em>R_{12} = 0</em> (see <a href="Complete-Orthogonal-Decomposition.html#Complete-Orthogonal-Decomposition">Complete Orthogonal Decomposition</a>).
</p>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005fdecomp"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_decomp</strong> <em>(gsl_matrix * <var>A</var>, gsl_vector * <var>tau</var>, gsl_permutation * <var>p</var>, int * <var>signum</var>, gsl_vector * <var>norm</var>)</em></dt>
<dd><p>This function factorizes the <em>M</em>-by-<em>N</em> matrix <var>A</var> into
the <em>QRP^T</em> decomposition <em>A = Q R P^T</em>. On output the
diagonal and upper triangular part of the input matrix contain the
matrix <em>R</em>. The permutation matrix <em>P</em> is stored in the
permutation <var>p</var>. The sign of the permutation is given by
<var>signum</var>. It has the value <em>(-1)^n</em>, where <em>n</em> is the
number of interchanges in the permutation. The vector <var>tau</var> and the
columns of the lower triangular part of the matrix <var>A</var> contain the
Householder coefficients and vectors which encode the orthogonal matrix
<var>Q</var>. The vector <var>tau</var> must be of length <em>k=\min(M,N)</em>. The
matrix <em>Q</em> is related to these components by, <em>Q = Q_k ... Q_2
Q_1</em> where <em>Q_i = I - \tau_i v_i v_i^T</em> and <em>v_i</em> is the
Householder vector <em>v_i =
(0,...,1,A(i+1,i),A(i+2,i),...,A(m,i))</em>. This is the same storage scheme
as used by <small>LAPACK</small>. The vector <var>norm</var> is a workspace of length
<var>N</var> used for column pivoting.
</p>
<p>The algorithm used to perform the decomposition is Householder QR with
column pivoting (Golub & Van Loan, <cite>Matrix Computations</cite>, Algorithm
5.4.1).
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005fdecomp2"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_decomp2</strong> <em>(const gsl_matrix * <var>A</var>, gsl_matrix * <var>q</var>, gsl_matrix * <var>r</var>, gsl_vector * <var>tau</var>, gsl_permutation * <var>p</var>, int * <var>signum</var>, gsl_vector * <var>norm</var>)</em></dt>
<dd><p>This function factorizes the matrix <var>A</var> into the decomposition
<em>A = Q R P^T</em> without modifying <var>A</var> itself and storing the
output in the separate matrices <var>q</var> and <var>r</var>.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005fsolve"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_solve</strong> <em>(const gsl_matrix * <var>QR</var>, const gsl_vector * <var>tau</var>, const gsl_permutation * <var>p</var>, const gsl_vector * <var>b</var>, gsl_vector * <var>x</var>)</em></dt>
<dd><p>This function solves the square system <em>A x = b</em> using the <em>QRP^T</em>
decomposition of <em>A</em> held in (<var>QR</var>, <var>tau</var>, <var>p</var>) which must
have been computed previously by <code>gsl_linalg_QRPT_decomp</code>.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005fsvx"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_svx</strong> <em>(const gsl_matrix * <var>QR</var>, const gsl_vector * <var>tau</var>, const gsl_permutation * <var>p</var>, gsl_vector * <var>x</var>)</em></dt>
<dd><p>This function solves the square system <em>A x = b</em> in-place using the
<em>QRP^T</em> decomposition of <em>A</em> held in
(<var>QR</var>,<var>tau</var>,<var>p</var>). On input <var>x</var> should contain the
right-hand side <em>b</em>, which is replaced by the solution on output.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005flssolve"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_lssolve</strong> <em>(const gsl_matrix * <var>QR</var>, const gsl_vector * <var>tau</var>, const gsl_permutation * <var>p</var>, const gsl_vector * <var>b</var>, gsl_vector * <var>x</var>, gsl_vector * <var>residual</var>)</em></dt>
<dd><p>This function finds the least squares solution to the overdetermined
system <em>A x = b</em> where the matrix <var>A</var> has more rows than
columns and is assumed to have full rank. The least squares solution minimizes
the Euclidean norm of the residual, <em>||b - A x||</em>. The routine requires as input
the <em>QR</em> decomposition of <em>A</em> into (<var>QR</var>, <var>tau</var>, <var>p</var>) given by
<code>gsl_linalg_QRPT_decomp</code>. The solution is returned in <var>x</var>. The
residual is computed as a by-product and stored in <var>residual</var>. For rank
deficient matrices, <code>gsl_linalg_QRPT_lssolve2</code> should be used instead.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005flssolve2"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_lssolve2</strong> <em>(const gsl_matrix * <var>QR</var>, const gsl_vector * <var>tau</var>, const gsl_permutation * <var>p</var>, const gsl_vector * <var>b</var>, const size_t <var>rank</var>, gsl_vector * <var>x</var>, gsl_vector * <var>residual</var>)</em></dt>
<dd><p>This function finds the least squares solution to the overdetermined
system <em>A x = b</em> where the matrix <var>A</var> has more rows than
columns and has rank given by the input <var>rank</var>. If the user does not
know the rank of <em>A</em>, the routine <code>gsl_linalg_QRPT_rank</code> can be
called to estimate it. The least squares solution is
the so-called “basic” solution discussed above and may not be the minimum
norm solution. The routine requires as input
the <em>QR</em> decomposition of <em>A</em> into (<var>QR</var>, <var>tau</var>, <var>p</var>) given by
<code>gsl_linalg_QRPT_decomp</code>. The solution is returned in <var>x</var>. The
residual is computed as a by-product and stored in <var>residual</var>.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005fQRsolve"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_QRsolve</strong> <em>(const gsl_matrix * <var>Q</var>, const gsl_matrix * <var>R</var>, const gsl_permutation * <var>p</var>, const gsl_vector * <var>b</var>, gsl_vector * <var>x</var>)</em></dt>
<dd><p>This function solves the square system <em>R P^T x = Q^T b</em> for
<var>x</var>. It can be used when the <em>QR</em> decomposition of a matrix is
available in unpacked form as (<var>Q</var>, <var>R</var>).
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005fupdate"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_update</strong> <em>(gsl_matrix * <var>Q</var>, gsl_matrix * <var>R</var>, const gsl_permutation * <var>p</var>, gsl_vector * <var>w</var>, const gsl_vector * <var>v</var>)</em></dt>
<dd><p>This function performs a rank-1 update <em>w v^T</em> of the <em>QRP^T</em>
decomposition (<var>Q</var>, <var>R</var>, <var>p</var>). The update is given by
<em>Q'R' = Q (R + w v^T P)</em> where the output matrices <em>Q'</em> and
<em>R'</em> are also orthogonal and right triangular. Note that <var>w</var> is
destroyed by the update. The permutation <var>p</var> is not changed.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005fRsolve"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_Rsolve</strong> <em>(const gsl_matrix * <var>QR</var>, const gsl_permutation * <var>p</var>, const gsl_vector * <var>b</var>, gsl_vector * <var>x</var>)</em></dt>
<dd><p>This function solves the triangular system <em>R P^T x = b</em> for the
<em>N</em>-by-<em>N</em> matrix <em>R</em> contained in <var>QR</var>.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005fRsvx"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_Rsvx</strong> <em>(const gsl_matrix * <var>QR</var>, const gsl_permutation * <var>p</var>, gsl_vector * <var>x</var>)</em></dt>
<dd><p>This function solves the triangular system <em>R P^T x = b</em> in-place
for the <em>N</em>-by-<em>N</em> matrix <em>R</em> contained in <var>QR</var>. On
input <var>x</var> should contain the right-hand side <em>b</em>, which is
replaced by the solution on output.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005frank"></a>Function: <em>size_t</em> <strong>gsl_linalg_QRPT_rank</strong> <em>(const gsl_matrix * <var>QR</var>, const double <var>tol</var>)</em></dt>
<dd><p>This function estimates the rank of the triangular matrix <em>R</em> contained in <var>QR</var>.
The algorithm simply counts the number of diagonal elements of <em>R</em> whose absolute value
is greater than the specified tolerance <var>tol</var>. If the input <var>tol</var> is negative,
a default value of <em>20 (M + N) eps(max(|diag(R)|))</em> is used.
</p></dd></dl>
<dl>
<dt><a name="index-gsl_005flinalg_005fQRPT_005frcond"></a>Function: <em>int</em> <strong>gsl_linalg_QRPT_rcond</strong> <em>(const gsl_matrix * <var>QR</var>, double * <var>rcond</var>, gsl_vector * <var>work</var>)</em></dt>
<dd><p>This function estimates the reciprocal condition number (using the 1-norm) of the <em>R</em> factor,
stored in the upper triangle of <var>QR</var>. The reciprocal condition number estimate, defined as
<em>1 / (||R||_1 \cdot ||R^{-1}||_1)</em>, is stored in <var>rcond</var>.
Additional workspace of size <em>3 N</em> is required in <var>work</var>.
</p></dd></dl>
<hr>
<div class="header">
<p>
Next: <a href="Complete-Orthogonal-Decomposition.html#Complete-Orthogonal-Decomposition" accesskey="n" rel="next">Complete Orthogonal Decomposition</a>, Previous: <a href="QR-Decomposition.html#QR-Decomposition" accesskey="p" rel="previous">QR Decomposition</a>, Up: <a href="Linear-Algebra.html#Linear-Algebra" accesskey="u" rel="up">Linear Algebra</a> [<a href="Function-Index.html#Function-Index" title="Index" rel="index">Index</a>]</p>
</div>
</body>
</html>
|
