UoG Numerical method (CEng2073) Worku.
M (Structural Engineer )
Chapter – four
Linear system of algebra equation
Consider a system of n-linear equation in n
unknowns
a 11 x 1 +a12 x 2 +a13 x 3 +...+a 1n x n =b1 where aij, i, j =1, 2, 3, …… n are the
a 21 x 1 +a 22 x 2 +a 23 x 3 +...+a2 n x n=b2 coefficients bi, i=1, 2, 3, …… n are known
on the right hand side values and xi; i = 1, 2,
a 31 x 1 +a 32 x 2 +a 33 x 3 +...+a3 n x n=b3 3 …. n are unknown to be determined.
...........................................................
............................................................
a n1 x 1 +an 2 x 2 +an 3 x 3 +...+ann x n =b n
In matrix notation, we write the system as AX=B.................(1)
[ ] {} {}
a11 a 12 . .. a1n x1 b1
A= a21 a 22 . .. a2n ; X= x 2 ; ∧ B= b2
. . . . . .
a n1 an 2 . .. a nn xn bn
Where
The matrix [A:B] obtained by appending the column B to the matrix A is called the augmented
matrix, that is
[ ]
a11 a 12 .. . a 1n : b1
a a 22 .. . a 2n : b2
[ A : B ]= 21
. . . . : .
an 1 an 2 .. . ann : bn
Remark: the system of equation 1 is consistent (has at least one solution) if rank (A) = rank (A:B)
= r.
a) if r = n (number of unknowns = Rank ) then, the system has a unique solution
b) if r < n then, the system has infinite solution
c) the system of equation 1 is in consistent (has no solution ), if rank of A¿ rank of (A:B)
Methods to solve system of linear equation
4.1. Matrix inversion
−1
From eq.(1) it can be written as AX=B⇒ X= A B that is A-1 the inverse of A.
1 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
Definition of inverse; Let A and B be two n × n matrices. Then matrix A is said to be invertible
and to have an associated inverse matrix B if AB = BA = I.
Where, I is identity matrix a square matrix, which the diagonal element is 1 and the other
element is zero.
Example-1. Solve the following equation using matrix inversion
2 x+ y−2 z=8
a , x+ y−z=2
−x− y+2 z=−4
It can be written in matrix form AX=B
[ ]{ } { } { } [ ] { }
−1
2 1 −2 x 8 x 2 1 −2 8
1 1 −1 y = 2 ⇒ y = 1 1 −1 2
−1 −1 2 z −4 z −1 −1 2 −4
[ ]
2 1 −2
A= 1 1 −1
Let determine the inverse of −1 −1 2
The augmented matrix
[ ] [ ]
2 1 −2 1 0
0 1 0 −1 1 −1 0
[ A : I ]= 1 1 −1 |0 0⃗
1 R 1=R 1−R 2 1 1 −1 |0 1 0
−1 −1 2 0 0
1 −1 −1 2 0 0 1
⃗R 2= R 2−R 1
¿
[ ][
1 0 −1 1 −1 0 ⃗R 3=R 3+ R 1
]
1 1 −1 | 0 1 0 1 0 −1 1 −1 0
−1 −1 2 0 0 1 ¿ 0 1 0 |−1 2 0
0 −1 1 1 −1 1
[ ][
R 3=R 3+ R 2
1 0 −1 1 −1 0
] [ ]
0 1 0 |−1 2 0 1 0 −1 1 −1 0 1 0 0 1 0 1 ¿
¿ 0 1 0 |−1 2 0⃗ R 1=R 3+ R 1 0 1 0 |−1 2 0
0 −1 1 1 −1 1
0 0 1 0 1 1 0 0 1 0 1 1
Thus the matrix A-1(inverse of matrix A)
[ ]
1 0 1
−1
A = −1 2 0
0 1 1
2 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
Let’s check,
[ ][ ] [ ][ ][ ]
2 1 −2 1 0 1 1 0 1 2 1 −2 1 0 0
−1 −1
AA = A A=I ⇒ 1 1 −1 −1 2 0 = −1 2 0 1 1 −1 = 0 1 0
−1 −1 2 0 1 1 0 1 1 −1 −1 2 0 0 1
The value of x, y, z is
{}[ ] { }[ ]{ } { }
−1
x 2 1 −2 8 1 0 1 8 12
y = 1 1 −1 2 = −1 2 0 2 = −4
z −1 −1 2 −4 0 1 1 −4 −2
Which is x=12 , y=−4 , z=−2
Exercise -1. Solve the following equation using matrix inversion (submitted)
x−z=12
a , −x+2 y=6
y+ z=−3
4.2. Gauss – seidle method
In this method, we use the updated values of x1, x2, x3, ….xi-1 in computing the values of xi.
We assume that aij¿ 0 for all i.
We write as
¿ 11 x 1 +a 12 x 2 +a13 x 3 =b1 a 11 x1 =b 1−(a12 x 2 +a13 x3 )
a 21 x 1 +a 22 x 2 +a 23 x 3 =b2 ⇒ a 22 x 2=b2 −( a21 x 1 +a 23 x 3 )
a 31 x 1 +a 32 x 2 +a 33 x 3 =b3 a 33 x 3 =b3 −(a31 x 1 +a32 x 2 )
The Gauss- seidle method is defined as
1
x ( k +1 )= (b −(a12 x k +a13 x k ))
1 a11 1 2 3
1
x ( k +1 )= (b2 −(a21 x ( k +1)+a23 x k ))
2 a22 1 3
1
x ( k +1)= (b 3−(a31 x ( k +1)+a32 x ( k +1) ))
3 a33 1 2
Where k= 0, 1, 2, 3…
Example 2: find a solution of a system of equation
3 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
45 x 1 +2 x 2 +3 x 3 =58
3 x 1 +22 x2 +32 x 3 =47
5 x 1 + x 2 +20 x 3 =67 Correct to 3-decimal places using Gauss- seidle method.
Solution: using G.S method
1
x ( k +1 )= (58−(2 x k +3 x k ))
1 45 2 3
1
x ( k +1 )= (47−(3 x (k +1 ) +32 x 3 k ))
2 22 1
1
x ( k +1)= 20 (67−(5 x (k +1) +x (k +1) ))
3 1 2 , k = 0,1 2, 3 …
0
Starting with x i =0 , i=1 ,2 , 3 .. ..
For this example, initialize with x2 = x3 = 0
Then the successive iteration becomes
k x1(k+1) x2(k+1) x3(k+1)
0 1.288889 1.960606 2.929747
1 1.006434 1.732782 3.011752
2 1.011093 1.724692 3.010992
3 1.011503 1.724705 3.010889
Exercise 2: find a solution of a system of equation
4 x1 −2 x 2 +3 x 3=50
3 x 1−22 x 2 +32 x 3=47
52 x1 +x 2−20 x 3 =6 Correct to 3-decimal places using Gauss-seidle method.
4.3. Gauss-elimination
The method is about based on the idea of reducing the given system of equations AX=b to upper
– triangular system of equation UX=Z using elementary operation, we know that this two system
are equivalent that is the solution of both system UX = Z is then solved by back substitution to
obtain the solution vector X.
We illustrate the method using the nxn system
Consider a system A(1)x = b(1)
4 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
[ ]{ } { }
a ( 1) a (1) ... a (1 ) x1 b ( 1)
11 12 1n 1
a a ... a x2 b ( 1)
21
( 1)
22
(1)
2n
(1 ) = 2
. . . . . .
a a ... a xn b ( 1)
(1 ) ( 1) ( 1)
n1 n2 nn n
a ( 1) ≠0
Step-1- 11 then defined a multiplier
a ( 1)
i1
M ( 1) = a11
where i=2 ,3 , 4 ,. .. . ..
i1
Then, we can eliminate X1 from the last (n-1) equation by
substituting the ith equations the Mi1 of the first equation. The first raw of A(1) and B(1) are unchanged
and are called the first pivot.
The new system becomes A(2)x = b(2)
[ ]{ } { }
a ( 1) a ( 1) .. . a ( 1) x1 b (1 )
11 12 1n 1
0 a .. . a x2 b (2 )
22
( 2)
2n
( 2) = 2
. . . . . .
0 a .. . a xn b (2 )
(2 ) (2 )
n2 nn n
When the new coefficients given by
a ( 2) =a ( 1) −M ( 1) a ( 1) i=1 , 2 ,3 , .. . .. ..
ij ij i1 1j
b ( 2)=b ( 1)−M ( 1) b ( 1)
i i i1 1
a ( 1)≠0
Step-2- 22 follow a similar way, eliminate X2 from the last (n-2) equation and we get A(3)x =b(3)
[ ]{ } { }
a ( 1) a ( 1) ... a ( 1) x1 b ( 1)
11 12 1n 1
0 a ... a x2 b
22
( 2)
2n
( 2) = 2( 2)
. . . . . .
0 0 ... a xn b ( 3)
( 3)
nn n
Where the new coefficient are given by
a ( 3 )=a (2 )−M ( 2) a (2 ) i=1 ,2 , 3 ,. . ..∧ j=1 , 2 ,3 , .. . .. .
ij ij i2 2j
b ( 3)=b ( 2)−M ( 2) b (2 )
i i i2 2 Where
a (2 )
i2
M ( 2) = a where i=3 , 4 ,. .. . ..
i2 ( 2)
22
5 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
By continuing in this manner after n-1 steps
[ ]{ } { }
a ( 1) a ( 1) ... a ( 1) b (1 )
11 12 1n x1 1
0 a ( 2) ... a ( 2) x2 b (2 )
22 2n = 2
. . a (3) a ( 3)
. .
33 3n
xn b (n )
0 0 0 .. . a (n ) n
nn
n
∑ aik x k
k =i+1
x i=bi− i=n−1 ,n−2 ,. . .. .. 1
Using back substitution a ij
Remark: when any one of the pivot is zero or its very small number as elimination progress the
gauss elimination fail.
If the pivot is zero then division by zero gives over flow error (not defined)
If the pivot is very small number then division by it introduce large around of error and a
solution may contain large error.
Example: 3. solve the following system of linear equation using Gauss elimination method
x 1−x 2 + 2 x 3 −x 4 =8
x 1 +2 x 2 +x 3 =0 2 x1 −2 x 2 +3 x 3−3 x 4 =20
2 x1 +2 x 2 +3 x 3 =3 x 1 + x 2 + x 3 =−2
a.
−x 1−2 x 2 =2 b. x 1−x 2 + 4 x 3 + 3 x 4 =4
Solution
The matrix form is written as AX=B
[ ]{ } { }
1 2 1 x1 0
2 2 3 x2 = 3
−1 −3 0 x2 2
a ( 1) a ( 1)
21 2 31 −1
M 21= = =2 ; M 31= = =−1
Step-1, a 11≠0 so a11 1 a11 1
Then the augmented matrix
[ ] [ ]
1 2 1 : 0 1 2 1 : 0
2 2 3 : 3⃗ R 2= R 2−2 R 1 ∧ R 3 =R 3+ R 1 0 −2 1 : 3
−1 −3 0 : 2 0 −1 1 : 2
6 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
a ( 1)
32 −1 1
M 32= = =
a22 −2 2
[ ]
1 2 1 : 0
[ ]
1 2 1 : 0
⃗ 1 0 −2 1 : 3
0 −2 1 : 3 R 3= R3 − 2 R 2
1 1
0 −1 1 : 2 0 0 :
2 2
The n th matrix form of UX = B’
[ ]{ } { }
1 2 1 0 x 1 +2 x 2 + x 3 =0
x1
0 −2 1 3 −2 x 2 + x 3=3
x2 = ⇒
0 0 1 1 x3 1
2 x2 2 =
2 2 By back substitution we gate x1=1, x2=-1, x3=1.
x 1−x 2 +2 x 3 −x 4 =−8
2 x1 −2 x 2 +3 x 3−3 x 4 =−20
x 1 +x 2 +x 3 =−2
b. x 1−x 2 +4 x 3 +3 x 4 =4
Solution
[ ] [ ]
1 −1 2 −1 : −8 ⃗ R 2=R 2 −2 R 1
1 −1 2 −1 : −8
2 −2 3 −3 : −20 ⃗ 0 0 −1 −1 : −4
R =R −R
1 1 1 0 : −2 3 3 1 0 2 −1 1 : 6
[A:B]=
1 −1 4 3 : 4 ⃗ R4 =R4 −R1 0 0 2 4 : 12
The procedure cannot be continued in its present form, since the pivot element is zero
(a ( 2)=0)
a ( 2) and a ( 2)
22 so exchange raw 3 and raw 2 in searching element 32 42 for the
first non-zero element.
[ ][ ] [ ]
1 −1 2 −1 : −8 1 −1 2 −1 : −8 1 −1 2 −1 : −8
0 0 −1 −1 : −4 0 2 −1 1 : 6 0 2 −1 1 : 6
⃗ ⃗
R 4 =R 4 +2 R 3
0 2 −1 1 : 6 0 0 −1 −1 : −4 0 0 −1 −1 : −4
0 0 2 4 : 12 0 0 2 4 : 12 0 0 0 2 : 4
Then 1
x =−4 , x =2 , x =0 , x 4 =2
2 3
Exercise 3: solve a system of equation
x 1−x 2 + 2 x 3 −x 4 =8
x 1 +10 x 2−x 3 =7 2 x1 −2 x 2 +3 x 3−3 x 4 =20
2 x1 +4 x 2 +3 x 3 =9 x 1 + x 2 + x 3 =−2
a.
10 x 1−x 2 +x 3 =4 b. x 1−x 2 + 4 x 3 + 3 x 4 =4
7 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
4.4. LU decomposition method (factorized method)
This method is based on the fact a square matrix a can be factored in to the form LU where
L unit lower triangular matrix (elements of the diagonal are one and above the diagonal are zero)
U unit upper triangular matrix (elements of the diagonal are one and below the diagonal are zero)
Provided all principal minor of A are none zero. If the factorization of A is exist then it is unique.
a 11 x 1 +a12 x 2 +a13 x 3 =b 1
a 21 x 1 +a 22 x 2 +a 23 x 3 =b2
Consider
a 31 x 1 +a 32 x 2 +a 33 x 3 =b3
It can be write as in matrix form AX=B...............................(1) where
[ ] {} {}
a11 a12 a13 x1 b1
A= a 21 a22 a23 X = x2 B= b 2
a 31 a32 a33 x b3
and
Let A=LU .........................(2) Where L is lower triangular matrix
[ ] [ ]
1 0 0 U 11 U 12 U 13
L= L21 1 0 U = 0 U 22 U 23
L31 L32 1 0 0 U 33
and the equation becomes LUX=b . .. .. .(3 ) & we write
V 1 =b 1
LV =b ⇔ L21 V 1 +V 2 =b2
1
UX=V .......( 4 ) Then eq. (3) becomes L31 V 1 + L32 V 2 +V 3 =b3
Can be solved for v1, v2, v3 by forward substitution. Then after v is determined UX =V gives
u11 x 1 +u12 x 2 +u13 x 3 =v 1
u22 x2 +u23 x 3 =v 2
u33 x 3=v 3 Which can be solved for x1, x2,& x3 by back substitution
Remark-1
To compute the matrix L and U we write equation 2 as
[ ][ ][ ]
1 0 0 U 11 U 12 U 13 a11 a12 a13
L21 1 0 0 U 22 U 23 = a21 a22 a23
L31 L32 1 0 0 U 33 a31 a32 a33
By multiplying LU we get
8 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
i. U 11 =a 11 ; U 12=a 12 ; U 13=a13
a a a a
ii . L21 U 11 =a21 ⇒ L21= 21 = 21 iii. L31 U 11=a31 ⇒ L31= 31 = 31
U 11 a 11 U 11 a11
a21
iv . L21 U 12 +U 22=a22 ⇒U 22=a22−L21 U 12=a 22− a12
a11
a21
v . L21 U 13 +U 23 =a23 ⇒U 23=a23−L21 U 13=a 23− a13
a11
1 1 a
vi. L31 U 12 + L32 U 22=a32 ⇒ L32= (a32−L31 U 12 )= ( a32− 31 a 12 )
U 22 U 22 a 11
vii. L31 U 13+ L32 U 23 +U 33=a32 ⇒U 33=a 32−( L31 U 13 + L32 U 23 )=a32−( L31 U 13+ L32 U 23 )
To compute the element of L and U in the following set order
1. First raw of U i.e. U11, U12, U13
2. The first column of L i.e. L21, L31
3. 2nd raw of U
4. 2nd columns of L
5. 3rd raw of U
Example 4. Solve the following equation by Appling LU decomposition
3 x+2 y+7 z=4
2 x+3 y+z=5
3 x+4 y+z=7
Solution.
Step-1, From eq. (2) LU=A
Let
[ ][ ][ ]
1 0 0 U 11 U 12 U 13 3 2 7
L21 1 0 0 U 22 U 23 = 2 3 1
L31 L32 1 0 0 U 33 3 4 1
Implies
9 2012 E.C
UoG Numerical method (CEng2073) Worku.M (Structural Engineer )
i . U 11 =3 ; U 12=2 ; U 13=7
2 2 3 3
ii . L21 U 11 =2⇒ L21= = iii . L31 U 11 =3 ⇒ L31= = =1
U 11 3 U 11 3
2 5
iv . L21 U 12 +U 22=3 ⇒U 22=3−L21 U 12=3− 2=
3 3
2 −11
v . L21 U 13 +U 23 =1⇒ U 23=1−L21 U 13=1− 7=
3 3
1 1 6
vi. L31 U 12 + L32 U 22=4 ⇒ L32= ( 4−L31 U 12 )= ( 4−1∗2 )=
U 22 5 5
3
6
∗−11
5 −8
vii. L31 U 13+ L32 U 23 +U 33=1 ⇒U 33 =1−( L31 U 13+ L32 U 23 )=1−(1∗7+ )=
3 5
Thus
[ ][ ]
1 0 0 3 2 7
2 5 −11
1 0 0
A= 3 3 3
6 −8
1 1 0 0
5 5
Step-2. LV=B
[ ]{ } { }
1 0 0 v1 = 4
2 v1 4 2 7
1 0 v +v 2 =5 ⇒ v 2=
LV =B ⇒ 3 v2 = 5 ⇒ 3 1 3
6 v3 7 6 1
1 1 v 1+ v + v 3 =7⇒ v 3 =
5 5 2 5
Step-3. Hence the original system becomes UX=V
[ ]{ } { }
3 2 7 4 3 x +2 y+ 7 z=4
5 −11 x 7 5 11 7
0 y− z=
UX =V ⇔ 3 3 y = 3 ⇒ 3 3 3
−8 z 1 −8 1
0 0 z=
5 5 5 5
by back substitution
−1 9 7
z= ; y= ; ∧ x=
8 8 8
Exercise 4. Solve the following equation by using LU decomposition method
x +3 y +6 z=2 3 x + y +2 z=3
a , 3 x− y +4 z =9 b , x−3 y +z=−3
x−4 y+2 z=7 x −2 y +z =−4
10 2012 E.C