function [reseni, chyba, iterace] = newtonova_metoda_G(funkce,derivace,x0,epsilon,max_iter)
%Implementace Newtonovy metody pro reseni nelinearni rovnice f(x)=0
%
%	vstupni parametry:
%	funkce			- odkaz na zadanou funkci
%	derivace		- odkaz na derivaci zadane funkce
%	x0			- pocatecni aproximace 
%	epsilon 		- nastaveni pozadovane presnosti 
%	max_iter		- nastaveni maximalniho poctu iteraci 
%
%	vystupni paramery:
%	reseni			- posledni aproximace reseni
%	chyba   		- rozdil dvou po sobe jdoucich aproximaci
%	iterace			- pocet iteracnich kroku metody

clc;
tisk=sprintf('Newtonova metoda pro reseni nelinearni rovnice f(x)=0');
disp(tisk);
tisk=sprintf('pro pocatecni aproximaci x0=%0.8g a zastavovaci podminku |x(k)-x(k-1)|<%0.8g',...
	     x0,epsilon);
disp(tisk);
disp(' ');
tisk=sprintf('Funkce f(x) je zadana v souboru %s.m takto:',funkce);
disp(tisk);
nazev=sprintf('%s.m',funkce);
type(nazev)
disp(' ');
tisk=sprintf('Derivace funkce f(x) je zadana v souboru %s.m takto:',derivace);
disp(tisk);
nazev=sprintf('%s.m',derivace);
type(nazev)

% Pouze pro potreby vykreslovani
  figure(1);
  clf;
  hold on;
  zoom on;

  pom=zeros(max_iter+2,1);
  pom(1)=x0;
  xmin=x0;
  xmax=x0;
  ymin=feval(funkce,x0);
  ymax=feval(funkce,x0);

  iterace   = 0;
  tmp2      = x0;
  x         = tmp2;
  chyba     = inf;
  
  while (iterace < max_iter) && (chyba > epsilon)
    tmp1    = tmp2;
    x       = tmp1;
    tmp2    = tmp1 - feval(funkce,x)/feval(derivace,x);
    chyba   = abs(tmp2 - tmp1);
    iterace = iterace + 1;
    pom(iterace+1)=tmp2;
    if xmin>tmp2, xmin=tmp2; end 
    if xmax<tmp2, xmax=tmp2; end 
    if ymin>feval(funkce,tmp2), ymin=feval(funkce,tmp2); end 
    if ymax<feval(funkce,tmp2), ymax=feval(funkce,tmp2); end 
  end
  rozsah_x=(xmax-xmin);   
  rozsah_y=(ymax-ymin);
  xmin=xmin-0.3*rozsah_x;
  xmax=xmax+0.3*rozsah_x;
  ymin=ymin-0.3*rozsah_y;
  ymax=ymax+0.3*rozsah_y;
  axis([xmin xmax ymin ymax]);

  
  plot([xmin xmax],[0 0],'k');
  plot([0 0],[ymin ymax],'k');
  i=0;
  krok=rozsah_x/1000;
  for x=xmin:krok:xmax
    i=i+1;
    y(i)=feval(funkce,x);
  end
  x=xmin:krok:xmax;
  plot(x,y,'b');
% konec pripravy obrazku  
  
  iterace   = 0;
  tmp2      = x0;
  x         = tmp2;
  chyba     = inf;

  plot(tmp2,feval(funkce,tmp2),'ro');
  plot([tmp2,tmp2],[0,feval(funkce,tmp2)],'k--');
  if pom(2)<pom(1), zn1=1; else zn1=-1; end
  plot([pom(2)-zn1*abs(pom(1)-pom(2))/5,pom(1)+zn1*abs(pom(1)-pom(2))/5],...
       [-feval(funkce,tmp2)/5,feval(funkce,tmp2)*6/5],'m');
  
  
  tisk=sprintf( ... 
  '| iterace |       x(k)       | dx(k)=x(k)-x(k-1) |      f(x(k))      | ');
  disp(tisk)
  delka=length(tisk);
  retezec='';
  
  if (delka>80), delka=80; end
   iter=0;
   while iter<delka
    retezec=['-',retezec];
    iter=iter+1;
   end
  disp(retezec);

  tisk=sprintf('| %7.0f | %16.8f |                   |', iterace, tmp2 );
  disp(tisk)
  pause;

  while (iterace < max_iter) && (chyba > epsilon)
    tmp1    = tmp2;
    x       = tmp1;
    tmp2    = tmp1 - feval(funkce,x)/feval(derivace,x);
    chyba   = abs(tmp2 - tmp1);
    iterace = iterace + 1;
    fprintf('| %7.0f | %16.8f | %17.8f | %17.8f |\n', iterace, tmp2, ...
                 chyba,feval(funkce,tmp2));
    plot(tmp2,feval(funkce,tmp2),'ro');
    plot([tmp2,tmp2],[0,feval(funkce,tmp2)],'k--');
    if pom(iterace+2)<pom(iterace+1), zn1=1; else zn1=-1; end
    if pom(iterace+2)~=0
       plot([pom(iterace+2)-zn1*abs(pom(iterace+1)-pom(iterace+2))/5,...
	     pom(iterace+1)+zn1*abs(pom(iterace+1)-pom(iterace+2))/5],...
            [-feval(funkce,tmp2)/5,feval(funkce,tmp2)*6/5],'m');
    end
    pause;  
  end
  reseni    = tmp2;








