3
<title>QR Decomposition - 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="LU-Decomposition.html#LU-Decomposition" title="LU Decomposition">
10
<link rel="next" href="QR-Decomposition-with-Column-Pivoting.html#QR-Decomposition-with-Column-Pivoting" title="QR Decomposition with Column Pivoting">
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"></a>
43
Next: <a rel="next" accesskey="n" href="QR-Decomposition-with-Column-Pivoting.html#QR-Decomposition-with-Column-Pivoting">QR Decomposition with Column Pivoting</a>,
44
Previous: <a rel="previous" accesskey="p" href="LU-Decomposition.html#LU-Decomposition">LU Decomposition</a>,
45
Up: <a rel="up" accesskey="u" href="Linear-Algebra.html#Linear-Algebra">Linear Algebra</a>
49
<h3 class="section">13.2 QR Decomposition</h3>
51
<p><a name="index-QR-decomposition-1209"></a>
52
A general rectangular M-by-N matrix A has a
53
QR decomposition into the product of an orthogonal
54
M-by-M square matrix Q (where Q^T Q = I) and
55
an M-by-N right-triangular matrix R,
56
This decomposition can be used to convert the linear system A x =
57
b into the triangular system R x = Q^T b, which can be solved by
58
back-substitution. Another use of the QR decomposition is to
59
compute an orthonormal basis for a set of vectors. The first N
60
columns of Q form an orthonormal basis for the range of A,
61
ran(A), when A has full column rank.
64
— Function: int <b>gsl_linalg_QR_decomp</b> (<var>gsl_matrix * A, gsl_vector * tau</var>)<var><a name="index-gsl_005flinalg_005fQR_005fdecomp-1210"></a></var><br>
65
<blockquote><p>This function factorizes the M-by-N matrix <var>A</var> into
66
the QR decomposition A = Q R. On output the diagonal and
67
upper triangular part of the input matrix contain the matrix
68
R. The vector <var>tau</var> and the columns of the lower triangular
69
part of the matrix <var>A</var> contain the Householder coefficients and
70
Householder vectors which encode the orthogonal matrix <var>Q</var>. The
71
vector <var>tau</var> must be of length k=\min(M,N). The matrix
72
Q is related to these components by, Q = Q_k ... Q_2 Q_1
73
where Q_i = I - \tau_i v_i v_i^T and v_i is the
74
Householder vector v_i =
75
(0,...,1,A(i+1,i),A(i+2,i),...,A(m,i)). This is the same storage scheme
76
as used by <span class="sc">lapack</span>.
78
<p>The algorithm used to perform the decomposition is Householder QR (Golub
79
& Van Loan, <cite>Matrix Computations</cite>, Algorithm 5.2.1).
80
</p></blockquote></div>
83
— Function: int <b>gsl_linalg_QR_solve</b> (<var>const gsl_matrix * QR, const gsl_vector * tau, const gsl_vector * b, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQR_005fsolve-1211"></a></var><br>
84
<blockquote><p>This function solves the square system A x = b using the QR
85
decomposition of A into (<var>QR</var>, <var>tau</var>) given by
86
<code>gsl_linalg_QR_decomp</code>. The least-squares solution for rectangular systems can
87
be found using <code>gsl_linalg_QR_lssolve</code>.
88
</p></blockquote></div>
91
— Function: int <b>gsl_linalg_QR_svx</b> (<var>const gsl_matrix * QR, const gsl_vector * tau, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQR_005fsvx-1212"></a></var><br>
92
<blockquote><p>This function solves the square system A x = b in-place using the
93
QR decomposition of A into (<var>QR</var>,<var>tau</var>) given by
94
<code>gsl_linalg_QR_decomp</code>. On input <var>x</var> should contain the
95
right-hand side b, which is replaced by the solution on output.
96
</p></blockquote></div>
99
— Function: int <b>gsl_linalg_QR_lssolve</b> (<var>const gsl_matrix * QR, const gsl_vector * tau, const gsl_vector * b, gsl_vector * x, gsl_vector * residual</var>)<var><a name="index-gsl_005flinalg_005fQR_005flssolve-1213"></a></var><br>
100
<blockquote><p>This function finds the least squares solution to the overdetermined
101
system A x = b where the matrix <var>A</var> has more rows than
102
columns. The least squares solution minimizes the Euclidean norm of the
103
residual, ||Ax - b||.The routine uses the QR decomposition
104
of A into (<var>QR</var>, <var>tau</var>) given by
105
<code>gsl_linalg_QR_decomp</code>. The solution is returned in <var>x</var>. The
106
residual is computed as a by-product and stored in <var>residual</var>.
107
</p></blockquote></div>
110
— Function: int <b>gsl_linalg_QR_QTvec</b> (<var>const gsl_matrix * QR, const gsl_vector * tau, gsl_vector * v</var>)<var><a name="index-gsl_005flinalg_005fQR_005fQTvec-1214"></a></var><br>
111
<blockquote><p>This function applies the matrix Q^T encoded in the decomposition
112
(<var>QR</var>,<var>tau</var>) to the vector <var>v</var>, storing the result Q^T
113
v in <var>v</var>. The matrix multiplication is carried out directly using
114
the encoding of the Householder vectors without needing to form the full
116
</p></blockquote></div>
119
— Function: int <b>gsl_linalg_QR_Qvec</b> (<var>const gsl_matrix * QR, const gsl_vector * tau, gsl_vector * v</var>)<var><a name="index-gsl_005flinalg_005fQR_005fQvec-1215"></a></var><br>
120
<blockquote><p>This function applies the matrix Q encoded in the decomposition
121
(<var>QR</var>,<var>tau</var>) to the vector <var>v</var>, storing the result Q
122
v in <var>v</var>. The matrix multiplication is carried out directly using
123
the encoding of the Householder vectors without needing to form the full
125
</p></blockquote></div>
128
— Function: int <b>gsl_linalg_QR_Rsolve</b> (<var>const gsl_matrix * QR, const gsl_vector * b, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQR_005fRsolve-1216"></a></var><br>
129
<blockquote><p>This function solves the triangular system R x = b for
130
<var>x</var>. It may be useful if the product b' = Q^T b has already
131
been computed using <code>gsl_linalg_QR_QTvec</code>.
132
</p></blockquote></div>
135
— Function: int <b>gsl_linalg_QR_Rsvx</b> (<var>const gsl_matrix * QR, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQR_005fRsvx-1217"></a></var><br>
136
<blockquote><p>This function solves the triangular system R x = b for <var>x</var>
137
in-place. On input <var>x</var> should contain the right-hand side b
138
and is replaced by the solution on output. This function may be useful if
139
the product b' = Q^T b has already been computed using
140
<code>gsl_linalg_QR_QTvec</code>.
141
</p></blockquote></div>
144
— Function: int <b>gsl_linalg_QR_unpack</b> (<var>const gsl_matrix * QR, const gsl_vector * tau, gsl_matrix * Q, gsl_matrix * R</var>)<var><a name="index-gsl_005flinalg_005fQR_005funpack-1218"></a></var><br>
145
<blockquote><p>This function unpacks the encoded QR decomposition
146
(<var>QR</var>,<var>tau</var>) into the matrices <var>Q</var> and <var>R</var>, where
147
<var>Q</var> is M-by-M and <var>R</var> is M-by-N.
148
</p></blockquote></div>
151
— Function: int <b>gsl_linalg_QR_QRsolve</b> (<var>gsl_matrix * Q, gsl_matrix * R, const gsl_vector * b, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fQR_005fQRsolve-1219"></a></var><br>
152
<blockquote><p>This function solves the system R x = Q^T b for <var>x</var>. It can
153
be used when the QR decomposition of a matrix is available in
154
unpacked form as (<var>Q</var>, <var>R</var>).
155
</p></blockquote></div>
158
— Function: int <b>gsl_linalg_QR_update</b> (<var>gsl_matrix * Q, gsl_matrix * R, gsl_vector * w, const gsl_vector * v</var>)<var><a name="index-gsl_005flinalg_005fQR_005fupdate-1220"></a></var><br>
159
<blockquote><p>This function performs a rank-1 update w v^T of the QR
160
decomposition (<var>Q</var>, <var>R</var>). The update is given by Q'R' = Q
161
R + w v^T where the output matrices Q' and R' are also
162
orthogonal and right triangular. Note that <var>w</var> is destroyed by the
164
</p></blockquote></div>
167
— Function: int <b>gsl_linalg_R_solve</b> (<var>const gsl_matrix * R, const gsl_vector * b, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fR_005fsolve-1221"></a></var><br>
168
<blockquote><p>This function solves the triangular system R x = b for the
169
N-by-N matrix <var>R</var>.
170
</p></blockquote></div>
173
— Function: int <b>gsl_linalg_R_svx</b> (<var>const gsl_matrix * R, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fR_005fsvx-1222"></a></var><br>
174
<blockquote><p>This function solves the triangular system R x = b in-place. On
175
input <var>x</var> should contain the right-hand side b, which is
176
replaced by the solution on output.
177
</p></blockquote></div>