MATLAB Gauss-Seidel method - diverging solutions

Vote count:

0

I have to implement the Gauss-Seidel method into MATLAB, here is my code..

function [ x ] = gs( A, b, x0, k )
%UNTITLED Summary of this function goes here
%   Detailed explanation goes here
n=length(A);
D = diag(diag(A));
L = tril(A-D);
U = triu(A-D);
B = L+D;
x=x0;     %ADD COMMENTS!%

for j = 1:k;
    for i = 1:n;
        yk(i) = b(i) - U(i,1:n)*x;

        %Now we solve B*x = yk by forward substitution%

        x(i) = yk(i);
        for m = 1:i-1
            x(i) = x(i) - B(i,m)*x(m);
        end
        x(i) = x(i)/B(i,i);

    end
end

A is an nxn matrix, b and nx1 vector, x0 is a starting 'guess' solution, and k is the number of iterations. If i set the number of iterations to k=1, i get a solution that is very close to the actual solution x = inv(A)*b. However, when i increase k, the solutions seem to just tend to infinity. The question also states that I am not allowed to calculated the inverse of (L+D) directly, so I have to use forward substitution to then find x.

I'm not too sure where my error lies, but i know that my algorithm should converge to x=inv(A)*b and I cannot see where the error is.

asked 32 secs ago

user65972

101

MATLAB Gauss-Seidel method - diverging solutions