function stpdsc(nrf, mu, err_goal, niter_max)
% Parametry wejsciowe:
% nrf       - numer minimalizowanej funkcji
% mu        - wspolczynnik kroku 
% err_goal  - dokladnosc osiagniecia minimum
% niter_max - maksymalna liczba krokow
mc = 0.0;		% STALA ODFILTROWYWANIA OSCYLACJI
switch nrf
   case 1
      % FUNKCJA TESTUJACA NR 1
      xmin = 0; xmax = 5; ymin = 0; ymax = 5;
   case 2
      % FUNKCJA TESTUJACA NR 2
      xmin = 0.5; xmax = 2; ymin = 0.0; ymax = 1.5;
   case 3   
      % FUNKCJA TESTUJACA NR 3
      xmin = 0.5; xmax = 1.5; ymin = 0.5; ymax = 1.5;
   case 4
      % FUNKCJA TESTUJACA NR 4
      xmin = -2.5; xmax = 2.5; ymin = -1.5; ymax = 1.5;
   otherwise
      error('Nieprawidlowy wybor funkcji');  
end      

% RYSOWANIE POWIERZCHNI MINIMALIZOWANEJ FUNKCJI
xrange = linspace(xmin, xmax, 50);
yrange = linspace(ymin, ymax, 50);
[X, Y] = meshgrid(xrange, yrange);
Z = fun(X, Y, nrf);
mesh(xrange, yrange, Z)
xlabel('x'); ylabel('y'); zlabel('z=f(x,y)')
title('Powierzchnia minimalizowanej funkcji')
pause


contour(xrange, yrange, Z, 50)
axis equal
xlabel('x'); ylabel('y');
title('Wykres poziomicowy minimalizowanej funkcji')
hold on

% OKRESLENIE PRZYBLIZENIA POCZATKOWEGO
disp('Wybierz mysza przyblizenie poczatkowe...')
[x, y] = ginput(1)
xconv = x; yconv = y; fconv = fun(x, y, nrf);
[gx, gy] = grad(x, y, nrf);
dx = 0; dy = 0;
for niter = 1: niter_max
      % SPRAWDZ KRYTERIUM STOPU
      if norm([gx gy], inf) < err_goal
            niter = niter - 1;
            break
         end

      x_old = x; y_old = y;

      % WYZNACZ KROK W KIERUNKU NAJSZYBSZEGO SPADKU
      dx = mc * dx - (1 - mc) * mu * gx; dy = mc * dy - (1 - mc) * mu * gy;       
      x = x_old + dx; y = y_old + dy;
      [gx, gy] = grad(x, y, nrf);
      xconv = [xconv x]; yconv = [yconv y]; fconv = [fconv fun(x, y, nrf)];

      % WYSWIETL POPRAWIONE ROZWIAZANIE
      plot( [x_old x], [y_old y], 'r-', 'EraseMode', 'none')
      pause2(1.0)
      drawnow
   end
hold off   

% WYSWIETL REZULTATY
if niter < niter_max
      disp('Minimalizacja zakonczyla sie sukcesem !')
      disp(sprintf('(x_opt, y_opt) = (%6.3f, %6.3f)', x, y))
      disp(sprintf('grad_opt = (%g, %g)', gx, gy))
      disp(sprintf('Minimalna wartoscia funkcji jest %g', fun(x, y, nrf)))
   else
      disp('Przekroczono maksymalna liczbe krokow')
      disp(sprintf('(x_last, y_last) = (%6.3f, %6.3f)', x, y))
      disp(sprintf('grad_last = (%g, %g)', gx, gy))
      disp(sprintf('Ostatnio otrzymana wartoscia funkcji jest %g', fun(x, y, nrf)))
   end
pause

% WYKRES ZBIEZNOSCI WARTOSCI FUNKCJI
plot(0:length(fconv)-1, fconv);
xlabel('Liczba krokow')
ylabel('Wartosc funkcji')
title('Ewolucja wartosci funkcji') 
pause

plot(0:length(xconv)-1, xconv);
xlabel('Liczba krokow')
ylabel('Wspolrzedna x')
title('Ewolucja wspolrzednej x') 
pause

plot(0:length(yconv)-1, yconv);
xlabel('Liczba krokow')
ylabel('Wspolrzedna y')
title('Ewolucja wspolrzednej y') 
pause

close all


function Z = fun(X, Y, nrf)
switch nrf
  case 1
     % FUNKCJA TESTOWA NR 1
     Z = (X - 2).^2 + 2 * (Y - 3).^2 + 3;
  case 2
     % FUNKCJA TESTOWA NR 2 (Engvall) Z MINIMUM W PUNKCIE (x, y) = (1.31, 0.95)
     Z = X.^4 + Y.^4 - 2 * X.^2 .* Y - 4 * X + 3;
  case 3
     % FUNKCJA TESTOWA NR 3 (Rosenbrock)	
     Z = (X - 1).^2 + 10 * (X - Y.^2).^2;
  case 4
     % FUNKCJA TESTOWA NR 4 (2-D six-hump camel back function) Z MINIMUM GLOBALNYM
     % W PUNKTACH (x, y) = (-0.0988, 0.7126) ORAZ (0.0898, -0.7126)
     Z = (4 - 2.1 * X.^2 + X.^4 / 3) .* X.^2 + X .* Y + (-4 + 4 * Y.^2) .* Y.^2;
  otherwise
     error('Nieprawidlowy wybor funkcji');  
end


function [gx, gy] = grad(x, y, nrf)
switch nrf
   case 1
      gx = 2 * (x - 2);
      gy = 4 * (y - 3);
   case 2
      gx = 4 * x^3 - 4 * x * y - 4;
      gy = 4 * y^3 - 2 * x^2;
   case 3
      p = 10;
      gx = 2 * x - 2 + 2 * p * (x - y^2);
      gy = -4 * p * (x - y^2) * y;
   case 4
      gx = (-4.2 * x + 4.0 / 3.0 * x^3) * x^2 + 2 * (4 - 2.1 * x^2 + x^4 / 3) * x + y;
      gy = x + 8 * y^3 + 2 * (-4 + 4 * y^2) * y;
  otherwise
     error('Nieprawidlowy wybor funkcji');  
end
   


   