lecture7 - Numerical Methods in Engineering ENGR 391 System...

Info iconThis preview shows page 1. Sign up to view the full content.

View Full Document Right Arrow Icon
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: Numerical Methods in Engineering ENGR 391 System of linear equations Faculty of Engineering and Computer Sciences Concordia University LU decomposition In Gauss elimination, we consider an augmented matrix combining [A] and {b} to solve the system a11 a12 a13 b1 a21 a22 a23 b2 a31 a32 a33 b3 Gauss elimination becomes inefficient when solving equations with the same coefficients for [A] but with different b’s. LU decomposition separates the time consuming elimination of [A] from the manipulation of {b} Hence, the decomposed [A] could be used with several {b} ’s in an efficient manner. What is LU decomposition? LU decomposition LU decomposition is based on the fact that any square matrix [A] can be written as a product of two matrices as: [A]=[L][U] where [L] is a lower triangular matrix and [U] is an upper triangular matrix. In other words, [A] is factored or “decomposed into lower [L] and upper [U] triangular matrices. a11 a12 a13 100 . . . a21 a22 a23 .10 0 . . a31 a32 a33 . 00. . 1 How can we get the LU decomposition of a Matrix? LU decomposition from Gauss elimination Gauss algorithm: replaces Ax b 0 Ux d 0 by where U is an upper diagonal matrix. Operator form: what is the matrix T which does the same ? T Ax b Ux d Gaussian Elimination a11 A a12 a13 a 21 a 22 a 23 a31 The matrix a32 a33 a11 is replaced by: ~ A a12 a13 a21 ca11 a22 ca12 a31 a32 a23 ca13 a33 We can write this operation in the form of ~ A T21 c A where T21 c is a matrix The Gauss operator • A Gauss elimination step can be written as: 1 a12 a13 c 1 0 a21 a22 a23 0 0 0 a11 0 1 a31 a32 We have found a matrix T21 a33 a11 a12 a21 ca11 a22 ca12 a31 a32 a13 a23 ca13 a33 c which, when multiplying the matrix does the operation subtracting c times lines one from line 2 A , Gaussian Elimination Let us apply the theorem to a 3X3 matrix. First step of the Gauss elimination algorithm is the elimination of the a21 coefficient 1 a21 a11 0 0 0 a11 a12 a13 1 0 a21 a22 a23 0 1 a31 a32 a33 a21 a11 a12 a21 a21 a11 a22 a12 a11 a11 a31 a32 a11 a12 a13 0 ~ a22 ~ a23 a31 a32 a33 a23 a13 a21 a13 a11 a33 Gaussian Elimination Next step is the elimination of the a31 coefficient. 1 0 0 a11 a12 a13 0 a31 a11 10 0 ~ a22 ~ a23 0 1 a31 a32 a33 a11 a31 a12 a13 ~ ~ 0 a22 a23 a31 a31 a31 a11 a32 a11 a33 a13 a11 a11 a11 a11 a12 a13 0 ~ a22 ~ a23 0 ~ a32 ~ a33 Gaussian Elimination And finally the elimination of the a32 coefficient 1 0 0 a11 0 1 ~ a32 ~ a22 0 0 1 0 0 a12 ~ a22 ~ a 32 a13 ~ a23 ~ a 33 a11 0 0 ~ a32 a12 ~ a22 ~ a32 ~ ~ a22 a 22 a11 a12 a13 0 ~ a22 0 0 ~ a23 ~ ~ a33 ~ a33 a13 ~ a23 ~ a32 ~ ~ a23 a 22 Gaussian Elimination This result can be generalized. The Gaussian elimination algorithm can always be represented by a product of single lower triangular matrices Tij(-c) with the appropriate coefficients in the lower triangular part of the matrix. These appropriate coefficients are nothing else than the coefficient used in the elimination algorithm. 1 0 0 1 00 0 1 ~ a32 ~ a 22 0 0 a31 a11 10 0 1 U 01 T32 1 a 21 a11 0 0 0 a11 a12 a13 a11 a12 a13 1 0 a 21 a 22 a 23 0 ~ a22 0 1 a31 a32 a33 0 0 ~ a23 ~ ~ a33 ~ a32 ~ T31 a22 a31 T21 a11 a21 A a11 An example 123 A 235 356 1 0 0 1 00 0 1 0 0 10 0 11 301 T 1 00123 1 210235 0 0 01356 A 0 = 2 1 0 U 3 1 2 What can we do with this ? • We can represent the Gauss algorithm with a product of matrixes representing the elimination steps. T1T2T3 A U • Using Matrix Algebra, we can easily invert the algorithm. Inverting the Gauss algorithm • From T1T2T3 A U , we can have : A T1T2T3 1 1 1 U 1 A T3 T2 T1 U A LU What is the inverse of Tij(-c)? It turns out that the computation of this inverse is very simple. Following theorem applies: Theorem: The inverse of Tij c is Tij 1 **The proof can be done directly by verification c Tij c What is the inverse of T ? 1 a12 a13 c 1 0 a21 a22 a23 0 0 0 a11 0 1 a31 a32 1 00 01 a12 a21 ca11 a22 ca12 a33 c10 0 a11 a31 1 1 a32 00 c10 0 01 a13 a23 ca13 a33 For our example 1 0 0 1 00 0 1 0 0 10 0 301 11 1 00123 1 210235 0 0 01356 T = 00 1 001 0 01 210 123 0 100 1 00 1 235 0 356 A A 0 = 01 3010 T-1 110 2 3 1 0 1 2 U 2 1 0 U 3 1 2 1 001 0 01 0 100 1 00 1 235 0 356 A 00 210 123 01 3010 110 T-1 = 123 1001 235 2100 356 3110 A = LU 2 1 0 U 2 1 0 3 1 2 3 1 2 1 21 T a31 T321 a11 a21 T311 a11 1 a 21 a11 0 ~ a32 ~ a22 00 1 001 0 0 0 a 0 1 31 a11 100 1 ~ a32 ~ a 22 0 10 010 1 Therefore: . 1 a 21 a11 a31 a11 0 0 a11 a12 a13 a11 a12 a13 1 0 0 ~ a 22 a 21 a 22 a 23 ~ a32 ~ a 22 1 0 0 ~ a 23 ~ ~ a33 a31 a32 a33 Note: be careful of the correct choice of the sign used in the coefficient of the matrix T 1 a 21 a11 a31 a11 0 0 1 0 ~ a32 ~ a 22 1 Example 2 A 5 1 6 11 1 4 23 2 5 1 1 0 0 2 6 11 1 3 1 0 0 21 0 4 23 2 5 4 0 1 2 1 Solution 100 A 2 5 1 010 6 11 1 001 4 Gauss elimination for column 1 23 1 A 002 3 100 2010 5 4 8 1 Gauss elimination for column 2 2 5 1 A 0 02 3 1 00 2 210 5 4 0 1 2 1 Crout’s method 0 0 1 u12 u13 l 21 l 22 l31 l32 0 l33 0 0 u 23 1 l11 a11 a21 a12 a22 a13 a23 a31 a32 a33 a11 a21 a12 a22 a13 a23 l11 l21 l11u12 l21u12 l22 a31 a32 a33 l31 l31u12 l32 1 0 l11u13 l21u13 l22u23 l31u13 l32u23 l33 Crout’s method We can find, therefore, the elements of the matrices [L] and [U] by equating the two above matrices: l11 l11u12 a11 ; l21 a21 ; l31 a12 , hence u12 a31 a12 l11 l21u12 l22 a22 , hence l22 l31u12 l32 a32 , hence l32 a12 a11 l11u13 a13 , hence u13 l21u13 l22u23 a22 l21u12 a32 l31u12 a13 l11 a23 , hence u 23 l31u13 l32u23 l33 a13 a11 a23 l21u13 l22 a33 , hence l33 a33 l31u13 l32u23 Crout’s method For a general n n matrix, you have to apply the following expressions to find the LU decomposition of a matrix [A]: j1 lij aij lik u kj ;i j; i=1,2,…,n k1 i1 aij lik u kj k1 uij ; i < j; j=2,3,…,n lii uii 1 ; i =1,2,…,n Solution of equations Now to solve our system of linear equations, we can express our initial system: Ax b Under the following form Ax LU x b To find the solution {x} , we define first a vector {z} z Ux Our initial system becomes, then: Lz b As [L] is a lower triangular matrix the {z} can be computed starting by z1 until zn. Then the value of {z} can be found using the equation Lz b As [U] is an upper triangular matrix, it is possible to compute {x} using a backward substitution process starting xn until x1. z Ux Solving equations with the LU decomposition In summary, following method can be applied to solve a system of linear equation • Step 1: Decompose the matrix A by LU decomposition: A=LU • Step 2: Compute z by forward substitution using Lz=b • Step 3: Compute x by backward substitution using Ux=z. An example x 2y z 3 2x y 2z 3 3x yz 6 1 2 1x 3 2 1 2y 3 31 1 z 6 Step 1: A=LU 1 2 1 1 0 2 1 2 2 1 00 7 10 3 31 1 3 01 2 3 0 1 0 2 Step 2: Compute z with Lz=b 1 0 2 1 0 z2 7 1 z3 3 3 0 z1 3 3 6 z1 z2 3 3 z3 4 Step 3: Compute x with Ux=z 1 0 0 2 1x 3 0 3 3 2z 0 y 4 x 3 y 1 z 2 ...
View Full Document

This note was uploaded on 07/11/2011 for the course ENGR 391 taught by Professor Hoidick during the Winter '09 term at Concordia Canada.

Ask a homework question - tutors are online