Download LU Factorization: Solving Matrix Equations using Lower and Upper Triangular Matrices - Pro and more Study notes Mathematical Methods for Numerical Analysis and Optimization in PDF only on Docsity!
Introduction to Scientific Computing^1 Math 3176 - Spring 2009 Instructor: Alexandros Sopasakis
Direct methods and matrix factorization
The main problem which is encountered in several research applications relating to matrix algebra is AX = B (1)
where A and B are known matrices while X is the matrix of the unknowns. A variety of problems can be transformed into this set up and therefore being able to solve such a system effectively is very important. There are several techniques that can be applied in solving (1). When the size of the matrices is sufficiently large the idea of simply finding the inverse of matrix A (if it exists) is no longer attractive for several reasons. Simply though the computer resources required for such a calculation are great! Thus alternate more effective techniques are used in practice. In early algebra classes we have already encountered elimination and pivoting methods such as Gauss-Jordan method or Gaussian elimination and backward substitution. These methods require n^3 operations for matrices of size n × n. Thus as the size of the matrix increases the cost in operation skyrockets. If instead we produce a factorized version of A then we can reduce the operational cost of solving for X from n^3 to n^2. This translates to almost 99% reduction in calculations assuming that the matrices are larger than 100×100 (not unusual for applications nowadays). However to factor out A we do require n^3 operation anyway. The only benefit though is that once the matrix A has been factored then we can solve any general system AX = B with n^2 operations only.
LU factorization - Doolittle’s version
In the first method which we examine we will factor matrix A into two other matrices. A lower triangular matrix L and an upper triangular matrix U such that
A = LU
The idea for solving system (1) will be as follows. Starting from (1)
AX = B
we let a new variable Y be, Y = U X so that,
AX = LU X = B becomes LY = B (^1) DEPARTMENT OF MATHEMATICS AND STATISTICS, University OF North Carolina at Charlotte
Since L is a lower triangular matrix this system is almost trivial to solve for the unknown Y ’s. Once we find all the values for Y then we can start solving the system U X = Y
Note that this system is also very easy to solve since U is an upper triangular matrix. Thus finding X with this method is no big deal. The only think left to do is actually find the matrices L and U with the required properties for which A = LU. This is accomplished by the usual Gaussian elimination method which is applied only up to the point of obtaining an upper triangular matrix (without the back substitution). Let us look at a simple example: Example Solve the following matrix system using an LU factorization.
x 1 x 2 x 3
Solution The main part will be to produce the LU factorization. Once this is done then solving the system will be easy. To produce the factorization we start with the usual Gaussian elimination method. For ease in notation we denote by R each row of A. Then to create zeros below the element a(1, 1) we simply do the following − 3 R 1 + R 2 → R 2 , − 1 R 1 + R 3 → R 3 ,
Last we create zero below a(2, 2) via − 1 R 2 + R 3 → R 3 ,
This procedure remarkably has already produced our required matrices L and U from A. In fact the matrices are,
A = LU =
Do the multiplication to check the result! Note that U is simply the resulting matrix from our Gaussian elimination while L is consisting of 1’s in the diagonal while the
- First solve LY = B.
- For i = 1, 2 ,... , n do
b(i) −
∑i− 1 j=1 l(i, j)y(j)
/l(i, i)
- Now solve U X = Y by back substitution in exactly the same way,
- For i = n, n − 1 ,... , 1 do
y(i) −
j=n u(i, j)x(j)
/u(i, i)
There are a couple of results which are interesting since they give us an incite at when these methods work. The following definition is necessary first,
Definition 1. The n × n matrix A is said to be strictly diagonally dominant if
|a(i, i)| >
∑^ n
j 6 =i
|a(i, j)| for all i = 1, 2 ,... , n
The results comes indirectly through Gaussian elimination:
Theorem 2. A strictly diagonally dominant matrix A is non-singular. Furthermore Gaussian elimination can be performed on any linear system of the form Ax = B to obtain its unique solution without row or column interchanges, and the computations are stable with respect to the growth of round-off errors.
When can we perform LU decomposition? The following theorem gives the answer,
Theorem 3. If Gaussian elimination can be performed on the linear system AX = B without row interchanges then the matrix A can be factored into the product of a lower-triangular matrix L and an upper-triangular matrix U , where A = LU.
More relevant results follow in the next section. There is another type of factorization which is in fact very similar to this LU or Doolittle’s decomposition. The alternate factorization method also produces an LU decomposition where U being a unit upper triangular matrix instead of L. This is called Crout’s factorization. Naturally either factorization will do the job and producing one or the other is a matter of taste than anything else. You can change the provided pseudo-code very easily in order to produce such a factorization.
LDLT^ and LLT^ or Choleski’s factorization
We continue here by presenting more methods for factoring A. All the techniques presented, similarly to the LU decomposition of A are of the same operation cost of n^3 overall. As the name denotes an LDLT^ type factorization takes the following form,
A = LDLT
where L as usual is lower triangular and D is a diagonal matrix with positive entries in the diagonal. Similarly the Choleski factorization
A = LLT
consists of a lower and upper triangular matrix where neither of which have 1’s in the diagonal (in contrast to either Doolittle’s or Crout’s factorizations). It is very easily to construct any of the above factorizations once you have an LU decomposition of A. Let us look at the equivalent factorizations for the following matrix
A =
Using our pseudo-code we obtain the following LU decomposition of A,
A = LU =
Now the equivalent LDLT^ decomposition consist of the following three matrices,
A = LDLT^ =
Note how the new upper triangular matrix has been obtained by simply dividing each row of the old upper triangular matrix with the respective diagonal element. Now that we have the LDLT^ factorization we can also easily obtain the equivalent Crout’s factorization,
A =
Note here that the new lower triangular matrix is constructed by simply multiplying out matrices L and D. Last the Choleski decomposition is also easily constructed
