function demo_maehly(P)

presnost=0.001;
n=length(P)-1;
A=max(abs(P(2:n+1)));
x0=1+A/abs(P(1));

% nalezeni prvniho korene
DP=polyder(P);
x=x0;
Px=polyval(P,x);
znam=sign(Px);
DPx=polyval(DP,x);
while znam*Px>0
   x=x-2*Px/DPx;
   Px=polyval(P,x);
   DPx=polyval(DP,x);
end
xo=x;
x=x-Px/DPx;
while abs(xo-x)>presnost
   Px=polyval(P,x);
   DPx=polyval(DP,x);
   xo=x;
   x=x-Px/DPx;
end

koreny=x;
k_M=x;

%vypocet dalsich korenu delenim
Q=P;

for ii=2:n
  Q=deconv(Q,[1, -x]);
  DQ=polyder(Q);
  x=x-0.001;
  xo=x;
  Qx=polyval(Q,x);
  DQx=polyval(DQ,x);
  x=x-Qx/DQx;
  while abs(xo-x)>presnost
    Qx=polyval(Q,x);
    DQx=polyval(DQ,x);
    xo=x;
    x=x-Qx/DQx;
  end
  koreny=[koreny,x];
end

%vypocet dalsich korenu Maehlyovou metodou

x=k_M-100*presnost;
for ii=1:n-1
  Px=polyval(P,x);
  DPx=polyval(DP,x);
  xo=x;
  pom=x-k_M(1:ii);
  pom=sum(1./pom);
  x=x-Px/(DPx-Px*pom);
  while abs(xo-x)>presnost
    Px=polyval(P,x);
    DPx=polyval(DP,x);
    xo=x;
    pom=x-k_M(1:ii);
    pom=sum(1./pom);
    x=x-Px/(DPx-Px*pom);
  end
  k_M=[k_M, x];
  x=x-100*presnost;
end

display('Koreny pomoci deleni polynomu:')
koreny
display('Koreny pomoci Maehlyovy metody:')
k_M
display('Skutecne koreny:')
roots(P)
