Happy New Year everyone! :)
I'm writing a Gauss-Seidel function in Matlab and I'm encountering some problems.
The iteration must stop when we reach 6 decimal digits of precision. It means that the infinite norm (asked to use it) of x-xprevious must be less than 0.5*10^(-6).
Firstly, here's my function:
function [x] = ex1_3(A,b)
format long
sizeA=size(A,1);
x=zeros(sizeA,1);
%Just a check for the conditions of the Gauss-Seidel Method
for i=1:sizeA
sum=0;
for j=1:sizeA
if i~=j
sum=sum+A(i,j);
end
end
if A(i,i)<sum
fprintf('\nGauss-Seidel''s conditions not met!\n');
return
end
end
%Actual Gauss-Seidel Method
max_temp=10^(-6); %Pass first iteration
while max_temp>(0.5*10^(-6))
xprevious=x;
for i=1:sizeA
x(i,1)=b(i,1);
for j=1:sizeA
if i~=j
x(i,1)=x(i,1)-A(i,j)*x(j,1);
end
end
x(i,1)=x(i,1)/A(i,i);
end
x
%Calculating infinite norm of vector x-xprevious
temp=x-xprevious;
max_temp=temp(1,1);
for i=2:sizeA
if abs(temp(i,1))>max_temp
max_temp=abs(temp(i,1));
end
end
end
And now the problems! When I call the function for a 3x3 array, I think it works. However, when I call it for a 10x10 array x becomes Inf (I guess it's out of machine numbers limits). Is there anything I can do to prevent this, except for changing the infinite norm and the 6 decimal digits precision (I must use these two, because my tutor told me so) ?
In the array I use (which was given to me) the entries outside the diagonal are -1 and the ones on the diagonal are 3. b is like this b=[2;1;1;1;1;1;1;1;1;2] (for n=10)
Your condition of the Gauss-Seidel Method is not correct:
D=diag(diag(A));
L=-tril(A,-1);U=-triu(A,1);
B=(D-L)\U;
R = max(abs(eig(B)));
if R>=1
fprintf('\nGauss-Seidel''s conditions not met!\n');
return
end
R is called spectral radius of iterative matrix B. It has to be less than 1 that Gauss-Seidel converges. Actually the matrix A in your test case has the R=1.8092, thus Gauss-Seidel method won't converge.
Check this slide from page 18 for more details.
EDIT
According to @LutzL's comment, you may use Gershgorin circle theorem to estimate the eigenvalue rather than calculate them with computational cost.
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With