Gauss-Seidel method is a popular iterative method of solving linear system of algebraic equations. It is applicable to any converging matrix with non-zero elements on diagonal. The method is named after two German mathematicians: Carl Friedrich Gauss and Philipp Ludwig von Seidel.

Gauss-Seidel is considered an improvement over Gauss Jacobi Method. In this method, just like any other iterative method, an approximate solution of the given equations is assumed, and iteration is done until the desired degree of accuracy is obtained.

In earlier tutorials, we’ve already gone through the C program and algorithm/flowchart for Gauss-Seidel method. Here, we’re going to write a program code for **Gauss-Seidel method in MATLAB**, discuss its theoretical background, and analyze the MATLAB program’s result with a numerical example.

## Gauss-Seidel Method Theory:

Let’s go through a brief theoretical/mathematical background of Gauss-Seidel method. The matrices, iterations, and the procedure explained below cover the basic guidelines to write the program code for Gauss-Seidel method in MATLAB.

Consider the following system of linear equations:

a_{11}x_{1 }+ a_{12}x_{2 }+ a_{13}x_{3} + a_{14}x_{4} + a_{15}x_{5 }+ a_{16}x_{6 ……. }+ a_{1n}x_{n }= b_{1 }a_{21}x_{1 }+ a_{22}x_{2 }+ a_{23}x_{3} + a_{24}x_{4} + a_{25}x_{5 }+ a_{26}x_{6 ……. }+ a_{2n}x_{n }= b_{2 }a_{31}x_{1 }+ a_{32}x_{2 }+ a_{33}x_{3} + a_{34}x_{4} + a_{35}x_{5 }+ a_{36}x_{6 ……. }+ a_{3n}x_{n }= b_{3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ….. . . . . . . . .. . . . . . . . . . . . . .. . . .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ….. . . . . . . . .. . . . . . . . . . . . . .. . . .. . . . }a_{n1}x_{1 }+ a_{n2}x_{2 }+ a_{n3}x_{3} + a_{n4}x_{4} + a_{n5}x_{5 }+ a_{n6}x_{6 ……. }+ a_{nn}x_{n }= b_{n}

where, a_{ij }represents the coefficient of unknown terms x_{i}.

The above equations can be presented in matrix form as follows:

Or simply, it can be written as: [A][X] = [B]

Now, decomposing the matrix A into its lower triangular component and upper triangular component, we get:

A = L x U

where,

Further, the system of linear equations can be expressed as:

L x X = B – UX —–(a)

In Gauss-Seidel method, the equation (a) is solved iteratively by solving the left hand value of x and then using previously found x on right hand side. Mathematically, the iteration process in Gauss-Seidel method can be expressed as:

X^{(k+1)} = L ^{-1}( B –UX^{(k) })

Applying forward substitution, the elements of X^{(k+1)} can be computed as follows:

The same procedure aforementioned is followed in the MATLAB program for this method. The process of iteration is continued till the values of unknowns are under the limit of desired tolerance.

## Gauss-Seidel Method in MATLAB:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 | % Gauss-Seidel Method in MATLAB function x = gauss_siedel( A ,B ) disp ( 'Enter the system of linear equations in the form of AX=B') %Inputting matrix A A = input ( 'Enter matrix A : \n') % check if the entered matrix is a square matrix [na , ma ] = size (A); if na ~= ma disp('ERROR: Matrix A must be a square matrix') return end % Inputting matrix B B = input ( 'Enter matrix B : ') % check if B is a column matrix [nb , mb ] = size (B); if nb ~= na || mb~=1 disp( 'ERROR: Matrix B must be a column matrix') return end % Separation of matrix A into lower triangular and upper triangular matrices % A = D + L + U D = diag(diag(A)); L = tril(A)- D; U = triu(A)- D % check for convergence condition for Gauss-Seidel method e= max(eig(-inv(D+L)*(U))); if abs(e) >= 1 disp ('Since the modulus of the largest Eigen value of iterative matrix is not less than 1') disp ('this process is not convergent.') return end % initial guess for X..? % default guess is [ 1 1 .... 1] r = input ( 'Any initial guess for X? (y/n): ','s'); switch r case 'y' % asking for initial guess X0 = input('Enter initial guess for X :\n') % check for initial guess [nx, mx] = size(X0); if nx ~= na || mx ~= 1 disp( 'ERROR: Check input') return end otherwise X0 = ones(na,1); end % allowable error in final answer t = input ( 'Enter the error allowed in final answer: '); tol = t*ones(na,1); k= 1; X( : , 1 ) = X0; err= 1000000000*rand(na,1);% initial error assumption for looping while sum(abs(err) >= tol) ~= zeros(na,1) X ( : ,k+ 1 ) = -inv(D+L)*(U)*X( : ,k) + inv(D+L)*B;% Gauss-Seidel formula err = X( :,k+1) - X( :, k);% finding error k = k + 1; end fprintf ('The final answer obtained after %g iterations is \n', k) X( : ,k) |

Also see,

Gauss-Seidel C Program

Gauss-Seidel Algorithm/Flowchart

In the above MATLAB program, a function, **x = gauss_siedel( A ,B )**, is initially defined. Here, A and B are the matrices generated with the coefficients used in the linear system of equations. The elements of A and B are input into the program following the basic syntax of MATLAB programming.

A and B are to be checked: A should be a square matrix and B must be a column matrix to satisfy the criteria of Gauss-Seidel method. Then, as explained in the theory, matrix A is split into its upper triangular and lower triangular parts to get the value of first iteration.

The value of variables obtained from the first iteration are used to start the second iteration, and the program keeps on iterating till the solution are in the desired limit of tolerance as provided by the user.

Here’s a sample output screen of the MATLAB program:

## Gauss-Seidel Method Example:

The above MATLAB program of Gauss-Seidel method in MATLAB is now solved here mathematically. The equations given are:

4x_{1} – x_{2} –x_{3} = 3

-2x_{1 }+ 6x_{2} + x_{3} = 9

-x_{1} + x_{2} – 7x_{3} = -6

In order to get the value of first iteration, express the given equations as follows:

4x_{1} – 0 –0 = 3

-2x_{1 }+ 6x_{2} + 0 = 9

-x_{1} + x_{2} – 7x_{3} = -6

From the first equation: *x*_{1} = 3/4 = 0.750

Substitute the value of x_{1 } in the second equation : *x*_{2} = [9 + 2(0.750)] / 6 = 1.750

Substitute the values of x_{1} and x_{2} in the third equation: *x*_{3} = [-6 + 0.750 − 1.750] / 7 = − 1.000

Thus, the result of first iteration is: ( 0.750, 1.750, -1.000 )

The whole iteration procedure that goes on in Gauss-Seidel method (and the above MATLAB program) is presented below:

where, k is the number of iteration.

The final solution obtained is (1.000, 2.000, -1.000).

If you have any questions regarding Gauss-Seidel method, its theory, or MATLAB program, drop them in the comments. You can find more Numerical methods tutorial using MATLAB here.