## Welcome to the World of Modelling and Simulation

### What is Modelling?

This blog is all about system dynamics modelling and simulation applied in the engineering field, especially mechanical, electrical, and ... ### Gauss-Seidel Iteration Method to Solve System of Algebraic Equations

The following matrix represents a system of linear algebraic equations. In this problem, we are going to solve these equations by applying Gauss-Seidel iteration approach.

The following program solves system of linear algebraic equations iteratively with successive approximation by using most recent solution vectors. However, there is a condition to work for this program which is strictly diagonal dominance. And, the given matrix in question is not diagonally dominant although matrix switching operations are carried out. For this circumstance, convergence is not guaranteed and the program shows the problem while calculating error as it increases with no of iterations.

MATLAB Program for Gauss-Seidel Iterative Method

% Gauss-Seidel Iterative Algorithm

function x = gauss_seidal(a,b,x0,e,N)

% a - (nxn) matrix
% b - column vector of length n
% x0 - initial solution
% e - error definition
% N - no. of iterations

format long;

m=size(a,1);   % get number of rows in matrix a

n=length(b);   % get length of b

p=length(x0);  % get length of x0

q=length(e);   % get length of e

if (m ~= n)

error('a and b do not have the same number of rows')

end

y0=x0;

for i=1:m       % Interchanging row/column to make diagonally dominant matrix

[~,f]=max(abs(a(1:m,i:i)));

% disp(f);

a([i,f],1:n)=a([f,i],1:n);

end

for
% while switching rows. This loop just interchange between rows to avoid '0'

if a(t,t)==0

[~,g]=max(abs(a(1:m,t:t)));

a([t,g],1:n)=a([g,t],1:n);

end

end

disp(a);

c=1;

while(c<N)

for j=1:m   % Successive iteration, most recent solution used for next calculation

for k=1:p

Z(1,k)=(a(j,1:n)*y0(1:p,1));

end

Y(j:j,1)=sum(Z(1,k));

% disp(Y);

R(j:j,1)=b(j:j,1)-Y(j:j,1);

% disp (R);

x(j:j,1)=y0(j:j,1)+(R(j:j,1)/a(j:j,j:j));

y0(j:j,1)= x(j:j,1);

% disp (y0);
% disp (x);
end

if (abs(x(1:p,1)-x0(1:p,1)) < e(1:q,1)) % Error checking

break;

end

disp (x);
disp(abs(x(1:p,1)-x0(1:p,1)));

x0=y0;

c=c+1;

end

x=x(1:p,1); % Returns solution

end

Program Outputs:

>> gauss_seidal([1 -1 4 0 2 9; 0 5 -2 7 8 4; 1 0 5 7 3 -2; 6 -1 2 3 0 8; -4 2 0 5 -5 3; 0 7 -1 5 4 -2], [19;2;13;-7;-9;2], [0; 0;0; 0; 0; 0],[0.05;0.05;0.05;0.05;0.05;0.05], 3)  % function calling

6    -1     2     3     0     8    % rearranged matrix
0     7    -1     5     4    -2
1     0     5     7     3    -2
-4     2     0     5    -5     3
0     5    -2     7     8     4
1    -1     4     0     2     9

3.166666666666667   % first improved solutions
0.285714285714286
1.966666666666667
1.019047619047619
-1.703571428571429
-0.593386243386244

3.166666666666667   % error
0.285714285714286
1.966666666666667
1.019047619047619
1.703571428571429
0.593386243386244

2.840388007054674   % second improved solutions
0.642705971277400
1.390044091710759
-0.732311665406903
-0.241714380196523
-0.586047738025251

0.326278659611992   % error
0.356991685563114
0.576622574955908
1.751359284454522
1.461857048374905
0.007338505360992

ans =                 % final solutions

2.840388007054674
0.642705971277400
1.390044091710759
-0.732311665406903
-0.241714380196523
-0.586047738025251

Percentage of error for the six solutions: 10.12%125%29%173%85%  and  1.1%

These errors may be the results of having a matrix, which is not diagonally dominant and the solution diverges as iteration increases.