Numerical methods - lecture 6
Jiří Zelinka
Autumn 2017
Jiří Zelinka
1/22
Repetition
Iteration methods for solving system of linear equations
Ax = b —>    x = 7~x + g Iteration process:
x/c+i =Txk + g = 0,1,...
Solution:
x = (E-T)~1g
2/22
Theorem
The sequence {x^}^_0 determined by the iteration process x = Tx + g converges for every initial iteration x° g Rn       p(T) < 1. In this case
lim x^ = x. x = 7~x + g
k—>oo
Jacobi iteration method
/-th equation:
a/ixi H-----h a//X/ H-----h a/nX,, = 6/
The component X; is expressed as k-th iteration:
□ s1
Jiří Zelinka Numerical methods - lecture 6
Matrix notation
Ax = b,      A = D + L+U,
Ax = (D + L + U)x = b
(an 0 \
\ 0 ann J
Jiří Zelinka
□ iS1 Numerical methods - lecture 6
L =
í o
■
V anl ■ ( O 3i2
U =
Vo
a„,„-i  O /
3ln \
3/7—1, n
O
= -D"1(/.+ Ľ)x+ D_ib.
-i
■k+i _ n-i
= -D~l(L+ U)xk + D^b.
□ iS1
Jiří Zelinka Numerical methods - lecture 6
5/22
x^1 = Tjxk + D~1b:
Tj = -D~\L + U), ty = -f for i^j, ta = 0
3,7
a
Tj =
0
^21
^22
V
3nl
]nn
3l2
0
3n2
]nn
3l„ \
au
32 n 322
0
D_1b =
au
Ó2_ 322
V ann J
Jiří Zelinka
□ [31 Numerical methods - lecture 6
6/22
Gauss-Seidel iteration method
/c+l 1 / i k k k
*1 = - [Dl- 3l2X2 ~ ai3X3 - a14*4 ~
all
/c+l 1 /c+l /c k
+      ~~ (^2 ~~ ^21*1+   ~~ a23*3 ~~ ^24*4
x2
*3
322
*+1  "    1 (fc - a3ixr+1 - a32x2fc+1 - a34x4fc
^33
i-l n
7=1 y=/+l
Jiří Zelinka
□ [31 Numerical methods - lecture 6
{D + L+U)x = {D + L)x = x =
Ux + b
(D + Z.)"1 Ux + (D + /.)-
TG = (D + Z.)"1^,      xfc+1 = T^x* + (D + L)
k+l _
Theorem: If A is diagonally dominant matrix, i.e.
>
E
'U
or
>
E
"JI
then Jacobi and Gauss-Seidel methods converge.
□ [31
Jiří Zelinka Numerical methods - lecture 6
Relaxation (Succesive over-relaxation (SOR)) method
xk - k-th iteration
x^1 - the following iteration aquired by the Gauss-Seidel metod
u e (0,2) - relaxation parameter
xk+1 = (1 - u)xk + uxkGf
Jiří Zelinka
Numerical methods - lecture 6
9/22
Direct methods for solving system of linear equations
Ax = b, A - non-singular Gaussian elimination method:
reduction [A\b] to the system with upper triangular matrix [U\b].
Operations:
• Row exchange
• Adding c-multiple of the i-th row to the j-th row
□ S
Jiří Zelinka Numerical methods - lecture 6
10/22
Corresponding matrices:
/ 1 0 0 1
P* =
1
0 0
0 0 0
0
0 0 0 0 0 1
0
1 0 0 0
0 0
1
0 0
0 0
1 0 0 0
0 0 0 0 0 1
0
0 0 0
0 0
1
P-hk - permutation matrixm Pik = P-k = P,^
T
Jiří Zelinka
Numerical methods - lecture 6
/ 1 o
O 1
Li,k,c —
0\ O
1
1
V o
1/
Jiří Zelinka
12/22
/ 1 o
O 1
Vo
1
—c
1
o
1/
Jiří Zelinka
13/22
Gaussian elimination method without row exchange
(1) zeroing the first column under the diagonal
A™ I b(1)) = U ■ {A I b)
Li =
( 1 • hi 1
V Inl
0
lkl = ~
3kl
an
0
1
□ iS1
Jiří Zelinka Numerical methods - lecture 6
(i) zeroing the i-th column under diagonal
(A® I bW) = L; ■ I b^'"1))
í I
0
1
'n,
Iki = "
a(,'-1} (/-I) ?
4
/ = 2,..., n — 1
□ iS1
Jiří Zelinka Numerical methods - lecture 6
15/22
U I b) =Ln_1-...L2-L1-(A\b)
SO
U = Ln_ľ • ... L2- U - A
then
1—-^ /_2 • • • _LJ — A.
-i
Matrices /_/ are lower triangular so matrices /_11 are lower triangular, too. Then
A = L-U,      L = L^Lž1 ...C-i-
-i
/_ - lower triangular matrix.
□ [31
Jiří Zelinka Numerical methods - lecture 6
16/22
Jiří Zelinka
17/22
LU decomposition
The product A = L • U is called LU (LR) decomposition of the matrix A
Solving linear equations:
substitution Ux = y, and we get system Ly = b with the lower triangular matrix, then we solve Ux = y with the upper triangular matrix.
LU decomposition with row exchange
P-A = L-U for suitable permutation matrix P. We must do row exchange if a,-,- = 0.
Pradcically, we find out the row exchanges during calculation.
18/22
Pivoting (partial)
To improve the numerical stability the element with maximal absolut value is chosen in each column in the rest of the matrix to be modify. Then we exchange appropriate rows.
Example
2xi + 4x2 ~  x3 — —5
xi +  x2 — 3x3 = -9
4xi +  x2 + 2x3 = 9
Jiří Zelinka
Numerical methods - lecture 6
19/22
/ 2 4  -1 \ / 1  0 0 \
4=1  1-3,/.=     * 10, \ 4  1    2 J \ *  *  1 /
pivot
The first line is multiplied by —^ and added to the second line. The first line is multiplied by — | and added to the third line.
20/22
1 2
L= [  h   1  0 I, A-+ I 0  ?   -I  \      [ 0  I -2
3 _7
4 2
pivot
L-> I i  1 0 |,        I 1
The second line is multiplied by —^ and added to the third line.
' 4  1    2   \ / 1   0 0
\ 0  I   -2   \=U,L^\ \ 10
Jiří Zelinka Numerical methods - lecture 6
Theorem
If all main minors of the matrix A are non-zero then it is possible to do Gaussian elimination method without row exchange and the LU decomposition has form
A = L-U
Jiří Zelinka
Numerical methods - lecture 6
22/22