function [X,iter]=nestacionarna(A,b,P,maxiter,k,tol)
% A - matrica sistema
% b - vektor desne strane
% P - pocetno resenje
% maxiter - maksimalni broj iteracija
% k - koliko razlicitih parametara tau odredjujemo
% tol - tolerancija

% [X,iter]=nestacionarna(toeplitz(4:-1:1),[0.1 2 0.1 2]',[1 1 1 1]',100,10,1e-4)

% Hermiteova
if ~isequal(A,A')
    error('Mora da vazi A*=A')
end

% Pozitivno definitna
[R,p] = chol(A);
if p~=0
    error('Matrica nije pozitivno definitna');
end

M=max(abs(eig(A)));
m=min(abs(eig(A)));

for j=1:k
    tau(j)=2/(M+m+(M-m)*cos((2*j-1)/(2*k)*pi));
end

% Ukoliko zelimo da promenimo redosled tau-ova (prvi, poslednji, drugi,
% pretposlednji...)
j=1;
for i=1:k/2
    tau2(j)=tau(i);
    tau2(j+1)=tau(k-i+1);
    j=j+2;
end

% Ako zelimo drugaciji redosled, skinuti komentar sa narednog reda
%tau=tau2;

greska=[];

% Ako zelimo ciklicno da se vrte tau-ovi onda while petlja.
% Ako zelimo da uporedimo resenje dobijeno za fiksiranih k iteracija u
% zavisnosti od razlicitog rasporeda tau-ova onda for petlja
%while 1    
for i=1:maxiter
    r=A*P-b;
    greska=[greska norm(r,inf)];
    if norm(r,inf)<tol
        plot(1:i,greska)
        iter=i;
        greska
        return 
    end
    X=P-tau(mod(i-1,k)+1)*r; % ako je k<maxiter ciklicno vrtimo tau-ove dobijene na pocetku
    P=X;       
end
