function [x,k,kstep]=durand(p,x0,delta,maxsteps)

% [y,k]=DURAND(p,x0,eps,maxsteps) izvaja Durand-Kernerjevo 
% metodo za racunanje nicel polinoma p z zacetnimi priblizki x0
% Racunanje se konca, ko so vse razlike med starimi in novimi priblizki
% manjse od delta ali ko presezemo maksimalno stevilo korakov maxsteps.
% 
% Primer uporabe:
% p = poly([1 2 3 4 5])
% x0 = rand(5,1)
% [x,k,step] = durand(p,x0,1e-12,50)

n=length(x0); % stevilo priblizkov

if nargin<3, delta=1e-12; end
if nargin<4, maxsteps=50*n; end

kon = zeros(n,1);  % tu bomo oznacili elemente, ki so ze skonvergirali
stkon = 0; % stevilo ze skonvergiranih vrednosti
iter = 0; % stevilo iteracij posameznih elementov
x = x0;
kstep = zeros(n,1);
korak = 0;

alfa = p(1); % vodilni koeficient
while (korak<maxsteps) && (stkon<n)
   korak = korak + 1; 
   for j = 1:n
     if kon(j)==0
        iter = iter+1;
        kstep(j) = kstep(j)+1;
        f = polyval(p,x(j));
        if f ~= 0
          dif = 1;
          for k = 1:n
             if k~=j
                dif = dif*(x(j)-x(k));
              end
          end
          novo = x(j)-f/(alfa*dif);
       else
          novo = x(j);
       end
       razlika(j) = abs(novo-x(j))/abs(novo);
       x(j) = novo;
       if razlika(j)<delta
         stkon = stkon+1;
         kon(j) = 1;
       end
     end
   end
   fprintf('Korak %3d:\n',[korak])  
   for j=1:n
     if (imag(x(j))~=0)
        disp(sprintf('%3d: %15.15e %+15.15ei  %0.1e',[j,real(x(j)),imag(x(j)),razlika(j)]));  
     else
        disp(sprintf('%3d: %15.15e  %0.1e',[j,x(j),razlika(j)]));  
     end
   end
   pause(0.3)
end
k=iter;


   