function y = laguerre(p,x0,verzija,delta,maxsteps)

% y=LAGUERRE(p,x0,verzija,delta,maxsteps) dela Laguerrovo iteracijo za 
% računanje ničle polinoma p z zacetnim priblizkom x0. Verzija je lahko 1 
% (stabilna), 2 (fiksno v desno), 3 (fiskno v levo), 
% ce jo ne podamo je privzeto 1.
% 
% iteracija se konca, ko se zadnja dva priblizka razlikujeta za manj kot 
% delta ali ko je presezeno stevilo iteracij maxsteps. Ce vrednosti ne 
% podamo v klicu metode, so privzete vrednosti delta=1e-12 in maxsteps=50.
%
% Primer uporabe
% p=poly([1 2 3 4 5])
% laguerre(p,2.2,1,1e-12,50)
% laguerre(p,5.2,2,1e-12,50)
% laguerre(p,5.2,3,1e-12,50)

% BP 2004/2005

if nargin<3, verzija=1; end
if nargin<4, delta=1e-8; end
if nargin<5, maxsteps=50; end

p1 = polyder(p);
p2 = polyder(p1);
m = length(p)-1;
n = m; 

xn = x0; x0=Inf; korak=0;                    						 

while (abs(xn-x0)>delta) && (korak<maxsteps)
   
   x0 = xn;
   f = polyval(p,x0); f1 = polyval(p1,x0); f2 = polyval(p2,x0);
   
   predznak = sign(f1); 
   if predznak==0, predznak=1; end 
   imen1 = f1+predznak*sqrt((n-1)*((n-1)*f1^2-n*f*f2));
   imen2 = f1-predznak*sqrt((n-1)*((n-1)*f1^2-n*f*f2));
   xn = x0-n*f/imen1;
   xn1 = x0-n*f/imen2;
   disp(sprintf('                                    x1=%15.15f   x2=%15.15f',[xn,xn1]))
   if verzija==1
      if abs(imen2)>abs(imen1), xn=xn1; end
   end
   if verzija==2
      if xn<x0, xn=xn1; end
   end
   if verzija==3
      if xn>x0, xn=xn1; end
   end
   
   korak = korak+1;           						 
   if (imag(xn)~=0)
      disp(sprintf('%3d: %15.15e %+15.15ei  %0.1e',[korak,real(xn),imag(xn),abs(xn-x0)]));  
   else
      disp(sprintf('%3d: %15.15e  %0.1e',[korak,xn,abs(xn-x0)]));  
   end
   pause(0.2)
end   

y=xn;