KUBNA1

function [x,fx,n,rez] = kubna1( fun, dfun, x1, x2, tol )
% Kubna metoda - trazenje minimuma funkcije jedne promenljive FUN
% aproksimacijom polinomom 3. reda;
% Pocetni interval je [X1, X3], TOL je trazena tacnost.
% Minimum je u X i ima vrednost FX posle N iteracija.

X = [ x1 x2 ]';
Y = [ [1 1]' X X.*X X.*X.*X;  [0 0]' [1 1]' 2*X 3*X.*X ];
F = [ feval( fun,  X);  feval( dfun, X) ];
abcd = Y \ F;  b = abcd(2); c = abcd(3); d = abcd(4);
D = 4*c*c - 12*b*d;
xa = ( -2*c - sqrt(D) ) / 6 / d;
xb = ( -2*c + sqrt(D) ) / 6 / d;
n = 0;
rez = [n X' F' abcd' D xa xb ];

if (xa>X(1)) & (xa<X(2))
  x = xa;
else
  if (xb>X(1)) & (xb<X(2))
    x = xb;
  else
    disp('Greska u racunanju')
    return
  end
end
if feval(dfun,x) < 0
  X(1) = x;
else
  X(2) = x;
end
rez = [rez X' abs( [1 x x*x x*x*x] * abcd - feval(fun,x))];

while abs( [1 x x*x x*x*x] * abcd - feval(fun,x)) > tol
  Y = [ [1 1]' X X.*X X.*X.*X;  [0 0]' [1 1]' 2*X 3*X.*X ];
  F = [ feval( fun,  X);  feval( dfun, X) ];
  abcd = Y \ F;  b = abcd(2); c = abcd(3); d = abcd(4);
  D = 4*c*c - 12*b*d;
  xa = ( -2*c - sqrt(D) ) / 6 / d;
  xb = ( -2*c + sqrt(D) ) / 6 / d;
  n = n + 1;
  pom=[n X' F' abcd' D xa xb ];

  if (xa>X(1)) & (xa<X(2))
    x = xa;
  else
    if (xb>X(1)) & (xb<X(2))
      x = xb;
    else
      disp('Greska u racunanju')
      break
    end
  end
  if feval(dfun,x) < 0
    X(1) = x;
  else
    X(2) = x;
  end

  rez = [rez; [pom X' abs( [1 x x*x x*x*x] * abcd - feval(fun,x))] ];
end
fx = feval(fun, x);