Der Aufruf des Programms lautet x=solve(a,b,c,d,x0,tol,nmax)
Es soll Ax=d lösen

A soll nicht als n x n Matrix gespeichert werden sonder nur die Hauptdiagonale a und die Nebendiagonalen c und b!
x0 ist der Startvektor x soll alle Iterierten x_k Spaltenweise enthalten

Es soll Abgebrochen werden wenn
norm(Ax_k-d)_2 + norm(x_k-x_(k-1))_2 < tol oder sobald nmax überschritten wird

 

function x=solve1(a,b,c,d,x0,tol,nmax,par)


[n,m]=size(a);
A=zeros(n,n);
[s,t]=size(b);
[k,l]=size(c)
b(n)=0
c(n)=0
A=diag(a(1:n));
B=diag(b(1:n-1),-1);
C=diag(c(1:n-1),+1);
A=A+B+C

A*x0

x(1)=a(1)*x0(1)+(b(1)*x0(2))
%x(2)=c(1)*x0(1)+a(2)*x0(2)+b(2)*x0(3)
x(n)=c(n-1)*x0(n-1)+a(n)*x0(n)
for i=2:n
    x(i)=c(i-1)*x0(i)+a(i)*x0(i+1)+b(i)*x0(i+2)
    size(x(i))
end
x(n)=c(n-1)*x0(n-1)+a(n)*x0(n)


%entweder x=[x,neuervektor]
%oder x(:,i)=neuervektor

switch (par)

    case 1
        D = diag(a) ;            
        L = diag(b(1:n-1),-1) ;                
        U = diag(c(1:n-1),1) ;
        x=[n]; 
        i = 2;                                   
        while  i<=nmax
            x0=x(1);
            x(i)=inv(D)*(L+U)*x(i-1)+inv(D)*b;
            if sqrt((mvmult(a,b,c,x0))^2-d^2)+sqrt(x(i)^2+x(i-1)^2)<tol
                break;       
            end
            i=i+1;
        end  
       
    case 2
        D = diag(a) ;            
        L = diag(b(1:n-1),-1) ;                
        U = diag(c(1:n-1),1) ;   
        x=[n] ;
        i = 2;                                   
        while  i<=nmax 
            x0=x(1)
            x(i)=inv(D-L)*(U)*x(i-1)+inv(D-L)*b;
            if sqrt((mvmult(a,b,c,x0))^2-d^2)+sqrt(x(i)^2+x(i-1)^2)<tol
                break;       
            end
            i=i+1;
        end  
end