The Mathematical Garden

Iterative Method

If we are to solve

Ax = é ê ê ê ê ê ê ê ê ê ë  1 2  1 3 0  1 3 1  1 3 0  1 3  1 2 ù ú ú ú ú ú ú ú ú ú û é ê ê ê ë x1 x2 x3 ù ú ú ú û = é ê ê ê ë 5 10 5 ù ú ú ú û =b.

There are at least three different ways of splitting matrix A. They result in three types of iterative methods.

A = é ê ê ê ë 1 0 0 0 1 0 0 0 1 ù ú ú ú û + é ê ê ê ê ê ê ê ê ê ë  -1 2  1 3 0  1 3 0  1 3 0  1 3  -1 2 ù ú ú ú ú ú ú ú ú ú û     case 1 = é ê ê ê ê ê ê ê ë  1 2 0 0 0 1 0 0 0  1 2 ù ú ú ú ú ú ú ú û + é ê ê ê ê ê ê ê ê ê ë 0  1 3 0  1 3 0  1 3 0  1 3 0 ù ú ú ú ú ú ú ú ú ú û     case 2 = é ê ê ê ê ê ê ê ê ê ë  1 2 0 0  1 3 1 0 0  1 3  1 2 ù ú ú ú ú ú ú ú ú ú û + é ê ê ê ê ê ê ê ë 0  1 3 0 0 0  1 3 0 0 0 ù ú ú ú ú ú ú ú û     case3 = S + (A-S)

Now

Ax=(S + (A-S)) x = b

and therefore

S x + (A-S) x = b.

Hence we may write

x = S-1 b -S-1 (A - S) x

where we assume that S^-1 exists.

Given an initial guess x⁽⁰⁾ of the solution of Ax=b, one may consider the following iterative scheme:

x(k+1) = S-1 b - S-1 (A-S) x(k)

Clearly if x^(k) ® x as k ® ¥ then we have x = A^-1 b.

Theorem 2: If there exists a matrix ||.||_M such that ||S^-1 (A-S)||_M < 1 then the iterative scheme converges to the solution of Ax=b. The proof can be found in the 1st item in the reference section.

Remark: A matrix norm ||.||_M is similar to a vector norm. For any two square matrices A and B of the same size, we have (i) ||A||_M=0 if and only if A=0, (ii) ||cA||_M = |c|||A||_M and (iii) ||A+B||_M £ ||A||_M + ||B||_M.

Examples

Let A be the matrix and b be the right hand side vector and we use x⁽⁰⁾=[0, 0, 0]^t as the initial guess.

Case 1: When S=I, this is called the Richardson Method.

x(k+1) = b- (A-I) x(k) = é ê ê ê ë 5 10 5 ù ú ú ú û - é ê ê ê ê ê ê ê ê ê ë -  1 2  1 3 0  1 3 0  1 3 0  1 3 -  1 2 ù ú ú ú ú ú ú ú ú ú û x(k)

x(1) = [ 5, 10, 5]t x(2) = [4.1667, 6.6667, 4.1667]t x(3) = [4.8611, 7.2222, 4.8611]t x(4) = [5.0231, 6.7593, 5.0231]t : x(30) = [5.9983, 6.0014, 5.9983]t.

Case 2: When S = Diag [a₁₁,¼, a_nn] then this is called the Jacobi method.

Therefore

x(k+1) = S-1 b- S-1 (A-S) x(k) = é ê ê ê ë 10 10 10 ù ú ú ú û - é ê ê ê ê ê ê ê ë  1 2 1  1 2 ù ú ú ú ú ú ú ú û -1   é ê ê ê ê ê ê ê ê ê ë 0  1 3 0  1 3 0  1 3 0  1 3 0 ù ú ú ú ú ú ú ú ú ú û x(k) = [10 10 10]t - é ê ê ê ê ê ê ê ê ê ë 0  2 3 0  1 3 0  1 3 0  2 3 0 ù ú ú ú ú ú ú ú ú ú û x(k)

x(1) = [ 10, 10, 10]t x(2) = [3.3333, 3.3333, 3.3333]t x(3) = [7.7778, 7.7778, 7.7778]t : x(30) = [6.0000, 6.0000, 6.0000]t.

Case 3: When S = Lower triangular part of the matrix A. This method is called the Gauss-Seidel Method.

x(k+1) = S-1 b- S-1 (A-S) x(k) = é ê ê ê ê ê ê ê ë 10  20 3  50 9 ù ú ú ú ú ú ú ú û - é ê ê ê ê ê ê ê ê ê ë  1 2 0 0  1 3 1 0 0  1 3  1 2 ù ú ú ú ú ú ú ú ú ú û -1   é ê ê ê ê ê ê ê ë 0  1 3 0 0 0  1 3 0 0 0 ù ú ú ú ú ú ú ú û x(k)

x(1) = [ 10,  20 3 ,  50 9 ]t x(2) = [5.5556, 6.2963, 5.8025]t x(3) = [5.8025, 6.1317, 5.9122]t x(4) = [5.9122, 6.0585, 5.9610]t : x(14) = [6.0000, 6.0000, 6.0000]t.

Theorem 3: If A is diagonally dominant, then both the Jacobi method and Gauss-Seidel method converge to the solution of Ax=b. The proof can be found in the book listed in the reference section.