3
<title>LU 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="next" href="QR-Decomposition.html#QR-Decomposition" title="QR Decomposition">
10
<link href="http://www.gnu.org/software/texinfo/" rel="generator-home" title="Texinfo Homepage">
12
Copyright (C) 1996, 1997, 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 The GSL Team.
14
Permission is granted to copy, distribute and/or modify this document
15
under the terms of the GNU Free Documentation License, Version 1.2 or
16
any later version published by the Free Software Foundation; with the
17
Invariant Sections being ``GNU General Public License'' and ``Free Software
18
Needs Free Documentation'', the Front-Cover text being ``A GNU Manual'',
19
and with the Back-Cover Text being (a) (see below). A copy of the
20
license is included in the section entitled ``GNU Free Documentation
23
(a) The Back-Cover Text is: ``You have freedom to copy and modify this
24
GNU Manual, like GNU software.''-->
25
<meta http-equiv="Content-Style-Type" content="text/css">
26
<style type="text/css"><!--
27
pre.display { font-family:inherit }
28
pre.format { font-family:inherit }
29
pre.smalldisplay { font-family:inherit; font-size:smaller }
30
pre.smallformat { font-family:inherit; font-size:smaller }
31
pre.smallexample { font-size:smaller }
32
pre.smalllisp { font-size:smaller }
33
span.sc { font-variant:small-caps }
34
span.roman { font-family:serif; font-weight:normal; }
35
span.sansserif { font-family:sans-serif; font-weight:normal; }
41
<a name="LU-Decomposition"></a>
42
Next: <a rel="next" accesskey="n" href="QR-Decomposition.html#QR-Decomposition">QR Decomposition</a>,
43
Up: <a rel="up" accesskey="u" href="Linear-Algebra.html#Linear-Algebra">Linear Algebra</a>
47
<h3 class="section">13.1 LU Decomposition</h3>
49
<p><a name="index-LU-decomposition-1182"></a>
50
A general square matrix A has an LU decomposition into
51
upper and lower triangular matrices,
52
where P is a permutation matrix, L is unit lower
53
triangular matrix and U is upper triangular matrix. For square
54
matrices this decomposition can be used to convert the linear system
55
A x = b into a pair of triangular systems (L y = P b,
56
U x = y), which can be solved by forward and back-substitution.
57
Note that the LU decomposition is valid for singular matrices.
60
— Function: int <b>gsl_linalg_LU_decomp</b> (<var>gsl_matrix * A, gsl_permutation * p, int * signum</var>)<var><a name="index-gsl_005flinalg_005fLU_005fdecomp-1183"></a></var><br>
61
— Function: int <b>gsl_linalg_complex_LU_decomp</b> (<var>gsl_matrix_complex * A, gsl_permutation * p, int * signum</var>)<var><a name="index-gsl_005flinalg_005fcomplex_005fLU_005fdecomp-1184"></a></var><br>
62
<blockquote><p>These functions factorize the square matrix <var>A</var> into the LU
63
decomposition PA = LU. On output the diagonal and upper
64
triangular part of the input matrix <var>A</var> contain the matrix
65
U. The lower triangular part of the input matrix (excluding the
66
diagonal) contains L. The diagonal elements of L are
67
unity, and are not stored.
69
<p>The permutation matrix P is encoded in the permutation
70
<var>p</var>. The j-th column of the matrix P is given by the
71
k-th column of the identity matrix, where k = p_j the
72
j-th element of the permutation vector. The sign of the
73
permutation is given by <var>signum</var>. It has the value (-1)^n,
74
where n is the number of interchanges in the permutation.
76
<p>The algorithm used in the decomposition is Gaussian Elimination with
77
partial pivoting (Golub & Van Loan, <cite>Matrix Computations</cite>,
79
</p></blockquote></div>
81
<p><a name="index-linear-systems_002c-solution-of-1185"></a>
84
— Function: int <b>gsl_linalg_LU_solve</b> (<var>const gsl_matrix * LU, const gsl_permutation * p, const gsl_vector * b, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fLU_005fsolve-1186"></a></var><br>
85
— Function: int <b>gsl_linalg_complex_LU_solve</b> (<var>const gsl_matrix_complex * LU, const gsl_permutation * p, const gsl_vector_complex * b, gsl_vector_complex * x</var>)<var><a name="index-gsl_005flinalg_005fcomplex_005fLU_005fsolve-1187"></a></var><br>
86
<blockquote><p>These functions solve the square system A x = b using the LU
87
decomposition of A into (<var>LU</var>, <var>p</var>) given by
88
<code>gsl_linalg_LU_decomp</code> or <code>gsl_linalg_complex_LU_decomp</code>.
89
</p></blockquote></div>
92
— Function: int <b>gsl_linalg_LU_svx</b> (<var>const gsl_matrix * LU, const gsl_permutation * p, gsl_vector * x</var>)<var><a name="index-gsl_005flinalg_005fLU_005fsvx-1188"></a></var><br>
93
— Function: int <b>gsl_linalg_complex_LU_svx</b> (<var>const gsl_matrix_complex * LU, const gsl_permutation * p, gsl_vector_complex * x</var>)<var><a name="index-gsl_005flinalg_005fcomplex_005fLU_005fsvx-1189"></a></var><br>
94
<blockquote><p>These functions solve the square system A x = b in-place using the
95
LU decomposition of A into (<var>LU</var>,<var>p</var>). On input
96
<var>x</var> should contain the right-hand side b, which is replaced
97
by the solution on output.
98
</p></blockquote></div>
100
<p><a name="index-refinement-of-solutions-in-linear-systems-1190"></a><a name="index-iterative-refinement-of-solutions-in-linear-systems-1191"></a><a name="index-linear-systems_002c-refinement-of-solutions-1192"></a>
103
— Function: int <b>gsl_linalg_LU_refine</b> (<var>const gsl_matrix * A, const gsl_matrix * LU, const gsl_permutation * p, const gsl_vector * b, gsl_vector * x, gsl_vector * residual</var>)<var><a name="index-gsl_005flinalg_005fLU_005frefine-1193"></a></var><br>
104
— Function: int <b>gsl_linalg_complex_LU_refine</b> (<var>const gsl_matrix_complex * A, const gsl_matrix_complex * LU, const gsl_permutation * p, const gsl_vector_complex * b, gsl_vector_complex * x, gsl_vector_complex * residual</var>)<var><a name="index-gsl_005flinalg_005fcomplex_005fLU_005frefine-1194"></a></var><br>
105
<blockquote><p>These functions apply an iterative improvement to <var>x</var>, the solution
106
of A x = b, using the LU decomposition of A into
107
(<var>LU</var>,<var>p</var>). The initial residual r = A x - b is also
108
computed and stored in <var>residual</var>.
109
</p></blockquote></div>
111
<p><a name="index-inverse-of-a-matrix_002c-by-LU-decomposition-1195"></a><a name="index-matrix-inverse-1196"></a>
114
— Function: int <b>gsl_linalg_LU_invert</b> (<var>const gsl_matrix * LU, const gsl_permutation * p, gsl_matrix * inverse</var>)<var><a name="index-gsl_005flinalg_005fLU_005finvert-1197"></a></var><br>
115
— Function: int <b>gsl_linalg_complex_LU_invert</b> (<var>const gsl_matrix_complex * LU, const gsl_permutation * p, gsl_matrix_complex * inverse</var>)<var><a name="index-gsl_005flinalg_005fcomplex_005fLU_005finvert-1198"></a></var><br>
116
<blockquote><p>These functions compute the inverse of a matrix A from its
117
LU decomposition (<var>LU</var>,<var>p</var>), storing the result in the
118
matrix <var>inverse</var>. The inverse is computed by solving the system
119
A x = b for each column of the identity matrix. It is preferable
120
to avoid direct use of the inverse whenever possible, as the linear
121
solver functions can obtain the same result more efficiently and
122
reliably (consult any introductory textbook on numerical linear algebra
124
</p></blockquote></div>
126
<p><a name="index-determinant-of-a-matrix_002c-by-LU-decomposition-1199"></a><a name="index-matrix-determinant-1200"></a>
129
— Function: double <b>gsl_linalg_LU_det</b> (<var>gsl_matrix * LU, int signum</var>)<var><a name="index-gsl_005flinalg_005fLU_005fdet-1201"></a></var><br>
130
— Function: gsl_complex <b>gsl_linalg_complex_LU_det</b> (<var>gsl_matrix_complex * LU, int signum</var>)<var><a name="index-gsl_005flinalg_005fcomplex_005fLU_005fdet-1202"></a></var><br>
131
<blockquote><p>These functions compute the determinant of a matrix A from its
132
LU decomposition, <var>LU</var>. The determinant is computed as the
133
product of the diagonal elements of U and the sign of the row
134
permutation <var>signum</var>.
135
</p></blockquote></div>
137
<p><a name="index-logarithm-of-the-determinant-of-a-matrix-1203"></a>
140
— Function: double <b>gsl_linalg_LU_lndet</b> (<var>gsl_matrix * LU</var>)<var><a name="index-gsl_005flinalg_005fLU_005flndet-1204"></a></var><br>
141
— Function: double <b>gsl_linalg_complex_LU_lndet</b> (<var>gsl_matrix_complex * LU</var>)<var><a name="index-gsl_005flinalg_005fcomplex_005fLU_005flndet-1205"></a></var><br>
142
<blockquote><p>These functions compute the logarithm of the absolute value of the
143
determinant of a matrix A, \ln|\det(A)|, from its LU
144
decomposition, <var>LU</var>. This function may be useful if the direct
145
computation of the determinant would overflow or underflow.
146
</p></blockquote></div>
148
<p><a name="index-sign-of-the-determinant-of-a-matrix-1206"></a>
151
— Function: int <b>gsl_linalg_LU_sgndet</b> (<var>gsl_matrix * LU, int signum</var>)<var><a name="index-gsl_005flinalg_005fLU_005fsgndet-1207"></a></var><br>
152
— Function: gsl_complex <b>gsl_linalg_complex_LU_sgndet</b> (<var>gsl_matrix_complex * LU, int signum</var>)<var><a name="index-gsl_005flinalg_005fcomplex_005fLU_005fsgndet-1208"></a></var><br>
153
<blockquote><p>These functions compute the sign or phase factor of the determinant of a
154
matrix A, \det(A)/|\det(A)|, from its LU decomposition,
156
</p></blockquote></div>