% Newtonova metoda 04

funkce='x + 3/(x^2)';       % ucelova funkce        
DER1='1-6/(x^3)';           % 1. derivace ucelove funkce
DER2='18/(x^4)';            % 2. derivace ucelove funkce
a=0.5; b=3;                 % krajni hodnoty intervalu I=[a,b] 
x=(a+b)/2;                  % 1. aproximace bude stred intervalu [a,b]

p = linspace(a,b);                   % vykreslovany interval
q = (p.^(1)) + 3*(p.^(-2));          % predpis funkce
plot(p,q,'r')                        % vykresleni funkce   
line([1.817,1.817],[2.426,3.026])    % vyznaceni skutecneho minima

fx1=eval(funkce)-0.1;
fx2=eval(funkce)+0.1;
line([x,x],[fx1,fx2])    % vyznaceni aproximace minima

i=0;                  % promenna urcena k pocitani iteraci 
chyba=1;              % "mira presnosti" - jak presny bude vysledek 

disp(' '); % vlozi prazdny radek na vystupni obrazovku

% vypis na obrazovku

disp(sprintf('Je zadana funkce f=%s a interval [a,b]=[%f,%f].',funkce,a,b));
disp(sprintf('Nasim cilem je najit minimum funkce f s presnosti 0,0001.'));

disp(' '); % vlozi prazdny radek na vystupni obrazovku


% zakladni WHILE-cyklus (cyklus jednotlivych kroku) konci ve chvili
% dosazeni pozadovane presnosti, ci pokud citac iteraci dosahne hodnoty 15

while ((chyba>0.0001) & (i<15))
    
   i=i+1;  % zvyseni poctu iteraci     
   
   
   DERx=eval(DER1);      % hodnota 1. derivace funkce f v bode x
   DERxx=eval(DER2);     % hodnota 2. derivace funkce f v bode x
     
   disp(sprintf('%d. krok: a=%f, b=%f, x=%f, df(x)=%f, ddf(x)=%f',i,a,b,x,DERx,DERxx));         
      
   xx=x;   % ulozeni stare hodnoty x (posledni aproximace)   

   x=xx - (DERx/DERxx); 
  
   chyba = abs(xx-x);     % chyba jako rozdil poslednich dvou aproximaci
   
   disp(sprintf('        => nova aproximace: x=%f, chyba=%f',x,chyba));     
  
  
   disp(' '); % vlozi prazdny radek na vystupni obrazovku
  
   fx1=eval(funkce)-0.1*(i+1);
   fx2=eval(funkce)+0.1*(i+1);
   line([x,x],[fx1,fx2])    % vyznaceni aproximace minima
 
   if (chyba==0) i=17;       % "trefili" jsme presne minimum
   end    
  keyboard;
end % konec while-cyklu (cyklu kroku vypoctu, radky 34-62)


% Vypis informaci o vypocitanem reseni v zavislosti na tom,
% jestli jsme nahodou "netrefili" minimum presne

if (i==17)
      disp(sprintf('Presne reseni: x=%f.',x));
      disp(sprintf('Hodnota minima: f(x)=%f.' ,eval(funkce)));
      
else  
   
disp(sprintf('Priblizne reseni: x=%f s presnosti %f.',x,chyba));
disp(sprintf('Hodnota minima: f(x)=%f. ',eval(funkce)));
   
end